Skip to content

Commit df47cfa

Browse files
rgknoxrgknox
committed
Added some functions for numerical checks.
1 parent 800e502 commit df47cfa

File tree

3 files changed

+267
-29
lines changed

3 files changed

+267
-29
lines changed

‎biogeochem/EDCanopyStructureMod.F90‎

Lines changed: 14 additions & 25 deletions
Original file line numberDiff line numberDiff line change
@@ -87,7 +87,9 @@ subroutine canopy_structure( currentSite , bc_in )
8787
! !USES:
8888

8989
use EDParamsMod, only : ED_val_comp_excln
90-
use EDtypesMod , only : ncwd, min_patch_area
90+
use EDtypesMod , only : ncwd
91+
use EDTypesMod , only : min_patch_area
92+
use EDTypesMod , only : val_check_ed_vars
9193
use FatesInterfaceMod, only : bc_in_type
9294
!
9395
! !ARGUMENTS
@@ -104,7 +106,9 @@ subroutine canopy_structure( currentSite , bc_in )
104106
integer :: patch_area_counter ! count iterations used to solve canopy areas
105107
logical :: area_not_balanced ! logical controlling if the patch layer areas
106108
! have successfully been redistributed
109+
integer :: return_code ! math checks on variables will return>0 if problems exist
107110
integer, parameter :: max_patch_iterations = 100
111+
108112

109113
!----------------------------------------------------------------------
110114

@@ -115,47 +119,32 @@ subroutine canopy_structure( currentSite , bc_in )
115119
currentSite%promotion_rate(:) = 0._r8
116120
currentSite%demotion_carbonflux = 0._r8
117121
currentSite%promotion_carbonflux = 0._r8
122+
118123
!
119124
! Section 1: Check total canopy area.
120125
!
121126
do while (associated(currentPatch)) ! Patch loop
122127

123128

124-
! Perform a numerical check on input data structures
125-
currentCohort => currentPatch%tallest
126-
do while (associated(currentCohort))
127-
if( currentCohort%dbh .ne. currentCohort%dbh ) then
128-
write(fates_log(),*) 'NAN COHORT DBH:'
129-
write(fates_log(),*) 'lat:',currentpatch%siteptr%lat
130-
write(fates_log(),*) 'lon:',currentpatch%siteptr%lon
131-
call endrun(msg=errMsg(sourcefile, __LINE__))
132-
end if
129+
! Perform numerical checks on some cohort and patch structures
130+
! ------------------------------------------------------------------------------
133131

134-
if( currentCohort%n .ne. currentCohort%n ) then
135-
write(fates_log(),*) 'NAN COHORT N'
136-
write(fates_log(),*) 'lat:',currentpatch%siteptr%lat
137-
write(fates_log(),*) 'lon:',currentpatch%siteptr%lon
138-
call endrun(msg=errMsg(sourcefile, __LINE__))
139-
end if
132+
call val_check_ed_vars(currentPatch,'co_n:co_dbh:pa_area',return_code)
133+
! No need to make error message, already generated in math_check_ed_vars
134+
if(return_code>0) call endrun(msg=errMsg(sourcefile, __LINE__))
140135

136+
! canopy layer has a special bounds check
137+
currentCohort => currentPatch%tallest
138+
do while (associated(currentCohort))
141139
if( currentCohort%canopy_layer < 1 .or. currentCohort%canopy_layer > nclmax+1 ) then
142140
write(fates_log(),*) 'lat:',currentpatch%siteptr%lat
143141
write(fates_log(),*) 'lon:',currentpatch%siteptr%lon
144142
write(fates_log(),*) 'BOGUS CANOPY LAYER: ',currentCohort%canopy_layer
145143
call endrun(msg=errMsg(sourcefile, __LINE__))
146144
end if
147-
148145
currentCohort => currentCohort%shorter
149146
enddo
150147

151-
if( currentPatch%area .ne. currentPatch%area )then
152-
write(fates_log(),*) 'NAN PATCH AREA'
153-
write(fates_log(),*) 'lat:',currentpatch%siteptr%lat
154-
write(fates_log(),*) 'lon:',currentpatch%siteptr%lon
155-
call endrun(msg=errMsg(sourcefile, __LINE__))
156-
end if
157-
158-
159148
if (currentPatch%area .gt. min_patch_area) then ! avoid numerical weirdness that shouldn't be happening anyway
160149

161150
! Does any layer have excess area in it? Keep going until it does not...

‎main/EDTypesMod.F90‎

Lines changed: 214 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -1,7 +1,8 @@
11
module EDTypesMod
22

3-
use FatesConstantsMod , only : r8 => fates_r8
4-
use shr_infnan_mod, only : nan => shr_infnan_nan, assignment(=)
3+
use FatesConstantsMod, only : r8 => fates_r8
4+
use FatesGlobals, only : fates_log
5+
use shr_infnan_mod, only : nan => shr_infnan_nan, assignment(=)
56

67
use FatesHydraulicsMemMod, only : ed_cohort_hydr_type
78
use FatesHydraulicsMemMod, only : ed_patch_hydr_type
@@ -602,5 +603,215 @@ function get_size_class_index(dbh) result(cohort_size_class)
602603
cohort_size_class = count(dbh-sclass_ed.ge.0.0_r8)
603604

604605
end function get_size_class_index
605-
606+
607+
! =====================================================================================
608+
609+
subroutine val_check_ed_vars(currentPatch,var_aliases,return_code)
610+
611+
! ----------------------------------------------------------------------------------
612+
! Perform numerical checks on variables of interest.
613+
! The input string is of the form: 'VAR1_NAME:VAR2_NAME:VAR3_NAME'
614+
! ----------------------------------------------------------------------------------
615+
616+
617+
use FatesUtilsMod,only : check_hlm_list
618+
use FatesUtilsMod,only : check_var_real
619+
620+
! Arguments
621+
type(ed_patch_type),intent(in), target :: currentPatch
622+
character(len=*),intent(in) :: var_aliases
623+
integer,intent(out) :: return_code ! return 0 for all fine
624+
! return 1 if a nan detected
625+
! return 10+ if an overflow
626+
! return 100% if an underflow
627+
! Locals
628+
type(ed_cohort_type), pointer :: currentCohort
629+
630+
631+
! Check through a registry of variables to check
632+
633+
if ( check_hlm_list(trim(var_aliases),'co_n') ) then
634+
635+
currentCohort => currentPatch%shortest
636+
do while(associated(currentCohort))
637+
call check_var_real(currentCohort%n,'cohort%n',return_code)
638+
if(.not.(return_code.eq.0)) then
639+
call dump_site(currentPatch%siteptr)
640+
call dump_patch(currentPatch)
641+
call dump_cohort(currentCohort)
642+
return
643+
end if
644+
currentCohort => currentCohort%taller
645+
end do
646+
end if
647+
648+
if ( check_hlm_list(trim(var_aliases),'co_dbh') ) then
649+
650+
currentCohort => currentPatch%shortest
651+
do while(associated(currentCohort))
652+
call check_var_real(currentCohort%dbh,'cohort%dbh',return_code)
653+
if(.not.(return_code.eq.0)) then
654+
call dump_site(currentPatch%siteptr)
655+
call dump_patch(currentPatch)
656+
call dump_cohort(currentCohort)
657+
return
658+
end if
659+
currentCohort => currentCohort%taller
660+
end do
661+
end if
662+
663+
if ( check_hlm_list(trim(var_aliases),'pa_area') ) then
664+
665+
call check_var_real(currentPatch%area,'patch%area',return_code)
666+
if(.not.(return_code.eq.0)) then
667+
call dump_site(currentPatch%siteptr)
668+
call dump_patch(currentPatch)
669+
return
670+
end if
671+
end if
672+
673+
674+
675+
return
676+
end subroutine val_check_ed_vars
677+
678+
! =====================================================================================
679+
680+
subroutine dump_site(csite)
681+
682+
type(ed_site_type),intent(in),target :: csite
683+
684+
685+
! EDTypes is
686+
687+
write(fates_log(),*) '----------------------------------------'
688+
write(fates_log(),*) ' Site Coordinates '
689+
write(fates_log(),*) '----------------------------------------'
690+
write(fates_log(),*) 'latitude = ', csite%lat
691+
write(fates_log(),*) 'longitude = ', csite%lon
692+
write(fates_log(),*) '----------------------------------------'
693+
return
694+
695+
end subroutine dump_site
696+
697+
! =====================================================================================
698+
699+
700+
subroutine dump_patch(cpatch)
701+
702+
type(ed_patch_type),intent(in),target :: cpatch
703+
704+
write(fates_log(),*) '----------------------------------------'
705+
write(fates_log(),*) ' Dumping Patch Information '
706+
write(fates_log(),*) ' (omitting arrays) '
707+
write(fates_log(),*) '----------------------------------------'
708+
write(fates_log(),*) 'pa%patchno = ',cpatch%patchno
709+
write(fates_log(),*) 'pa%age = ',cpatch%age
710+
write(fates_log(),*) 'pa%age_class = ',cpatch%age_class
711+
write(fates_log(),*) 'pa%area = ',cpatch%area
712+
write(fates_log(),*) 'pa%countcohorts = ',cpatch%countcohorts
713+
write(fates_log(),*) 'pa%ncl_p = ',cpatch%ncl_p
714+
write(fates_log(),*) 'pa%total_canopy_area = ',cpatch%total_canopy_area
715+
write(fates_log(),*) 'pa%total_tree_area = ',cpatch%total_tree_area
716+
write(fates_log(),*) 'pa%canopy_area = ',cpatch%canopy_area
717+
write(fates_log(),*) 'pa%bare_frac_area = ',cpatch%bare_frac_area
718+
write(fates_log(),*) 'pa%lai = ',cpatch%lai
719+
write(fates_log(),*) 'pa%zstar = ',cpatch%zstar
720+
write(fates_log(),*) 'pa%disturbance_rate = ',cpatch%disturbance_rate
721+
write(fates_log(),*) '----------------------------------------'
722+
return
723+
724+
end subroutine dump_patch
725+
726+
! =====================================================================================
727+
728+
subroutine dump_cohort(ccohort)
729+
730+
731+
type(ed_cohort_type),intent(in),target :: ccohort
732+
733+
write(fates_log(),*) '----------------------------------------'
734+
write(fates_log(),*) ' Dumping Cohort Information '
735+
write(fates_log(),*) '----------------------------------------'
736+
write(fates_log(),*) 'co%pft = ', ccohort%pft
737+
write(fates_log(),*) 'co%n = ', ccohort%n
738+
write(fates_log(),*) 'co%dbh = ', ccohort%dbh
739+
write(fates_log(),*) 'co%hite = ', ccohort%hite
740+
write(fates_log(),*) 'co%b = ', ccohort%b
741+
write(fates_log(),*) 'co%balive = ', ccohort%balive
742+
write(fates_log(),*) 'co%bdead = ', ccohort%bdead
743+
write(fates_log(),*) 'co%bstore = ', ccohort%bstore
744+
write(fates_log(),*) 'co%laimemory = ', ccohort%laimemory
745+
write(fates_log(),*) 'co%bsw = ', ccohort%bsw
746+
write(fates_log(),*) 'co%bl = ', ccohort%bl
747+
write(fates_log(),*) 'co%br = ', ccohort%br
748+
write(fates_log(),*) 'co%lai = ', ccohort%lai
749+
write(fates_log(),*) 'co%sai = ', ccohort%sai
750+
write(fates_log(),*) 'co%gscan = ', ccohort%gscan
751+
write(fates_log(),*) 'co%leaf_cost = ', ccohort%leaf_cost
752+
write(fates_log(),*) 'co%canopy_layer = ', ccohort%canopy_layer
753+
write(fates_log(),*) 'co%canopy_layer_yesterday = ', ccohort%canopy_layer_yesterday
754+
write(fates_log(),*) 'co%nv = ', ccohort%nv
755+
write(fates_log(),*) 'co%status_coh = ', ccohort%status_coh
756+
write(fates_log(),*) 'co%canopy_trim = ', ccohort%canopy_trim
757+
write(fates_log(),*) 'co%status_coh = ', ccohort%status_coh
758+
write(fates_log(),*) 'co%excl_weight = ', ccohort%excl_weight
759+
write(fates_log(),*) 'co%prom_weight = ', ccohort%prom_weight
760+
write(fates_log(),*) 'co%size_class = ', ccohort%size_class
761+
write(fates_log(),*) 'co%size_by_pft_class = ', ccohort%size_by_pft_class
762+
write(fates_log(),*) 'co%gpp_acc_hold = ', ccohort%gpp_acc_hold
763+
write(fates_log(),*) 'co%gpp_acc = ', ccohort%gpp_acc
764+
write(fates_log(),*) 'co%gpp_tstep = ', ccohort%gpp_tstep
765+
write(fates_log(),*) 'co%npp_acc_hold = ', ccohort%npp_acc_hold
766+
write(fates_log(),*) 'co%npp_tstep = ', ccohort%npp_tstep
767+
write(fates_log(),*) 'co%npp_acc = ', ccohort%npp_acc
768+
write(fates_log(),*) 'co%resp_tstep = ', ccohort%resp_tstep
769+
write(fates_log(),*) 'co%resp_acc = ', ccohort%resp_acc
770+
write(fates_log(),*) 'co%resp_acc_hold = ', ccohort%resp_acc_hold
771+
write(fates_log(),*) 'co%npp_leaf = ', ccohort%npp_leaf
772+
write(fates_log(),*) 'co%npp_froot = ', ccohort%npp_froot
773+
write(fates_log(),*) 'co%npp_bsw = ', ccohort%npp_bsw
774+
write(fates_log(),*) 'co%npp_bdead = ', ccohort%npp_bdead
775+
write(fates_log(),*) 'co%npp_bseed = ', ccohort%npp_bseed
776+
write(fates_log(),*) 'co%npp_store = ', ccohort%npp_store
777+
write(fates_log(),*) 'co%rdark = ', ccohort%rdark
778+
write(fates_log(),*) 'co%resp_m = ', ccohort%resp_m
779+
write(fates_log(),*) 'co%resp_g = ', ccohort%resp_g
780+
write(fates_log(),*) 'co%livestem_mr = ', ccohort%livestem_mr
781+
write(fates_log(),*) 'co%livecroot_mr = ', ccohort%livecroot_mr
782+
write(fates_log(),*) 'co%froot_mr = ', ccohort%froot_mr
783+
write(fates_log(),*) 'co%md = ', ccohort%md
784+
write(fates_log(),*) 'co%leaf_md = ', ccohort%leaf_md
785+
write(fates_log(),*) 'co%root_md = ', ccohort%root_md
786+
write(fates_log(),*) 'co%carbon_balance = ', ccohort%carbon_balance
787+
write(fates_log(),*) 'co%dmort = ', ccohort%dmort
788+
write(fates_log(),*) 'co%seed_prod = ', ccohort%seed_prod
789+
write(fates_log(),*) 'co%treelai = ', ccohort%treelai
790+
write(fates_log(),*) 'co%treesai = ', ccohort%treesai
791+
write(fates_log(),*) 'co%leaf_litter = ', ccohort%leaf_litter
792+
write(fates_log(),*) 'co%c_area = ', ccohort%c_area
793+
write(fates_log(),*) 'co%woody_turnover = ', ccohort%woody_turnover
794+
write(fates_log(),*) 'co%cmort = ', ccohort%cmort
795+
write(fates_log(),*) 'co%bmort = ', ccohort%bmort
796+
write(fates_log(),*) 'co%imort = ', ccohort%imort
797+
write(fates_log(),*) 'co%fmort = ', ccohort%fmort
798+
write(fates_log(),*) 'co%hmort = ', ccohort%hmort
799+
write(fates_log(),*) 'co%isnew = ', ccohort%isnew
800+
write(fates_log(),*) 'co%dndt = ', ccohort%dndt
801+
write(fates_log(),*) 'co%dhdt = ', ccohort%dhdt
802+
write(fates_log(),*) 'co%ddbhdt = ', ccohort%ddbhdt
803+
write(fates_log(),*) 'co%dbalivedt = ', ccohort%dbalivedt
804+
write(fates_log(),*) 'co%dbdeaddt = ', ccohort%dbdeaddt
805+
write(fates_log(),*) 'co%dbstoredt = ', ccohort%dbstoredt
806+
write(fates_log(),*) 'co%storage_flux = ', ccohort%storage_flux
807+
write(fates_log(),*) 'co%cfa = ', ccohort%cfa
808+
write(fates_log(),*) 'co%fire_mort = ', ccohort%fire_mort
809+
write(fates_log(),*) 'co%crownfire_mort = ', ccohort%crownfire_mort
810+
write(fates_log(),*) 'co%cambial_mort = ', ccohort%cambial_mort
811+
write(fates_log(),*) 'co%size_class = ', ccohort%size_class
812+
write(fates_log(),*) 'co%size_by_pft_class = ', ccohort%size_by_pft_class
813+
write(fates_log(),*) '----------------------------------------'
814+
return
815+
end subroutine dump_cohort
816+
606817
end module EDTypesMod

‎main/FatesUtilsMod.F90‎

Lines changed: 39 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -3,6 +3,9 @@ module FatesUtilsMod
33
! This module contains helper functions and subroutines which are general in nature.
44
! Think string parsing, timing, maybe numerics, etc.
55

6+
use FatesConstantsMod, only : r8 => fates_r8
7+
use FatesGlobals, only : fates_log
8+
69
contains
710

811

@@ -30,5 +33,40 @@ function check_hlm_list(hlms,hlm_name) result(astatus)
3033

3134
end function check_hlm_list
3235

33-
36+
! =====================================================================================
37+
38+
subroutine check_var_real(r8_var, var_name, return_code)
39+
40+
real(r8),intent(in) :: r8_var
41+
character(len=*),intent(in) :: var_name
42+
integer,intent(out) :: return_code
43+
44+
real(r8), parameter :: r8_type = 1.0
45+
real(r8), parameter :: overflow = huge(r8_type)
46+
real(r8), parameter :: underflow = tiny(r8_type)
47+
48+
return_code = 0
49+
50+
! NaN check
51+
if (r8_var /= r8_var) then
52+
write(fates_log(),*) 'NaN detected, ',trim(var_name),': ',r8_var
53+
return_code = 1
54+
end if
55+
56+
! Overflow check (within 100th of max precision)
57+
if (abs(r8_var) > 0.01*overflow) then
58+
write(fates_log(),*) 'Nigh overflow detected, ',trim(var_name),': ',r8_var
59+
return_code = return_code + 10
60+
end if
61+
62+
! Underflow check (within 100x of min precision)
63+
if (abs(r8_var) < 100.0_r8*underflow) then
64+
write(fates_log(),*) 'Nigh underflow detected, ',trim(var_name),': ',r8_var
65+
return_code = return_code + 100
66+
end if
67+
68+
69+
end subroutine check_var_real
70+
71+
3472
end module FatesUtilsMod

0 commit comments

Comments
 (0)