removing large 3D radiation arrays from restarts #428

rgknox · 2018-10-17T01:28:13Z

Description:

This PR is a refactor that removes the need to remember large 3D arrays in restart files. These arrays were storing the fraction of absorbed radiation, as dimensioned by crown layer x leaf layer x pft. These variables were forcing the restart files to maintain very large cohort dataspaces which were acting as surrogates to hold this data.

For example, with 14 PFTs, 4 crown layers and 30 leaf layers, this would necessitate space for 1680 values per each column. Since we were using the cohort dimension to hold these, this forces massive restart files.

The refactor instead saves 4 smaller variables related to ground albedo and solar zenith angle, and then reproduces these 3D arrays from existing information, immediately after the restart information is read in.

Fixes: #12 !

Collaborators:

@jenniferholm has identified that runs with multiple PFTs are really slow, which motivated this work.

Expectation of Answer Changes:

There should not be answer changes.

Checklist:

My change requires a change to the documentation.
I have updated the in-code documentation .AND. (the technical note .OR. the wiki) accordingly.
I have read the CONTRIBUTING document.
FATES PASS/FAIL regression tests were run
If answers were expected to change, evaluation was performed and provided

Test Results:

CTSM (or) E3SM (specify which) test hash-tag:

CTSM (or) E3SM (specify which) baseline hash-tag:

FATES baseline hash-tag:

Test Output:

All PASS: cheyenne, intel, fates test suite:

/gpfs/fs1/scratch/rgknox/clmed-tests/fates.cheyenne.intel.C46895e9-F3c87f2a-pft-arrays-vf3

rosiealice

That makes sense to me. I always thought that these things had made it into the restart file as some sort of unpleasant consequence of the timing of the albedo calls relative to everything else, so it's a good idea to find a work around for that.

… the 3D radiation patch variables, thus decreasing restart file sizes.

rgknox · 2018-10-26T19:17:31Z

exactly, the albedo calls come before the calls to Norman radiation, so the initial values needed for the albedo calculations are not available. I wonder if cold-starts have un-initialized values for that first albedo call, I should check that...

rgknox · 2018-10-26T19:42:20Z

For cold-starts, I don't believe we go through the process of "perfectly" initializing the canopy albedo for the first time-step. We initialize all the absorption fractions as 0, which I'm guessing makes the vegetation perfectly reflective for that first time step.

Part of me thinks this is not-a-big-deal because:

this only happens for the very first of the short-time steps.
This is is fine if the model is started at night (but half of the planet doesn't start at night).

ckoven · 2018-10-26T21:17:56Z

main/EDTypesMod.F90

@@ -329,6 +329,13 @@ module EDTypesMod
     integer  ::  ncan(nclmax,maxpft)                           ! number of total   leaf layers for each canopy layer and pft

     !RADIATION FLUXES      
+
+     integer  ::  solar_zenith_flag                           ! integer flag specifying daylight (based on zenith angle)
+     real(r8) ::  solar_zenith_angle                          ! solar zenith angle (radians)


why are solar_zenith_flag and solar_zenith_angle defined as patch variables? Couldn't that info be held at the site level?

it could. I was just following the precedent of the flag that is passed in from the HLM.

ckoven · 2018-10-26T21:28:34Z

biogeophys/EDSurfaceAlbedoMod.F90

@@ -195,756 +125,918 @@ subroutine ED_Norman_Radiation (nsites, sites, bc_in, bc_out )
                    bc_out(s)%fabd_parb(ifp,:) = 0.0_r8
                    bc_out(s)%fabi_parb(ifp,:) = 0.0_r8
                    do ib = 1,hlm_numSWb
+
+                       ! Requires a fix here, abld vs albi


so are we going to fix this here, or in https://github.com/NGEET/fates/pull/431/files#diff-68b652d39e00b970fbcf3d5b8d6f1165L199 ?

i was going to fix this, but didn't realize @jkshuman got to it first until I saw that PR. I think it depends on which one gets integrated first. There will be a conflict between the two either way (so they both will?)
For testing purposes, I will first run tests maintaining the status quo to test b4b. Then I will implement that fix.

add in the fix since b4b testing is complete.

ckoven · 2018-10-26T21:30:59Z

biogeophys/EDSurfaceAlbedoMod.F90

+    real(r8) :: angle
+
+    real(r8),parameter :: tolerance = 0.000000001_r8
+    real(r8), parameter :: pi   = 3.141592654                 ! PI


use FatesConstantsMod, only : pi => pi_const?

yeah, that would be better, will do

ckoven · 2018-10-26T21:31:45Z

biogeophys/EDSurfaceAlbedoMod.F90

+      !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++!
+      cosz = max(0.001_r8, currentPatch%solar_zenith_angle ) !copied from previous radiation code...
+      do ft = 1,numpft
+         sb = (90._r8 - (acos(cosz)*180/pi)) * (pi / 180._r8)


Suggested change

sb = (90._r8 - (acos(cosz)*180/pi)) * (pi / 180._r8)

sb = (90._r8 - (acos(cosz)*180._r8/pi)) * (pi / 180._r8)

got that change in my current working version, I'm currently working through a bug so will issue more changes.

ckoven · 2018-10-26T21:32:28Z

biogeophys/EDSurfaceAlbedoMod.F90

+                  tr_dif_z(L,ft,:) = 0._r8
+                  do iv = 1,currentPatch%nrad(L,ft)
+                     do j = 1,9
+                        angle = (5._r8 + dble(j - 1) * 10._r8) * 3.142_r8 / 180._r8


what's 3.142_r8?

this is code is all copied and pasted, so I didn't do too much of a pass on numerics and constants and such.
I essentially just pulled all the code in the norman radiation routine that was inside the patch loop, and made it a subroutine. The only things that are different, are how the variables forc_dir and forc_dif are constructed, and the bc_ variables are processed before the patch level subroutine.

yeah i know, just giving a hard time. but if it is okay not to be b4b, no reason not to clean this stuff up while we're looking at it.

ok, good point, I will do a pass through targeting b4b, and then after I'm satisfied, I will do a cleaning pass

ckoven · 2018-10-26T21:33:45Z

biogeophys/EDSurfaceAlbedoMod.F90

+               if (L == currentPatch%NCL_p.and.currentPatch%NCL_p > 1)then !is this the (incomplete) understorey?
+                  !Add on the radiation coming up through the canopy gaps.
+                  !diffuse to diffuse
+                  weighted_dif_up(L) = weighted_dif_up(L) +(1.0-sum(ftweight(L,1:numpft,1))) * &


Suggested change

weighted_dif_up(L) = weighted_dif_up(L) +(1.0-sum(ftweight(L,1:numpft,1))) * &

weighted_dif_up(L) = weighted_dif_up(L) +(1.0_r8-sum(ftweight(L,1:numpft,1))) * &

ckoven · 2018-10-26T21:33:58Z

biogeophys/EDSurfaceAlbedoMod.F90

+                       weighted_dif_down(L-1) * currentPatch%gnd_alb_dif(ib) 
+                  !direct to diffuse
+                  weighted_dif_up(L) = weighted_dif_up(L) + forc_dir(radtype) * &
+                       weighted_dir_tr(L-1) * (1.0-sum(ftweight(L,1:numpft,1))) * currentPatch%gnd_alb_dir(ib)


Suggested change

weighted_dir_tr(L-1) * (1.0-sum(ftweight(L,1:numpft,1))) * currentPatch%gnd_alb_dir(ib)

weighted_dir_tr(L-1) * (1.0_r8-sum(ftweight(L,1:numpft,1))) * currentPatch%gnd_alb_dir(ib)

ckoven · 2018-10-26T21:34:28Z

biogeophys/EDSurfaceAlbedoMod.F90

+                        if (L == 1)then !is this the top layer?
+                           Dif_dn(L,ft,1) = forc_dif(radtype)
+                        else
+                           Dif_dn(L,ft,1) = weighted_di


Suggested change

Dif_dn(L,ft,1) = weighted_di

(1.0_r8-sum(ftweight(L,1:numpft,1)))

that suggestion got messed up somehow i think? beta...

think I got that in, done

ckoven · 2018-10-26T21:37:15Z

biogeophys/EDSurfaceAlbedoMod.F90

+                        if (L == 1)then !is this the top layer?
+                           Dif_dn(L,ft,1) = forc_dif(radtype)
+                        else
+                           Dif_dn(L,ft,1) = weighted_di


Suggested change

Dif_dn(L,ft,1) = weighted_di

if ( abs(error) > 0.0001_r8)then

another suggestion got messed up?

actually should this .0001 be a constant, and used elsewhere?, e.g. https://github.com/NGEET/fates/pull/428/files#diff-68b652d39e00b970fbcf3d5b8d6f1165R941

I would rather punt on that number for now and keeping scope limited

ckoven · 2018-10-26T22:07:08Z

main/FatesRestartInterfaceMod.F90

@@ -1416,6 +1410,7 @@ subroutine create_patchcohort_structure(this, nc, nsites, sites, bc_in)
     use EDTypesMod,           only : ed_cohort_type
     use EDTypesMod,           only : ed_patch_type
     use EDTypesMod,           only : ncwd
+     use EDTypesMod,           only : maxSWb
     use EDTypesMod,           only : nlevleaf


Suggested change

use EDTypesMod, only : nlevleaf

note -- I'm suggesting we pull nlevleaf out of the restart infrastructure entirely as it should no longer be needed in this module. is that correct?

also there's a call to it in

fates/main/FatesRestartInterfaceMod.F90

Line 1031 in 44f66d6

use EDTypesMod, only : nlevleaf

which i can't suggest removal of in the inline editor. could you delete it there too?

ckoven · 2018-10-26T22:07:42Z

main/FatesRestartInterfaceMod.F90

@@ -1613,6 +1608,7 @@ subroutine get_restart_vectors(this, nc, nsites, sites)
     use EDTypesMod, only : ed_cohort_type
     use EDTypesMod, only : ed_patch_type
     use EDTypesMod, only : ncwd
+     use EDTypesMod, only : maxSWb
     use EDTypesMod, only : nlevleaf


Suggested change

use EDTypesMod, only : nlevleaf

ckoven · 2018-10-26T22:09:26Z

main/FatesRestartInterfaceMod.F90


-             if ( debug ) write(fates_log(),*) 'CLTV io_idx_pa_sunz 2 ',io_idx_pa_sunz
+             if ( debug ) write(fates_log(),*) 'CLTV 1186 ',nlevleaf,numpft,nclmax


Suggested change

if ( debug ) write(fates_log(),*) 'CLTV 1186 ',nlevleaf,numpft,nclmax

if ( debug ) write(fates_log(),*) 'CLTV 1186 ',numpft,nclmax

removed refs

rosiealice · 2018-10-30T09:21:14Z

main/FatesRestartInterfaceMod.F90

@@ -1058,7 +1061,7 @@ subroutine set_restart_vectors(this,nc,nsites,sites)
    integer  :: io_idx_co      ! cohort index
    integer  :: io_idx_pa_pft  ! each pft within each patch (pa_pft)
    integer  :: io_idx_pa_cwd  ! each cwd class within each patch (pa_cwd)
-    integer  :: io_idx_pa_sunz ! index for the combined dimensions for radiation
+    integer  :: io_idx_pa_sb   ! each SW band (vis/ir) per patch (pa_sb)


I think the wavelength band has a terminology (IB?) in CTSM already, so maybe we should mirror that? Altternatively, I think '_SWb' is more understandable than '_sb'.

updated naming suffix to _ib

rgknox · 2018-11-03T04:48:46Z

I had to revert some changes I made to the precision of constants and some math that had marginally changed answers, but after reverting the changes, this PR is now B4B with master.
fates test suite, cheyenne, intel
base ctsm: 27f7fa fates: af76d76
/glade/u/home/rgknox/ctsm/clmed-tests/fates.cheyenne.intel.C46895e9-Fd003357-pft-arrays-v400
I will next do a pass to clean up precision on some of the constants and variables.

rgknox · 2018-11-04T02:58:10Z

all tests pass, will integrate early next week unless I hear other comments or requests.

ckoven

this looks great, thanks @rgknox!

rgknox added the PR status: Not Ready The author is signaling that this PR is a work in progress and not ready for integration. label Oct 17, 2018

rosiealice reviewed Oct 18, 2018

View reviewed changes

Added a wrapper call to norman radiation, which allows us to re-build…

ba1a087

… the 3D radiation patch variables, thus decreasing restart file sizes.

rgknox force-pushed the pft-arrays branch from 7685441 to ba1a087 Compare October 24, 2018 23:19

rgknox mentioned this pull request Oct 26, 2018

fates refactor of 3d arrays in restart - api changes ESCOMP/CTSM#548

Merged

Small bugfix initializing output absorption fraction arrays

44f66d6

ckoven reviewed Oct 26, 2018

View reviewed changes

changed zenith flag to logical

dd41af4

rosiealice reviewed Oct 30, 2018

View reviewed changes

rosiealice approved these changes Oct 30, 2018

View reviewed changes

Reverting some precisions to original to attempt closer b4b results.

5d0148f

rgknox added 3 commits November 3, 2018 17:26

Updated precision on constants in norman radiation routine.

3440cfb

updated naming suffix for band indices in restart -> _ib

e0495bb

Removed references ti nlevleaf and nclmax in the restart code

1f0db48

rgknox added the PR status: Ready This PR is ready for integration after having passed through peer review and final testing. label Nov 4, 2018

rgknox removed the PR status: Not Ready The author is signaling that this PR is a work in progress and not ready for integration. label Nov 4, 2018

ckoven approved these changes Nov 5, 2018

View reviewed changes

rgknox merged commit b54091d into NGEET:master Nov 6, 2018

mvdebolskiy mentioned this pull request Oct 3, 2023

failure in ERS 3 month nocomp test #1051

Closed

rgknox deleted the pft-arrays branch October 31, 2023 18:52

	sb = (90._r8 - (acos(cosz)180/pi)) (pi / 180._r8)
	sb = (90._r8 - (acos(cosz)180._r8/pi)) (pi / 180._r8)

	weighted_dif_up(L) = weighted_dif_up(L) +(1.0-sum(ftweight(L,1:numpft,1))) * &
	weighted_dif_up(L) = weighted_dif_up(L) +(1.0_r8-sum(ftweight(L,1:numpft,1))) * &

	weighted_dir_tr(L-1) * (1.0-sum(ftweight(L,1:numpft,1))) * currentPatch%gnd_alb_dir(ib)
	weighted_dir_tr(L-1) * (1.0_r8-sum(ftweight(L,1:numpft,1))) * currentPatch%gnd_alb_dir(ib)

	Dif_dn(L,ft,1) = weighted_di
	(1.0_r8-sum(ftweight(L,1:numpft,1)))

	Dif_dn(L,ft,1) = weighted_di
	if ( abs(error) > 0.0001_r8)then


		if ( debug ) write(fates_log(),*) 'CLTV io_idx_pa_sunz 2 ',io_idx_pa_sunz
		if ( debug ) write(fates_log(),*) 'CLTV 1186 ',nlevleaf,numpft,nclmax

removing large 3D radiation arrays from restarts #428

removing large 3D radiation arrays from restarts #428

Conversation

rgknox commented Oct 17, 2018 • edited Loading

Description:

Collaborators:

Expectation of Answer Changes:

Checklist:

Test Results:

rosiealice left a comment

Choose a reason for hiding this comment

rgknox commented Oct 26, 2018

rgknox commented Oct 26, 2018

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

rgknox commented Nov 3, 2018

rgknox commented Nov 4, 2018

ckoven left a comment

Choose a reason for hiding this comment

rgknox commented Oct 17, 2018 •

edited

Loading