First-day host model SMP cause 100% seedling mortality when using the Tree Recruitment Scheme with seedling dynamics

[r-ward]

Describe the issue

Host model soil matric potential (SMP) spikes to ~negative 115 trillion on model day 1, which is remembered by the exponential moving averages that track accumulated moisture deficit new_seedling_mdd and new_seedling_layer_smp. The result is >100% moisture stress mortality rate, crashing fates when using the tree recruitment scheme with seedling dynamics.

Moisture deficit is accumulated in the UpdateFatesRMeansTStep subroutine. I'm pasting log output that shows the SPM spike (see log for c_index=2, by c_index=9 or 10 SMP returns to what seems to be a more reasonable value). Moisture deficit (new_seedling_mdd) is accumulated when deficit is beyond a critical threshold seedling_psi_crit and then added to the running mean Current EMA c_mean:

https://github.com/NGEET/fates/blob/3b58472a513e39fcb8d0dafc7bc0ae31b91f6e95/main/FatesInterfaceMod.F90#L2300C1-L2344C1

The moisture stress mortality rate is calculated using the running mean of seedling_mdd here:
https://github.com/NGEET/fates/blob/3b58472a513e39fcb8d0dafc7bc0ae31b91f6e95/biogeochem/EDPhysiologyMod.F90#L2316C1-L2344C32

On my branch I've avoided the problem by skipping the running mean updates on model day 1:

https://github.com/r-ward/fates/blob/f2c2233661320102a9eb46aed331e33a2c12936f/main/FatesInterfaceMod.F90#L2338C1-L2343C23

...but I'm not sure if that's the best solution, or if the SMP spike is concerning in itself!

Relevant log output

------------------------------------------------------------
 Successfully initialized the land model
 begin initial run at: 
    nstep=            1  year=         1985  month=            1  day=            1  seconds=         1800

************************************************************

 === MDD DEBUG Day:   1.0208333333333333     
 PFT:           1 Soil layer:           2
 Raw SMP from host:   0.0000000000000000     
 new_seedling_layer_smp:   0.0000000000000000     
 seedling_psi_crit:  -251995.70000000001     
 Deficit check val for check   251995.70000000001     
 new_seedling_mdd (after zero check):   0.0000000000000000     
 Current EMA c_mean:   0.0000000000000000     
 Current c_index:           1
 n_mem:        6048
 sdlng_mdd_timescale:   126.00000000000000     
 ======================== 
 --WARNING-- skipping CN balance check for first timesteps after startup or data assimilation
 === MDD DEBUG Day:   1.0208333333333333     
 PFT:           1 Soil layer:           2
 Raw SMP from host:  -115997489906016.48     
 new_seedling_layer_smp:  -115997489906016.48     
 seedling_psi_crit:  -251995.70000000001     
 Deficit check val for check  -115997489654020.78     
 new_seedling_mdd (after zero check):   14615683696406618.     
 Current EMA c_mean:   0.0000000000000000     
 Current c_index:           2
 n_mem:        6048
 sdlng_mdd_timescale:   126.00000000000000     
 ========================
 --WARNING-- skipping CN balance check for first timesteps after startup or data assimilation
 === MDD DEBUG Day:   1.0208333333333333     
 PFT:           1 Soil layer:           2
 Raw SMP from host:  -30257502454.039650     
 new_seedling_layer_smp:  -30257502454.039650     
 seedling_psi_crit:  -251995.70000000001     
 Deficit check val for check  -30257250458.339649     
 new_seedling_mdd (after zero check):   3812413557750.7959     
 Current EMA c_mean:   4871894565468872.0     
 Current c_index:           3
 n_mem:        6048
 sdlng_mdd_timescale:   126.00000000000000     
 ========================
 --WARNING-- skipping CN balance check for first timesteps after startup or data assimilation
 === MDD DEBUG Day:   1.0208333333333333     
 PFT:           1 Soil layer:           2
 Raw SMP from host:  -332845219.35745889     
 new_seedling_layer_smp:  -332845219.35745889     
 seedling_psi_crit:  -251995.70000000001     
 Deficit check val for check  -332593223.65745890     
 new_seedling_mdd (after zero check):   41906746180.839821     
 Current EMA c_mean:   3654874027491091.5     
 Current c_index:           4
 n_mem:        6048
 sdlng_mdd_timescale:   126.00000000000000     
 ========================
 === MDD DEBUG Day:   1.0208333333333333     
 PFT:           1 Soil layer:           2
 Raw SMP from host:  -38398293.236233741     
 new_seedling_layer_smp:  -38398293.236233741     
 seedling_psi_crit:  -251995.70000000001     
 Deficit check val for check  -38146297.536233738     
 new_seedling_mdd (after zero check):   4806433489.5654507     
 Current EMA c_mean:   2923907603342109.5     
 Current c_index:           5
 n_mem:        6048
 sdlng_mdd_timescale:   126.00000000000000     
 ========================
 === MDD DEBUG Day:   1.0208333333333333     
 PFT:           1 Soil layer:           2
 Raw SMP from host:  -8543591.0642576218     
 new_seedling_layer_smp:  -8543591.0642576218     
 seedling_psi_crit:  -251995.70000000001     
 Deficit check val for check  -8291595.3642576216     
 new_seedling_mdd (after zero check):   1044741015.8964603     
 Current EMA c_mean:   2436590470524006.5     
 Current c_index:           6
 n_mem:        6048
 sdlng_mdd_timescale:   126.00000000000000     
 ========================
 === MDD DEBUG Day:   1.0208333333333333     
 PFT:           1 Soil layer:           2
 Raw SMP from host:  -2752552.4450957361     
 new_seedling_layer_smp:  -2752552.4450957361     
 seedling_psi_crit:  -251995.70000000001     
 Deficit check val for check  -2500556.7450957359     
 new_seedling_mdd (after zero check):   315070149.88206273     
 Current EMA c_mean:   2088506266840722.2     
 Current c_index:           7
 n_mem:        6048
 sdlng_mdd_timescale:   126.00000000000000     
 ========================
 === MDD DEBUG Day:   1.0208333333333333     
 PFT:           1 Soil layer:           2
 Raw SMP from host:  -1023203.0811228864     
 new_seedling_layer_smp:  -1023203.0811228864     
 seedling_psi_crit:  -251995.70000000001     
 Deficit check val for check  -771207.38112288644     
 new_seedling_mdd (after zero check):   97172130.021483690     
 Current EMA c_mean:   1827443022869400.8     
 Current c_index:           8
 n_mem:        6048
 sdlng_mdd_timescale:   126.00000000000000     
 ========================
 === MDD DEBUG Day:   1.0208333333333333     
 PFT:           1 Soil layer:           2
 Raw SMP from host:  -440411.51312879863     
 new_seedling_layer_smp:  -440411.51312879863     
 seedling_psi_crit:  -251995.70000000001     
 Deficit check val for check  -188415.81312879862     
 new_seedling_mdd (after zero check):   23740392.454228625     
 Current EMA c_mean:   1624393808903037.2     
 Current c_index:           9
 n_mem:        6048
 sdlng_mdd_timescale:   126.00000000000000     
 ========================
 === MDD DEBUG Day:   1.0208333333333333     
 PFT:           1 Soil layer:           2
 Raw SMP from host:  -213703.35025217672     
 new_seedling_layer_smp:  -213703.35025217672     
 seedling_psi_crit:  -251995.70000000001     
 Deficit check val for check   38292.349747823289     
 new_seedling_mdd (after zero check):   0.0000000000000000     
 Current EMA c_mean:   1461954430386772.8     
 Current c_index:          10
 n_mem:        6048
 sdlng_mdd_timescale:   126.00000000000000     
 ========================
 === MDD DEBUG Day:   1.0208333333333333     
 PFT:           1 Soil layer:           2
 Raw SMP from host:  -135027.13194443836     
 new_seedling_layer_smp:  -135027.13194443836     
 seedling_psi_crit:  -251995.70000000001     
 Deficit check val for check   116968.56805556166     
 new_seedling_mdd (after zero check):   0.0000000000000000     
 Current EMA c_mean:   1329049482169793.2     
 Current c_index:          11
 n_mem:        6048
 sdlng_mdd_timescale:   126.00000000000000

FATES tag

sci.1.87.5_api.41.0.0-12-gf2c22336

Host land model tag

ctsm5.3.077

Machine

derecho

Other supported machine name

No response

Additional context

I have fixes for this and couple of additional bugs (which I'll make separate issues), on my branch debug_trs

[rosiealice]

Negative 115 trillion! That's maybe a record stupid-large number... Good job tracking it down. I feel like it might be good to fix this in the HLM(s), fundementally... Are you using ELM?

[r-ward]

Hahah...and then that number is essentially squared to calculate the h2o mort rate!

I'm using CTSM (5.3.077)

[rgknox]

Another nice find @r-ward . My gut tells me that we should be able to calculate a meaningful SMP on that first value. If we can specifically omit updating that variable, and not the others in the group of running means, that could work too.

[r-ward]

Here's the issue on the CTSM GitHub repo: ESCOMP/CTSM#3584 (comment)

[rgknox]

I think there are couple of issues here. Since this is the first timestep, I believe part of the issue is about the input boundary condition just simply not being ready (ie filled), bc_in(s)%smp_sl. But there is also a secondary vulnerability that we need to account for. If all or some of the soil is frozen, then there will also be wonky water potentials.
In our btran calculations we run checks to see if the layer has liquid water. See:

https://github.com/NGEET/fates/blob/sci.1.88.1_api.42.0.0/biogeophys/EDBtranMod.F90#L42

I think the solution here is to use that same check, and if there is no liquid water, do something in the trivial case.

[rgknox]

Note;
bc_in(s)%smp_sl(:) is set in CanopyFluxes (called on line 755 of clm_driver), which is indeed called earlier in the call sequence than WrapUpdateFatesRmean(nc,temperature_inst) (called on line 1170 of clm_driver).

So the issue is not about fates updating bc_in%smp_sl early enough in the call sequence, the issue either has to do with frozen soils or with the host side update of:
waterstatebulk_inst%h2osoi_vol_col(c,j) (and elm equivalent).

If it is the latter, then FATES btran calculations on the first step are using uninitialized water contents, but we probably just buffer our calculations using rails, ie keeping things between allowable upper and lower bounds... so it has minimal effect since its just one timestep.