From 02d2cf7927327bf703a88247ddc99cbc85bbdd9c Mon Sep 17 00:00:00 2001 From: Sam Rabin Date: Wed, 24 Jan 2024 15:44:32 -0700 Subject: [PATCH] Use collapse_to_dominant() in surfrd_hillslope(). For the Uniform and DominantPftLowland values of hillslope_pft_distribution_method --- .../include_user_mods | 1 + .../Hillslope_DominantPftLowland/user_nl_clm | 1 + src/main/surfrdMod.F90 | 70 +++++++------------ 3 files changed, 29 insertions(+), 43 deletions(-) create mode 100644 cime_config/testdefs/testmods_dirs/clm/Hillslope_DominantPftLowland/include_user_mods create mode 100644 cime_config/testdefs/testmods_dirs/clm/Hillslope_DominantPftLowland/user_nl_clm diff --git a/cime_config/testdefs/testmods_dirs/clm/Hillslope_DominantPftLowland/include_user_mods b/cime_config/testdefs/testmods_dirs/clm/Hillslope_DominantPftLowland/include_user_mods new file mode 100644 index 0000000000..fa2e50a80d --- /dev/null +++ b/cime_config/testdefs/testmods_dirs/clm/Hillslope_DominantPftLowland/include_user_mods @@ -0,0 +1 @@ +../Hillslope diff --git a/cime_config/testdefs/testmods_dirs/clm/Hillslope_DominantPftLowland/user_nl_clm b/cime_config/testdefs/testmods_dirs/clm/Hillslope_DominantPftLowland/user_nl_clm new file mode 100644 index 0000000000..4ab5c6cd56 --- /dev/null +++ b/cime_config/testdefs/testmods_dirs/clm/Hillslope_DominantPftLowland/user_nl_clm @@ -0,0 +1 @@ +hillslope_pft_distribution_method = 'DominantPftLowland' diff --git a/src/main/surfrdMod.F90 b/src/main/surfrdMod.F90 index 42534d418c..03a0082e97 100644 --- a/src/main/surfrdMod.F90 +++ b/src/main/surfrdMod.F90 @@ -896,11 +896,12 @@ subroutine surfrd_hillslope(begg, endg, ncid, ns) ! !USES: use clm_instur, only : ncolumns_hillslope, wt_nat_patch use clm_varctl, only : nhillslope,max_columns_hillslope - use clm_varpar, only : natpft_size, natpft_lb + use clm_varpar, only : natpft_size, natpft_lb, natpft_ub use ncdio_pio, only : ncd_inqdid, ncd_inqdlen use pftconMod , only : noveg use HillslopeHydrologyMod, only : pft_distribution_method, pft_standard, pft_from_file, pft_uniform_dominant_pft, pft_lowland_dominant_pft, pft_lowland_upland use array_utils, only: find_k_max_indices + use surfrdUtilsMod, only: collapse_to_dominant ! ! !ARGUMENTS: @@ -916,6 +917,8 @@ subroutine surfrd_hillslope(begg, endg, ncid, ns) logical :: readvar ! is variable on dataset integer,pointer :: arrayl(:) ! local array (needed because ncd_io expects a pointer) character(len=32) :: subname = 'surfrd_hillslope' ! subroutine name + logical, allocatable :: do_not_collapse(:) + integer :: n_dominant !----------------------------------------------------------------------- ! number of hillslopes per landunit @@ -960,52 +963,33 @@ subroutine surfrd_hillslope(begg, endg, ncid, ns) endif enddo - ! pft_uniform_dominant_pft uses the patch with the - ! largest weight for all hillslope columns in the gridcell - else if (pft_distribution_method == pft_uniform_dominant_pft) then - allocate(max_indices(1)) - do g = begg, endg - ! If hillslopes will be used in a gridcell, modify wt_nat_patch, - ! otherwise use original patch distribution - if(ncolumns_hillslope(g) > 0) then - - call find_k_max_indices(wt_nat_patch(g,:),natpft_lb,1,max_indices) - wt_nat_patch(g,:) = 0._r8 - wt_nat_patch(g,max_indices(1)) = 100._r8 + else if (pft_distribution_method == pft_uniform_dominant_pft & + .or. pft_distribution_method == pft_lowland_dominant_pft) then - endif - enddo - deallocate(max_indices) - - ! pft_lowland_dominant_pft uses the two patches with the - ! largest weights for the hillslope columns in the gridcell - else if (pft_distribution_method == pft_lowland_dominant_pft) then - allocate(max_indices(2)) + ! If hillslopes will be used in a gridcell, modify wt_nat_patch, + ! otherwise use original patch distribution + allocate(do_not_collapse(begg:endg)) + do_not_collapse(begg:endg) = .false. do g = begg, endg - ! If hillslopes will be used in a gridcell, modify wt_nat_patch, otherwise use original patch distribution - if(ncolumns_hillslope(g) > 0) then + if (ncolumns_hillslope(g) == 0) then + do_not_collapse(g) = .true. + end if + end do - ! Preserve the relative weights of the largest and - ! next largest weights using arbitrarily chosen values - ! (i.e. 1 should be larger than 2) that sum to 100. - ! This will minimize memory usage while still allowing - ! HillslopeDominantLowlandPft to pick out the two largest patch types. - - call find_k_max_indices(wt_nat_patch(g,:),natpft_lb,2,max_indices) - ! check that 2nd index weight is non-zero - if (wt_nat_patch(g,max_indices(2)) > 0._r8) then - wt_nat_patch(g,:) = 0._r8 - wt_nat_patch(g,max_indices(1)) = 75._r8 - wt_nat_patch(g,max_indices(2)) = 25._r8 - else - ! if only one pft exists, set its weight to 100 per cent - wt_nat_patch(g,:) = 0._r8 - wt_nat_patch(g,max_indices(1)) = 100._r8 - endif + if (pft_distribution_method == pft_uniform_dominant_pft) then + ! pft_uniform_dominant_pft uses the patch with the + ! largest weight for all hillslope columns in the gridcell + n_dominant = 1 + else if (pft_distribution_method == pft_lowland_dominant_pft) then + ! pft_lowland_dominant_pft uses the two patches with the + ! largest weights for the hillslope columns in the gridcell + n_dominant = 2 + else + call endrun( msg=' ERROR: unrecognized hillslope_pft_distribution_method'//errMsg(sourcefile, __LINE__)) + end if - endif - enddo - deallocate(max_indices) + call collapse_to_dominant(wt_nat_patch(begg:endg,:), natpft_lb, natpft_ub, begg, endg, n_dominant, do_not_collapse) + deallocate(do_not_collapse) else if (pft_distribution_method /= pft_standard) then call endrun( msg=' ERROR: unrecognized hillslope_pft_distribution_method'//errMsg(sourcefile, __LINE__))