migrate FATES to use the normal BGC call sequence

[rgknox]

Description of changes

FATES previously used its own reduced form of the BGC call sequence. These changes get FATES to interact with the same BGC call sequence used in CN. This step is necessary for nutrient dynamics and wood product fluxes. This will also likely help in the maintenance of the FATES and BGC calling code, as we won't have to manage each seperately.

Should be synchronized (api change) with FATES PR 999

Specific notes

Grid-scale mass balances for FATES are bypassed.

The cnveg_state and cnveg_flux data structures are not allocated when FATES is active. Two new filters have been added bgc_soilc and bgc_vegp, these are the filters that are passed to the BGC sequence for columns and patches respectively. The bgc_vegp filter currently will have no elements in it if FATES is active, but will have the same patch list as the soilp filter. The bgc_soilc filter is the same as the soilc filter, but is zero for an SP run or a FATES SP run.

Contributors other than yourself, if any:

@adrifoster @billsacks @ckoven @ekluzek @glemieux @wwieder

CTSM Issues Fixed (include github issue #):

May largely take care of #1879, but we need to add testing for that and verify
Fixes #2112

Are answers expected to change (and if so in what way)?

Non-FATES are B4B, other than the TOTAL C/N fields which were moved from BGC to Soil-BGC and as such are
roundoff different (because of change in ordering of the sums). FATES runs... I don't expect B4B with base, but I do expect qualitatively similar results with a base.

Any User Interface Changes (namelist or namelist defaults changes)? pending further discussion, FATES will default to using supplemental nitrogen in c-only mode, and no supplemental N when PRT2 is active (ie N cycling on).

Testing performed, if any:

Single site smoke tests.

Regular testing: on cheyenne and izumi now completed

[rgknox]

Here is a comparison of a single site FATES simulation with the the new call-sequence versus the old call sequence:

https://drive.google.com/file/d/1na1j_-4ZLu0dqBc3EhYqeKlCTV9L_mug/view?usp=sharing

Its a tropical site, BCI, with one tropical broadleaf evergreen. With the CLM BGC call sequence, the heterotrophic respiration is higher, leading to diminished SOM and litter pools. Both have accelerated spin-up turned off.

[rgknox]

I realized that in the previous version I had only enabled the BGC calls, but many of them still had if-clauses that would exit if it was a FATES run. I went through and removed those. Following, I ran a simulation comparing the new coding under two contexts, with and without supplemental N.

TLDR: What I am seeing makes sense, with some caveats. Firstly, we just shouldn't run a carbon-only FATES simulation without supplemental N turned on, the litter carbon just builds up because there is no N to accompany it, and therefore it just can't decompose. Second, N deposition still needs to be enabled. Third, the initial NO3 is super high, not sure why that is. Maybe it is a way of kick-starting a spin-up so that there is no need of an initial N supplementation period?

In the first figure I compare against the base version of the model.

And in the second figure, I just show the two new versions (the base model has no N-cycle, so its trivial).

[rgknox]

Added N deposition when FATES is on. Also, cleaned up use_fates filters by including a new use_fates_bgc filter. This is mutually exclusive with use_fates_sp, they both can't be true. They are both false though if use_fates is false:

These plots show that ndep is applied to the soil N balance in the new version, yet total N is conserved. The version without N deposition relied on a higher supplementation flux to achieve the same immobilization rates.

[rgknox]

@samsrabin thanks for mentioning. Yeah, I think I modified the logic on the n-deposition clause just above this block of code, and I had not merged in the PR that added this block with the cropcal_interp yet. it should be ok to just accept the HEAD in this case.

[ekluzek]

@rgknox can you write some notes as a comment in terms of things that you've noticed and figured out about the restart problem? That could be helpful in preparation of the brainstorming session tomorrow. Sorry I didn't think of it until now. I think there are notes in other places that I'll copy to here as well...

[ekluzek]

FATES restart tests are showing roundoff level differences in a few fields such as FATES_HET_RESP and FATES_NEP.

for example the test:

ERS_Ld60.f45_f45_mg37.I2000Clm50FatesCruRsGs.cheyenne_intel.clm-FatesColdPPhys

Notes from the Jul/27th CTSM-Software meeting:

But FATES tests are still not passing muster. The issue is not around bit-for-bit with baseline (no way to get around answer changes). But the issue is with reproducibility of restarts. The differences are incredibly small.
Visually it looks like most diffs are in Antarctica, but the diffs are actually for the vegetated landunit.
Weirdly, the diffs don’t appear until multiple days into the run.
Could it be that diffs first appear in fields that aren’t output (or maybe output at too-low precision)?
Yes, maybe: many variables don’t appear in output... but Ryan is now turning on a bunch of output variables. He will double-check that this works prior to this branch.
It seems that the differences originate in the soil.
Could it be something related to running means?
The one thing Ryan thought of with that was something he addressed by setting it to 0 temporarily.

[ekluzek]

In CLM the field for total heterotrophic respiration term is hr_col in SoilBGCCarbonFlux.

The sequence for FATES is currently in clm_driver as...

          call clm_fates%WrapUpdateFatesRmean(nc,temperature_inst)

          call EDBGCDyn(bounds_clump,                                                              &
               filter(nc)%num_soilc, filter(nc)%soilc,                                             &
               filter(nc)%num_soilp, filter(nc)%soilp,                                             &
               filter(nc)%num_pcropp, filter(nc)%pcropp, doalb,                                    &
               bgc_vegetation_inst%cnveg_carbonflux_inst, &
               bgc_vegetation_inst%cnveg_carbonstate_inst, &
               soilbiogeochem_carbonflux_inst, soilbiogeochem_carbonstate_inst,                    &
               soilbiogeochem_state_inst, clm_fates,                                               &
               soilbiogeochem_nitrogenflux_inst, soilbiogeochem_nitrogenstate_inst,                &
               c13_soilbiogeochem_carbonstate_inst, c13_soilbiogeochem_carbonflux_inst,            &
               c14_soilbiogeochem_carbonstate_inst, c14_soilbiogeochem_carbonflux_inst,            &
               active_layer_inst, atm2lnd_inst, water_inst%waterfluxbulk_inst,                     &
               canopystate_inst, soilstate_inst, temperature_inst, crop_inst, ch4_inst)

          if ( decomp_method /= no_soil_decomp )then
             call EDBGCDynSummary(bounds_clump,                                             &
                   filter(nc)%num_soilc, filter(nc)%soilc,                                  &
                   filter(nc)%num_soilp, filter(nc)%soilp,                                  &
                   soilbiogeochem_carbonflux_inst, soilbiogeochem_carbonstate_inst,         &
                   c13_soilbiogeochem_carbonflux_inst, c13_soilbiogeochem_carbonstate_inst, &
                   c14_soilbiogeochem_carbonflux_inst, c14_soilbiogeochem_carbonstate_inst, &
                   soilbiogeochem_nitrogenflux_inst, soilbiogeochem_nitrogenstate_inst,     &
                   nc)
          end if

          call clm_fates%wrap_update_hifrq_hist(bounds_clump, &
               soilbiogeochem_carbonflux_inst, &
               soilbiogeochem_carbonstate_inst)


          if( is_beg_curr_day() ) then

             ! --------------------------------------------------------------------------
             ! This is the main call to FATES dynamics
             ! --------------------------------------------------------------------------

             if ( masterproc ) then
                write(iulog,*)  'clm: calling FATES model ', get_nstep()
             end if
             call clm_fates%dynamics_driv( nc, bounds_clump,                        &
                  atm2lnd_inst, soilstate_inst, temperature_inst, active_layer_inst, &
                  water_inst%waterstatebulk_inst, water_inst%waterdiagnosticbulk_inst, &
                  water_inst%wateratm2lndbulk_inst, canopystate_inst, soilbiogeochem_carbonflux_inst, &
                  frictionvel_inst, soil_water_retention_curve)

             ! TODO(wjs, 2016-04-01) I think this setFilters call should be replaced by a
             ! call to reweight_wrapup, if it's needed at all.
             call setFilters( bounds_clump, glc_behavior )

with hr_col calculated in the Summary step above.

In the reworked sequence the Summary calls is part of EcosystemDynamicsPostDrainage

          call bgc_vegetation_inst%EcosystemDynamicsPostDrainage(bounds_clump, &
               filter(nc)%num_allc, filter(nc)%allc, &
               filter(nc)%num_bgc_soilc, filter(nc)%bgc_soilc, &
               filter(nc)%num_bgc_vegp, filter(nc)%bgc_vegp, &
               filter(nc)%num_actfirec, filter(nc)%actfirec,                 &
               filter(nc)%num_actfirep, filter(nc)%actfirep,                 &
               doalb, crop_inst, &
               soilstate_inst, soilbiogeochem_state_inst, &
               water_inst%waterstatebulk_inst, water_inst%waterdiagnosticbulk_inst, &
               water_inst%waterfluxbulk_inst, frictionvel_inst, canopystate_inst, &
               soilbiogeochem_carbonflux_inst, soilbiogeochem_carbonstate_inst, &
               c13_soilbiogeochem_carbonflux_inst, c13_soilbiogeochem_carbonstate_inst, &
               c14_soilbiogeochem_carbonflux_inst, c14_soilbiogeochem_carbonstate_inst, &
               soilbiogeochem_nitrogenflux_inst, soilbiogeochem_nitrogenstate_inst)
!...
          call clm_fates%WrapUpdateFatesRmean(nc,temperature_inst)


          call clm_fates%wrap_update_hifrq_hist(bounds_clump, &
               soilbiogeochem_carbonflux_inst, &
               soilbiogeochem_carbonstate_inst)


          if( is_beg_curr_day() ) then

             ! --------------------------------------------------------------------------
             ! This is the main call to FATES dynamics
             ! --------------------------------------------------------------------------

             call clm_fates%dynamics_driv( nc, bounds_clump,                        &
                  atm2lnd_inst, soilstate_inst, temperature_inst, active_layer_inst, &
                  water_inst%waterstatebulk_inst, water_inst%waterdiagnosticbulk_inst, &
                  water_inst%wateratm2lndbulk_inst, canopystate_inst, soilbiogeochem_carbonflux_inst, &
                  frictionvel_inst, soil_water_retention_curve)

             ! TODO(wjs, 2016-04-01) I think this setFilters call should be replaced by a
             ! call to reweight_wrapup, if it's needed at all.
             call setFilters( bounds_clump, glc_behavior )

[ekluzek]

So te sequence is a little different. Specifically I wonder if the call to

WrapUpdateFatesRmean

Needs to be moved to before EcosystemDynamicsPostDrainage, since that's how it's done in the baseline version?

[rgknox]

fates tests are now passing restarts, see: /glade/scratch/rgknox/tests_0810-143226ch

aux_clm tests running on cheyenne and izumi:
/glade/scratch/rgknox/tests_0810-210341ch
/scratch/cluster/rgknox/tests_0810-211127iz

[billsacks]

@ekluzek and others - we had talked in the standup meeting about making a comment here about the new logic around the allocations, but I see that @rgknox has changed this from what we were looking at, and now the logic is more in line with how we thought it should be – using an alloc_full_veg logical that is based on use_fates_bgc. (@rgknox the potential issue we discussed was where one processor had a filter size of 0 simply due to not having any vegetated points due to the decomposition, and so wanting to use a globally-consistent logical flag to control the allocation; we were especially concerned about the handling of history variables with the previous code. However, the changes you have made address this.)

[ekluzek]

I had a bit of an offline discussion with @rgknox about this that I think it's worth letting the group think about...

Ahh, OK, and the new name is alloc_full_veg. It's set in CNVegFacade and then passed down. I think this makes sense.
We tend to use global use_* settings for a ton of uses, and it might be better to make what they do more localized like this.
So this makes sense to me. And it sounds like @billsacks is good also from his comment he added to the PR.

My first reaction was that we should consistently use the global "use_fates*" type variables to control this type of thing, since you always knows what it means since it's a code global. But, in thinking more about this sometimes that covers up the specific reason one of those globals are used. So this local use is maybe a better paradigm, because you also know why it's being done. Here the real reason is "How should these arrays be allocated?". Right now it's one way for BGC and one for FATES, but there could be an expansion of that in the future. The name of the variable makes it obvious that that's what this is actually about, and it's local to just the subroutines that need it. So perhaps that is a better pattern to use for our code? Note for example this issue #942 where "use_cn" is used all over the code, but the reasons for its use vary, so its hard to know what the actual reason is. And perhaps using something that would be more localized would be better.

A programming adage that I like to follow especially for data arrays (and we have endorsed in CTSM) is "make everything as local as it can be". We've ignored that for control variables and parameters -- but maybe we shouldn't if possible?

Should this be something we should have a short group discussion on?

[ekluzek]

We did see the #2109 issue in the standard test list on Cheyenne for the FSURDATMODIFYCTSM test.

[samsrabin]

It actually looks like it's a different issue, which I've opened as #2111.