Bug fix: in copy cohort, npp_acc was not being copied! Now it is. Thi…

…s was leading to NaNs later on in runs and causing lots of problems. I added some numerical inquiry intrinsics to catch bad values, but I am considering putting them in a different place.

components/clm/src/ED/biogeochem/EDCohortDynamicsMod.F90

@@ -1038,13 +1038,15 @@ subroutine copy_cohort( currentCohort,copyc )
     n%gpp_acc_hold    = o%gpp_acc_hold
     n%gpp_acc         = o%gpp_acc
     n%gpp_tstep       = o%gpp_tstep
+
     n%npp_acc_hold    = o%npp_acc_hold
     n%npp_tstep       = o%npp_tstep
+    n%npp_acc         = o%npp_acc
 
     if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn Ia ',o%npp_acc
     if ( DEBUG ) write(fates_log(),*) 'EDcohortDyn Ib ',o%resp_acc
 
-    n%npp_acc_hold    = o%npp_acc_hold
+
     n%resp_tstep      = o%resp_tstep
     n%resp_acc        = o%resp_acc
     n%resp_acc_hold   = o%resp_acc_hold

components/clm/src/ED/biogeochem/EDPhysiologyMod.F90

@@ -873,6 +873,7 @@ subroutine Growth_Derivatives( currentSite, currentCohort)
           !what is the tax on the carbon available for growth? 
           currentCohort%carbon_balance = currentCohort%carbon_balance * (1.0_r8 - f_store)  
        else  !cbalance is negative. Take C out of store to pay for maintenance respn.
+
           currentCohort%storage_flux = currentCohort%carbon_balance 
 
           ! Note that npp_store only tracks the flux between NPP and storage.  Storage can 
@@ -887,14 +888,6 @@ subroutine Growth_Derivatives( currentSite, currentCohort)
        write(fates_log(),*) 'No target leaf area in GrowthDerivs? Bleaf(cohort) <= 0?'
        call endrun(msg=errMsg(sourcefile, __LINE__))
 
-!       currentCohort%storage_flux = currentCohort%carbon_balance
-!       ! Note that npp_store only tracks the flux between NPP and storage.  Storage can 
-!       ! also be drawn down to support some turnover demand.
-!       currentCohort%npp_store = min(0.0_r8,currentCohort%npp_acc_hold)
-!       currentCohort%carbon_balance = 0._r8
-!       write(iulog,*) 'ED: no leaf area in gd', currentCohort%indexnumber,currentCohort%n,currentCohort%bdead, &
-!             currentCohort%dbh,currentCohort%balive
-
     endif
 
     !Do we have enough carbon left over to make up the rest of the turnover demand? 

components/clm/src/ED/biogeophys/EDAccumulateFluxesMod.F90

@@ -9,14 +9,19 @@ module EDAccumulateFluxesMod
   ! Rosie Fisher. March 2014. 
   !
   ! !USES:
+  use abortutils, only        : endrun
+  use shr_log_mod , only      : errMsg => shr_log_errMsg
+  use FatesConstantsMod, only : fates_huge, fates_tiny
   implicit none
   private
   !
   public :: AccumulateFluxes_ED
 
   logical :: DEBUG = .false.  ! for debugging this module
-  !------------------------------------------------------------------------------
 
+  character(len=*), parameter, private :: sourcefile = &
+        __FILE__
+
 contains
 
   !------------------------------------------------------------------------------
@@ -33,6 +38,8 @@ subroutine AccumulateFluxes_ED(nsites, sites, bc_in, bc_out, dt_time)
     use EDTypesMod        , only : ed_patch_type, ed_cohort_type, &
                                    ed_site_type, AREA
     use FatesInterfaceMod , only : bc_in_type,bc_out_type
+    use, intrinsic :: IEEE_ARITHMETIC
+
     !
     ! !ARGUMENTS    
     integer,            intent(in)            :: nsites
@@ -67,12 +74,45 @@ subroutine AccumulateFluxes_ED(nsites, sites, bc_in, bc_out, dt_time)
                 ! _tstep fluxes are KgC/indiv/timestep _acc are KgC/indiv/day
 
                 if ( DEBUG ) then
-                   write(iulog,*) 'EDAccumFlux 64 ',ccohort%npp_acc, &
-                         ccohort%npp_tstep
+                   write(iulog,*) 'EDAccumFlux 64 ',ccohort%npp_tstep
                    write(iulog,*) 'EDAccumFlux 66 ',ccohort%gpp_tstep
                    write(iulog,*) 'EDAccumFlux 67 ',ccohort%resp_tstep
                 endif
+
+                ! Trap invalid values from photosynthesis and resp
+                ! -----------------------------------------------------------------------
+
+                if(ieee_is_nan(ccohort%gpp_tstep))then
+                   write(iulog,*)'GPP NaN Trap Triggered',s,ifp,ccohort%gpp_tstep
+                   call endrun(msg=errMsg(sourcefile, __LINE__))
+                end if
+
+                if(ieee_is_nan(ccohort%resp_tstep))then
+                   write(iulog,*)'RESP NaN Trap Triggered',s,ifp,ccohort%resp_tstep
+                   call endrun(msg=errMsg(sourcefile, __LINE__))
+                end if
 
+                if(ieee_is_nan(ccohort%npp_tstep))then
+                   write(iulog,*)'NPP NaN Trap Triggered',s,ifp,ccohort%npp_tstep
+                   call endrun(msg=errMsg(sourcefile, __LINE__))
+                end if
+
+                if(.not.ieee_is_finite(ccohort%gpp_tstep))then
+                   write(iulog,*)'GPP Infinite Trap Triggered',s,ifp,ccohort%gpp_tstep
+                   call endrun(msg=errMsg(sourcefile, __LINE__))
+                end if
+
+                if(.not.ieee_is_finite(ccohort%resp_tstep))then
+                   write(iulog,*)'RESP Infinite Trap Triggered',s,ifp,ccohort%resp_tstep
+                   call endrun(msg=errMsg(sourcefile, __LINE__))
+                end if
+
+                if(.not.ieee_is_finite(ccohort%npp_tstep))then
+                   write(iulog,*)'NPP Infinite Trap Triggered',s,ifp,ccohort%npp_tstep
+                   call endrun(msg=errMsg(sourcefile, __LINE__))
+                end if
+
+
                 ccohort%npp_acc  = ccohort%npp_acc  + ccohort%npp_tstep 
                 ccohort%gpp_acc  = ccohort%gpp_acc  + ccohort%gpp_tstep 
                 ccohort%resp_acc = ccohort%resp_acc + ccohort%resp_tstep

components/clm/src/ED/main/ChecksBalancesMod.F90

@@ -119,14 +119,6 @@ subroutine SummarizeNetFluxes( nsites, sites, bc_in, is_beg_day )
                      sum(currentPatch%seed_decay) + sum(currentPatch%leaf_litter_out) + &
                      sum(currentPatch%root_litter_out)) * &
                      ( currentPatch%area/AREA ) * 1.e3_r8 / ( 365.0_r8*SHR_CONST_CDAY )
-               print*,"TS TERMS:", sites(s)%fates_to_bgc_this_ts, &
-                     sum(currentPatch%CWD_AG_out), &
-                     sum(currentPatch%CWD_BG_out), &
-                     sum(currentPatch%seed_decay), & 
-                     sum(currentPatch%leaf_litter_out), &
-                     sum(currentPatch%root_litter_out), &
-                     ( currentPatch%area/AREA ) * 1.e3_r8 / ( 365.0_r8*SHR_CONST_CDAY )
-
                !
                sites(s)%tot_seed_rain_flux = sites(s)%tot_seed_rain_flux + &
                      sum(sites(s)%seed_rain_flux) * 1.e3_r8 / ( 365.0_r8*SHR_CONST_CDAY )
@@ -214,15 +206,6 @@ subroutine FATES_BGC_Carbon_Balancecheck(nsites, sites, bc_in, is_beg_day, dtime
                   sites(s)%fire_c_to_atm*SHR_CONST_CDAY)
             sites(s)%cbal_err_fates = sites(s)%cbal_err_fates / SHR_CONST_CDAY
 
-             print*,"ERR TERMS:",sites(s)%totfatesc, & 
-                  sites(s)%totfatesc_old, &
-                  sites(s)%npp_timeintegrated,  &
-                  sites(s)%tot_seed_rain_flux*SHR_CONST_CDAY,  &
-                  sites(s)%fates_to_bgc_this_ts*SHR_CONST_CDAY, &
-                  sites(s)%fire_c_to_atm*SHR_CONST_CDAY, &
-                  sites(s)%cbal_err_fates
-
-
             sites(s)%cbal_err_bgc = sites(s)%totbgcc - &
                   sites(s)%totbgcc_old - &
                   (sites(s)%fates_to_bgc_last_ts*SHR_CONST_CDAY - &

components/clm/src/ED/main/FatesConstantsMod.F90

@@ -55,4 +55,13 @@ module FatesConstantsMod
   real(fates_r8), parameter :: t_water_freeze_k_triple = 273.16_fates_r8      
 
 
+  ! For numerical inquiry
+  real(fates_r8), parameter :: fates_huge = huge(g_per_kg)
+
+  real(fates_r8), parameter :: fates_tiny = tiny(g_per_kg)
+
+
+
+
+
 end module FatesConstantsMod