diff --git a/bld/CLMBuildNamelist.pm b/bld/CLMBuildNamelist.pm index 42e0b35b92..1ca5bebdc7 100755 --- a/bld/CLMBuildNamelist.pm +++ b/bld/CLMBuildNamelist.pm @@ -2089,10 +2089,31 @@ sub setup_logic_soilstate { my ($opts, $nl_flags, $definition, $defaults, $nl) = @_; add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'organic_frac_squared' ); - add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'soil_layerstruct', - 'structure'=>$nl_flags->{'structure'}); add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'use_bedrock', 'use_fates'=>$nl_flags->{'use_fates'}, 'vichydro'=>$nl_flags->{'vichydro'} ); + + my $var1 = "soil_layerstruct_predefined"; + my $var2 = "soil_layerstruct_userdefined"; + my $var3 = "soil_layerstruct_userdefined_nlevsoi"; + my $soil_layerstruct_predefined = $nl->get_value($var1); + my $soil_layerstruct_userdefined = $nl->get_value($var2); + my $soil_layerstruct_userdefined_nlevsoi = $nl->get_value($var3); + + if (defined($soil_layerstruct_userdefined)) { + if (defined($soil_layerstruct_predefined)) { + $log->fatal_error("You have set both soil_layerstruct_userdefined and soil_layerstruct_predefined in your namelist; model cannot determine which to use"); + } + if (not defined($soil_layerstruct_userdefined_nlevsoi)) { + $log->fatal_error("You have set soil_layerstruct_userdefined and NOT set soil_layerstruct_userdefined_nlevsoi in your namelist; both MUST be set"); + } + } else { + if (defined($soil_layerstruct_userdefined_nlevsoi)) { + $log->fatal_error("You have set soil_layerstruct_userdefined_nlevsoi and NOT set soil_layerstruct_userdefined in your namelist; EITHER set both OR neither; in the latter case soil_layerstruct_predefined will be assigned a default value"); + } else { + add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'soil_layerstruct_predefined', + 'structure'=>$nl_flags->{'structure'}); + } + } } #------------------------------------------------------------------------------- diff --git a/bld/namelist_files/namelist_defaults_ctsm.xml b/bld/namelist_files/namelist_defaults_ctsm.xml index c9da63db74..04d4412629 100644 --- a/bld/namelist_files/namelist_defaults_ctsm.xml +++ b/bld/namelist_files/namelist_defaults_ctsm.xml @@ -125,9 +125,9 @@ attributes from the config_cache.xml file (with keys converted to upper-case). .true. .false. -5SL_3m -20SL_8.5m -10SL_3.5m +4SL_2m +20SL_8.5m +10SL_3.5m .false. .false. diff --git a/bld/namelist_files/namelist_definition_ctsm.xml b/bld/namelist_files/namelist_definition_ctsm.xml index 33947ebe4d..6e85a95fd4 100644 --- a/bld/namelist_files/namelist_definition_ctsm.xml +++ b/bld/namelist_files/namelist_definition_ctsm.xml @@ -114,13 +114,25 @@ If TRUE, square the organic fraction when it's used (as was done in CLM4.5) Otherwise use the fraction straight up (the default for CLM5.0) - + 10SL_3.5m = standard CLM4 and CLM4.5 version 23SL_3.5m = more vertical layers for permafrost simulations 49SL_10m = 49 layer soil column, 10m of soil, 5 bedrock layers 20SL_8.5m = 20 layer soil column, 8m of soil, 5 bedrock layers -5SL_3m = 5 layer soil column, 3m of soil, 0 bedrock layers +4SL_2m = 4 layer soil column, 2m of soil, 0 bedrock layers + + + +User-defined vector of dzsoi. The length of this vector determines nlevgrnd. When the user sets this vector, they have to set soil_layerstruct_userdefined_nlevsoi in the namelist, too; soil_layerstruct_userdefined_nlevsoi must be less than nlevgrnd in this version of the model, even though ideally soil_layerstruct_userdefined_nlevsoi could also equal nlevgrnd. +Default: rundef + + + +User-defined number of soil layers required to be set in the namelist when the user sets soil_layerstruct_userdefined in the namelist. +Default: iundef clm4.5: - clm5.0: + clm5.0: Satellite phenology: CN: Carbon Nitrogen model CNDV: CN with Dynamic Vegetation @@ -32,6 +32,7 @@ BGC (vert. resol. CN and methane) with prognostic crop, with modifications appropriate for CMIP6 WACCM DECK experiments: NWP configuration with satellite phenology: + NWP configuration with BGC and CROP: diff --git a/cime_config/config_compsets.xml b/cime_config/config_compsets.xml index e0acc4f4f3..f807de2d3f 100644 --- a/cime_config/config_compsets.xml +++ b/cime_config/config_compsets.xml @@ -248,6 +248,11 @@ 2000_DATM%GSWP3v1_CLM50%NWP-SP_SICE_SOCN_MOSART_SGLC_SWAV + + I2000Ctsm50NwpBgcCropGswpGs + 2000_DATM%GSWP3v1_CLM50%NWP-BGC-CROP_SICE_SOCN_MOSART_SGLC_SWAV + + diff --git a/cime_config/config_tests.xml b/cime_config/config_tests.xml index dbcd59eb6c..2d86e3b71d 100644 --- a/cime_config/config_tests.xml +++ b/cime_config/config_tests.xml @@ -87,4 +87,14 @@ SSP smoke CLM spinup test (only valid for CLM compsets with CLM45) $STOP_N + + CLM user-defined soil structure test + 1 + FALSE + FALSE + never + $STOP_OPTION + $STOP_N + + diff --git a/cime_config/testdefs/testlist_clm.xml b/cime_config/testdefs/testlist_clm.xml index e88924d3d5..f68c877ed2 100644 --- a/cime_config/testdefs/testlist_clm.xml +++ b/cime_config/testdefs/testlist_clm.xml @@ -787,13 +787,13 @@ - + - + @@ -1674,6 +1674,18 @@ + + + + + + + + + + + + diff --git a/cime_config/testdefs/testmods_dirs/clm/deepsoil_bedrock/user_nl_clm b/cime_config/testdefs/testmods_dirs/clm/deepsoil_bedrock/user_nl_clm index 0c98d81ff0..daaf704226 100644 --- a/cime_config/testdefs/testmods_dirs/clm/deepsoil_bedrock/user_nl_clm +++ b/cime_config/testdefs/testmods_dirs/clm/deepsoil_bedrock/user_nl_clm @@ -1,3 +1,3 @@ lower_boundary_condition = 2 -soil_layerstruct = '49SL_10m' +soil_layerstruct_predefined = '49SL_10m' use_bedrock = .true. diff --git a/cime_config/testdefs/testmods_dirs/clm/rm_indiv_lunits_and_collapse_to_dom/README b/cime_config/testdefs/testmods_dirs/clm/rm_indiv_lunits_and_collapse_to_dom/README deleted file mode 100644 index 15d6992b29..0000000000 --- a/cime_config/testdefs/testmods_dirs/clm/rm_indiv_lunits_and_collapse_to_dom/README +++ /dev/null @@ -1,10 +0,0 @@ -This test -1) Keeps landunits with areas exceeding the correspondng "toosmall_*" namelist - thresholds and removes the rest -2) Collapses unmanaged pfts to the 2 dominant ones -3) Collapses urban landunits to the dominant landunit -4) Collapses all landunits to the two dominant ones - -All this simplifies the gridcell representation so as to run faster. - -NB: This test is testing the 10x15 resolution only. diff --git a/cime_config/testdefs/testmods_dirs/clm/rm_indiv_lunits_and_collapse_to_dom/include_user_mods b/cime_config/testdefs/testmods_dirs/clm/rm_indiv_lunits_and_collapse_to_dom/include_user_mods deleted file mode 100644 index fe0e18cf88..0000000000 --- a/cime_config/testdefs/testmods_dirs/clm/rm_indiv_lunits_and_collapse_to_dom/include_user_mods +++ /dev/null @@ -1 +0,0 @@ -../default diff --git a/cime_config/testdefs/testmods_dirs/clm/rm_indiv_lunits_and_collapse_to_dom/user_nl_clm b/cime_config/testdefs/testmods_dirs/clm/rm_indiv_lunits_and_collapse_to_dom/user_nl_clm deleted file mode 100644 index 1ed05f191b..0000000000 --- a/cime_config/testdefs/testmods_dirs/clm/rm_indiv_lunits_and_collapse_to_dom/user_nl_clm +++ /dev/null @@ -1,9 +0,0 @@ -toosmall_soil = 5.0 -toosmall_crop = 5.0 -toosmall_glacier = 5.0 -toosmall_lake = 5.0 -toosmall_wetland = 5.0 -toosmall_urban = 5.0 -n_dom_pfts = 2 -n_dom_landunits = 2 -collapse_urban = .true. diff --git a/cime_config/testdefs/testmods_dirs/clm/vrtlay/user_nl_clm b/cime_config/testdefs/testmods_dirs/clm/vrtlay/user_nl_clm index cccaf2abf2..15a92be2d7 100644 --- a/cime_config/testdefs/testmods_dirs/clm/vrtlay/user_nl_clm +++ b/cime_config/testdefs/testmods_dirs/clm/vrtlay/user_nl_clm @@ -1 +1 @@ -soil_layerstruct = '23SL_3.5m' +soil_layerstruct_predefined = '23SL_3.5m' diff --git a/doc/ChangeLog b/doc/ChangeLog index b2a1e1d818..ea2eb10e9b 100644 --- a/doc/ChangeLog +++ b/doc/ChangeLog @@ -1,4 +1,145 @@ =============================================================== +Tag name: ctsm1.0.dev053 +Originator(s): slevis (Samuel Levis,SLevis Consulting LLC,303-665-1310) +Date: Thu Aug 1 16:56:09 MDT 2019 +One-line Summary: Soil layer definition clean-up and user-defined option + +Purpose of changes +------------------ + + Code clean-up clarifes that there are two types of soil layer + definition: the node-based and the thickness-based. + + User-defined option allows user to specify a soil layer profile in the + form of a dzsoi vector (values in meters) in the thickness-based + approach. + + Default nlevsoi for NWP configurations had to change from 5 to 4 for + consistency with the new error check described in known bugs below. + + Other code clean-up removes a couple of sections of repeating code. + +Bugs fixed or introduced +------------------------ + +Issues fixed (include CTSM Issue #): #279 #728 + +Known bugs found since the previous tag (include github issue ID): + #759 (this PR) bug causes model to abort when nlevsoi = nlevgrnd; + bug has been corrected with an error check + + +Significant changes to scientifically-supported configurations +-------------------------------------------------------------- + +Does this tag change answers significantly for any of the following physics configurations? +(Details of any changes will be given in the "Answer changes" section below.) + + [Put an [X] in the box for any configuration with significant answer changes.] + +[ ] clm5_0 + +[X] ctsm5_0-nwp + +[ ] clm4_5 + +Notes of particular relevance for users +--------------------------------------- + +Caveats for users (e.g., need to interpolate initial conditions): + if neither soil_layerstruct_predefined nor soil_layerstruct_userdefined + get specified in the namelist, then the model sets + soil_layerstruct_predefined to the old default setting for + soil_layerstruct (clm5: 20SL_8.5m, clm4.5: 10SL_3.5m) + +Changes to CTSM's user interface (e.g., new/renamed XML or namelist variables): + renamed soil_layerstruct to soil_layerstruct_predefined and added + soil_layerstruct_userdefined + +Changes made to namelist defaults (e.g., changed parameter values): + Default nlevsoi for NWP configurations had to change from 5 to 4 for + consistency with the new error check described in known bugs below + +Changes to the datasets (e.g., parameter, surface or initial files): none + +Substantial timing or memory changes: none + +Notes of particular relevance for developers +-------------------------------------------- +NOTE: Be sure to review the steps in README.CHECKLIST.master_tags as well as the coding style in the Developers Guide + +Caveats for developers (e.g., code that is duplicated that requires double maintenance): none + +Changes to tests or testing: + 1) New test... + ERP_P36x2_D_Ld5.f10_f10_musgs.I2000Ctsm50NwpBgcCropGswpGs.cheyenne_intel.clm-default + replaces existing test + ERS_D_Ld10.f10_f10_musgs.I2000Clm50BgcCropGs.cheyenne_intel.clm-rm_indiv_lunits_and_collapse_to_dom + to check the correction described in known bugs above; the new test + together with existing unit tests cover what the old test was testing + + 2) New test and new test type... + SOILSTRUCTUD_Ld5.f10_f10_musgs.I2000Clm50BgcCropGs.cheyenne_intel.clm-default + ensures that soil_layerstruct_userdefined gives bfb same answers as + soil_layerstruct_predefined = '4SL_2m' when set with the same dzsoi + values. The new test type was put together by: + - listing the test in (1) testlist_clm.xml (as all tests) and (2) config_tests.xml + - creating the file .../cime_config/SystemTests/soilstructud.py named after the test in lower case + +Code reviewed by: @billsacks + + +CTSM testing: + + [PASS means all tests PASS and OK means tests PASS other than expected fails.] + + build-namelist tests: + + cheyenne - + + tools-tests (test/tools): + + cheyenne - + + PTCLM testing (tools/shared/PTCLM/test): + + cheyenne - + + python testing (see instructions in python/README.md; document testing done): + + (any machine) - + + regular tests (aux_clm): + + cheyenne ---- OK + hobart ------ OK + +If the tag used for baseline comparisons was NOT the previous tag, note that here: + + +Answer changes +-------------- + +Changes answers relative to baseline: + + Summarize any changes to answers, i.e., + - what code configurations: nwp + - what platforms/compilers: all + - nature of change: larger than roundoff/same climate + This is due to the change of nlevsoi from 5 to 4 (more info above). + I confirmed that nwp does return bfb same answers when I revert this + change. + +Detailed list of changes +------------------------ + +List any externals directories updated (cime, rtm, mosart, cism, fates, etc.): none + +Pull Requests that document the changes (include PR ids): + https://github.com/ESCOMP/ctsm/pull/759 + +=============================================================== +=============================================================== Tag name: ctsm1.0.dev052 Originator(s): sacks (Bill Sacks) Date: Mon Jul 22 14:02:43 MDT 2019 diff --git a/doc/ChangeSum b/doc/ChangeSum index 701d7dda0d..ab379bb589 100644 --- a/doc/ChangeSum +++ b/doc/ChangeSum @@ -1,5 +1,6 @@ Tag Who Date Summary ============================================================================================================================ + ctsm1.0.dev053 slevis 08/01/2019 Soil layer definition clean-up and user-defined option ctsm1.0.dev052 sacks 07/22/2019 Fix rare soil color bug in mksurfdata_map ctsm1.0.dev051 sacks 07/19/2019 Update water tracers for remainder of first stage of hydrology ctsm1.0.dev050 slevis 07/15/2019 dz --> dz_lake bug-fix in LakeTemperatureMod.F90 line 960 diff --git a/src/biogeophys/SoilHydrologyInitTimeConstMod.F90 b/src/biogeophys/SoilHydrologyInitTimeConstMod.F90 index d10a842015..64ede52fc4 100644 --- a/src/biogeophys/SoilHydrologyInitTimeConstMod.F90 +++ b/src/biogeophys/SoilHydrologyInitTimeConstMod.F90 @@ -8,10 +8,12 @@ module SoilHydrologyInitTimeConstMod use shr_kind_mod , only : r8 => shr_kind_r8 use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type + use SoilStateType , only : soilstate_type use SoilHydrologyType , only : soilhydrology_type use WaterStateBulkType, only : waterstatebulk_type use LandunitType , only : lun use ColumnType , only : col + use SoilStateInitTimeConstMod, only: organic_max ! implicit none private @@ -31,27 +33,27 @@ module SoilHydrologyInitTimeConstMod contains !----------------------------------------------------------------------- - subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst) + subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst, soilstate_inst) ! ! !USES: use shr_const_mod , only : shr_const_pi use shr_spfn_mod , only : shr_spfn_erf use abortutils , only : endrun use spmdMod , only : masterproc - use clm_varctl , only : fsurdat, paramfile, iulog, use_vichydro, soil_layerstruct - use clm_varpar , only : nlevsoifl, toplev_equalspace - use clm_varpar , only : nlevsoi, nlevgrnd, nlevsno, nlevlak, nlevurb, nlayer, nlayert - use clm_varcon , only : zsoi, dzsoi, zisoi, spval, nlvic, dzvic, pc, grlnd + use clm_varctl , only : fsurdat, iulog, use_vichydro + use clm_varpar , only : toplev_equalspace + use clm_varpar , only : nlevsoi, nlevgrnd, nlayer, nlayert + use clm_varcon , only : dzsoi, spval, nlvic, dzvic, pc, grlnd use clm_varcon , only : aquifer_water_baseline - use landunit_varcon , only : istwet, istsoil, istdlak, istcrop, istice_mec + use landunit_varcon , only : istwet, istdlak, istice_mec use column_varcon , only : icol_shadewall, icol_road_perv, icol_road_imperv, icol_roof, icol_sunwall use fileutils , only : getfil - use organicFileMod , only : organicrd use ncdio_pio , only : file_desc_t, ncd_io, ncd_pio_openfile, ncd_pio_closefile ! ! !ARGUMENTS: type(bounds_type) , intent(in) :: bounds type(soilhydrology_type) , intent(inout) :: soilhydrology_inst + type(soilstate_type) , intent(in) :: soilstate_inst ! ! !LOCAL VARIABLES: integer :: p,c,j,l,g,lev,nlevs @@ -62,9 +64,6 @@ subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst) logical :: readvar type(file_desc_t) :: ncid character(len=256) :: locfn - real(r8) :: clay,sand ! temporaries - real(r8) :: om_frac ! organic matter fraction - real(r8) :: organic_max ! organic matter (kg/m3) where soil is assumed to act like peat real(r8) ,pointer :: b2d (:) ! read in - VIC b real(r8) ,pointer :: ds2d (:) ! read in - VIC Ds real(r8) ,pointer :: dsmax2d (:) ! read in - VIC Dsmax @@ -72,12 +71,6 @@ subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst) real(r8), pointer :: sandcol (:,:) ! column level sand fraction for calculating VIC parameters real(r8), pointer :: claycol (:,:) ! column level clay fraction for calculating VIC parameters real(r8), pointer :: om_fraccol (:,:) ! column level organic matter fraction for calculating VIC parameters - real(r8) ,pointer :: sand3d (:,:) ! read in - soil texture: percent sand - real(r8) ,pointer :: clay3d (:,:) ! read in - soil texture: percent clay - real(r8) ,pointer :: organic3d (:,:) ! read in - organic matter: kg/m3 - real(r8) ,pointer :: zisoifl (:) ! original soil interface depth - real(r8) ,pointer :: zsoifl (:) ! original soil midpoint - real(r8) ,pointer :: dzsoifl (:) ! original soil thickness !----------------------------------------------------------------------- ! Initialize VIC variables @@ -152,54 +145,6 @@ subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst) end do end do - allocate(sand3d(bounds%begg:bounds%endg,nlevsoifl)) - allocate(clay3d(bounds%begg:bounds%endg,nlevsoifl)) - allocate(organic3d(bounds%begg:bounds%endg,nlevsoifl)) - - call organicrd(organic3d) - - call getfil (fsurdat, locfn, 0) - call ncd_pio_openfile (ncid, locfn, 0) - call ncd_io(ncid=ncid, varname='PCT_SAND', flag='read', data=sand3d, dim1name=grlnd, readvar=readvar) - if (.not. readvar) then - call endrun(msg=' ERROR: PCT_SAND NOT on surfdata file'//errMsg(sourcefile, __LINE__)) - end if - call ncd_io(ncid=ncid, varname='PCT_CLAY', flag='read', data=clay3d, dim1name=grlnd, readvar=readvar) - if (.not. readvar) then - call endrun(msg=' ERROR: PCT_CLAY NOT on surfdata file'//errMsg(sourcefile, __LINE__)) - end if - call ncd_pio_closefile(ncid) - - ! Determine organic_max - call getfil (paramfile, locfn, 0) - call ncd_pio_openfile (ncid, trim(locfn), 0) - call ncd_io(ncid=ncid, varname='organic_max', flag='read', data=organic_max, readvar=readvar) - if ( .not. readvar ) then - call endrun(msg=' ERROR: organic_max not on param file'//errMsg(sourcefile, __LINE__)) - end if - call ncd_pio_closefile(ncid) - - ! get original soil depths to be used in interpolation of sand and clay - allocate(zsoifl(1:nlevsoifl), zisoifl(0:nlevsoifl), dzsoifl(1:nlevsoifl)) - do j = 1, nlevsoifl - zsoifl(j) = 0.025*(exp(0.5_r8*(j-0.5_r8))-1._r8) !node depths - enddo - - dzsoifl(1) = 0.5_r8*(zsoifl(1)+zsoifl(2)) !thickness b/n two interfaces - do j = 2,nlevsoifl-1 - dzsoifl(j)= 0.5_r8*(zsoifl(j+1)-zsoifl(j-1)) - enddo - dzsoifl(nlevsoifl) = zsoifl(nlevsoifl)-zsoifl(nlevsoifl-1) - - zisoifl(0) = 0._r8 - do j = 1, nlevsoifl-1 - zisoifl(j) = 0.5_r8*(zsoifl(j)+zsoifl(j+1)) !interface depths - enddo - zisoifl(nlevsoifl) = zsoifl(nlevsoifl) + 0.5_r8*dzsoifl(nlevsoifl) - - if ( masterproc )then - if ( soil_layerstruct /= '10SL_3.5m' ) write(iulog,*) 'Setting clay, sand, organic, in Soil Hydrology for VIC' - end if do c = bounds%begc, bounds%endc g = col%gridcell(c) l = col%landunit(c) @@ -211,41 +156,16 @@ subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst) ! do nothing else do lev = 1,nlevgrnd - if ( soil_layerstruct /= '10SL_3.5m' )then - if (lev .eq. 1) then - clay = clay3d(g,1) - sand = sand3d(g,1) - om_frac = organic3d(g,1)/organic_max - else if (lev <= nlevsoi) then - do j = 1,nlevsoifl-1 - if (zisoi(lev) >= zisoifl(j) .AND. zisoi(lev) < zisoifl(j+1)) then - clay = clay3d(g,j+1) - sand = sand3d(g,j+1) - om_frac = organic3d(g,j+1)/organic_max - endif - end do - else - clay = clay3d(g,nlevsoifl) - sand = sand3d(g,nlevsoifl) - om_frac = 0._r8 - endif + if ( lev <= nlevsoi )then + claycol(c,lev) = soilstate_inst%cellclay_col(c,lev) + sandcol(c,lev) = soilstate_inst%cellsand_col(c,lev) + om_fraccol(c,lev) = soilstate_inst%cellorg_col(c,lev) / organic_max else - ! duplicate clay and sand values from 10th soil layer - if (lev <= nlevsoi) then - clay = clay3d(g,lev) - sand = sand3d(g,lev) - om_frac = (organic3d(g,lev)/organic_max)**2._r8 - else - clay = clay3d(g,nlevsoi) - sand = sand3d(g,nlevsoi) - om_frac = 0._r8 - endif + claycol(c,lev) = soilstate_inst%cellclay_col(c,nlevsoi) + sandcol(c,lev) = soilstate_inst%cellsand_col(c,nlevsoi) + om_fraccol(c,lev) = 0.0_r8 end if - if (lun%urbpoi(l)) om_frac = 0._r8 - claycol(c,lev) = clay - sandcol(c,lev) = sand - om_fraccol(c,lev) = om_frac end do end if end if ! end of if not lake @@ -276,8 +196,6 @@ subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst) deallocate(b2d, ds2d, dsmax2d, ws2d) deallocate(sandcol, claycol, om_fraccol) - deallocate(sand3d, clay3d, organic3d) - deallocate(zisoifl, zsoifl, dzsoifl) end if ! end of if use_vichydro diff --git a/src/biogeophys/SoilStateInitTimeConstMod.F90 b/src/biogeophys/SoilStateInitTimeConstMod.F90 index 136c21ca1c..0a1e205c3a 100644 --- a/src/biogeophys/SoilStateInitTimeConstMod.F90 +++ b/src/biogeophys/SoilStateInitTimeConstMod.F90 @@ -5,6 +5,7 @@ module SoilStateInitTimeConstMod ! Set hydraulic and thermal properties ! ! !USES + use shr_kind_mod , only : r8 => shr_kind_r8 use SoilStateType , only : soilstate_type use LandunitType , only : lun use ColumnType , only : col @@ -19,6 +20,9 @@ module SoilStateInitTimeConstMod ! !PRIVATE MEMBER FUNCTIONS: private :: ReadNL ! + ! !PUBLIC DATA: + real(r8), public :: organic_max ! organic matter (kg/m3) where soil is assumed to act like peat + ! !PRIVATE DATA: ! Control variables (from namelist) logical, private :: organic_frac_squared ! If organic fraction should be squared (as in CLM4.5) @@ -87,7 +91,6 @@ end subroutine ReadNL subroutine SoilStateInitTimeConst(bounds, soilstate_inst, nlfilename) ! ! !USES: - use shr_kind_mod , only : r8 => shr_kind_r8 use shr_log_mod , only : errMsg => shr_log_errMsg use shr_infnan_mod , only : nan => shr_infnan_nan, assignment(=) use decompMod , only : bounds_type @@ -100,7 +103,7 @@ subroutine SoilStateInitTimeConst(bounds, soilstate_inst, nlfilename) use clm_varcon , only : zsoi, dzsoi, zisoi, spval use clm_varcon , only : secspday, pc, mu, denh2o, denice, grlnd use clm_varctl , only : use_cn, use_lch4, use_fates - use clm_varctl , only : iulog, fsurdat, paramfile, soil_layerstruct + use clm_varctl , only : iulog, fsurdat, paramfile, soil_layerstruct_predefined use landunit_varcon , only : istdlak, istwet, istsoil, istcrop, istice_mec use column_varcon , only : icol_roof, icol_sunwall, icol_shadewall, icol_road_perv, icol_road_imperv use fileutils , only : getfil @@ -140,7 +143,6 @@ subroutine SoilStateInitTimeConst(bounds, soilstate_inst, nlfilename) real(r8) :: tkm ! mineral conductivity real(r8) :: xksat ! maximum hydraulic conductivity of soil [mm/s] real(r8) :: clay,sand ! temporaries - real(r8) :: organic_max ! organic matter (kg/m3) where soil is assumed to act like peat integer :: dimid ! dimension id logical :: readvar type(file_desc_t) :: ncid ! netcdf id @@ -333,7 +335,11 @@ subroutine SoilStateInitTimeConst(bounds, soilstate_inst, nlfilename) g = col%gridcell(c) l = col%landunit(c) - if (lun%itype(l)==istwet .or. lun%itype(l)==istice_mec) then + ! istwet and istice_mec and + ! urban roof, sunwall, shadewall properties set to special value + if (lun%itype(l)==istwet .or. lun%itype(l)==istice_mec .or. & + (lun%urbpoi(l) .and. col%itype(c) /= icol_road_perv .and. & + col%itype(c) /= icol_road_imperv)) then do lev = 1,nlevgrnd soilstate_inst%bsw_col(c,lev) = spval @@ -358,48 +364,22 @@ subroutine SoilStateInitTimeConst(bounds, soilstate_inst, nlfilename) soilstate_inst%csol_col(c,lev)= spval end do - else if (lun%urbpoi(l) .and. (col%itype(c) /= icol_road_perv) .and. (col%itype(c) /= icol_road_imperv) )then - - ! Urban Roof, sunwall, shadewall properties set to special value - do lev = 1,nlevgrnd - soilstate_inst%watsat_col(c,lev) = spval - soilstate_inst%watfc_col(c,lev) = spval - soilstate_inst%bsw_col(c,lev) = spval - soilstate_inst%hksat_col(c,lev) = spval - soilstate_inst%sucsat_col(c,lev) = spval - soilstate_inst%watdry_col(c,lev) = spval - soilstate_inst%watopt_col(c,lev) = spval - soilstate_inst%bd_col(c,lev) = spval - if (lev <= nlevsoi) then - soilstate_inst%cellsand_col(c,lev) = spval - soilstate_inst%cellclay_col(c,lev) = spval - soilstate_inst%cellorg_col(c,lev) = spval - end if - end do - - do lev = 1,nlevgrnd - soilstate_inst%tkmg_col(c,lev) = spval - soilstate_inst%tksatu_col(c,lev) = spval - soilstate_inst%tkdry_col(c,lev) = spval - soilstate_inst%csol_col(c,lev) = spval - end do - else do lev = 1,nlevgrnd ! DML - this if statement could probably be removed and just the ! top part used for all soil layer structures - if ( soil_layerstruct /= '10SL_3.5m' )then ! apply soil texture from 10 layer input dataset + if ( soil_layerstruct_predefined /= '10SL_3.5m' )then ! apply soil texture from 10 layer input dataset if (lev .eq. 1) then clay = clay3d(g,1) sand = sand3d(g,1) - om_frac = organic3d(g,1)/organic_max + om_frac = organic3d(g,1)/organic_max else if (lev <= nlevsoi) then do j = 1,nlevsoifl-1 if (zisoi(lev) >= zisoifl(j) .AND. zisoi(lev) < zisoifl(j+1)) then clay = clay3d(g,j+1) sand = sand3d(g,j+1) - om_frac = organic3d(g,j+1)/organic_max + om_frac = organic3d(g,j+1)/organic_max endif end do else @@ -411,7 +391,7 @@ subroutine SoilStateInitTimeConst(bounds, soilstate_inst, nlfilename) if (lev <= nlevsoi) then ! duplicate clay and sand values from 10th soil layer clay = clay3d(g,lev) sand = sand3d(g,lev) - if ( organic_frac_squared )then + if (organic_frac_squared) then om_frac = (organic3d(g,lev)/organic_max)**2._r8 else om_frac = organic3d(g,lev)/organic_max @@ -423,25 +403,17 @@ subroutine SoilStateInitTimeConst(bounds, soilstate_inst, nlfilename) endif end if - if (lun%itype(l) == istdlak) then - - if (lev <= nlevsoi) then - soilstate_inst%cellsand_col(c,lev) = sand - soilstate_inst%cellclay_col(c,lev) = clay - soilstate_inst%cellorg_col(c,lev) = om_frac*organic_max - end if - - else if (lun%itype(l) /= istdlak) then ! soil columns of both urban and non-urban types + if (lun%urbpoi(l)) then + om_frac = 0._r8 ! No organic matter for urban + end if - if (lun%urbpoi(l)) then - om_frac = 0._r8 ! No organic matter for urban - end if + if (lev <= nlevsoi) then + soilstate_inst%cellsand_col(c,lev) = sand + soilstate_inst%cellclay_col(c,lev) = clay + soilstate_inst%cellorg_col(c,lev) = om_frac*organic_max + end if - if (lev <= nlevsoi) then - soilstate_inst%cellsand_col(c,lev) = sand - soilstate_inst%cellclay_col(c,lev) = clay - soilstate_inst%cellorg_col(c,lev) = om_frac*organic_max - end if + if (lun%itype(l) /= istdlak) then ! soil columns of both urban and non-urban types ! Note that the following properties are overwritten for urban impervious road ! layers that are not soil in SoilThermProp.F90 within SoilTemperatureMod.F90 diff --git a/src/main/clm_instMod.F90 b/src/main/clm_instMod.F90 index 7fa31ab4cd..406fa45224 100644 --- a/src/main/clm_instMod.F90 +++ b/src/main/clm_instMod.F90 @@ -306,7 +306,7 @@ subroutine clm_instInit(bounds) call soilhydrology_inst%Init(bounds, nlfilename, water_inst%waterstatebulk_inst, & use_aquifer_layer = use_aquifer_layer()) - call SoilHydrologyInitTimeConst(bounds, soilhydrology_inst) + call SoilHydrologyInitTimeConst(bounds, soilhydrology_inst, soilstate_inst) call saturated_excess_runoff_inst%Init(bounds) call infiltration_excess_runoff_inst%Init(bounds) diff --git a/src/main/clm_varctl.F90 b/src/main/clm_varctl.F90 index a6668a8148..0b100af1cf 100644 --- a/src/main/clm_varctl.F90 +++ b/src/main/clm_varctl.F90 @@ -264,7 +264,9 @@ module clm_varctl !---------------------------------------------------------- logical, public :: use_bedrock = .false. ! true => use spatially variable soil depth - character(len=16), public :: soil_layerstruct = '10SL_3.5m' + character(len=16), public :: soil_layerstruct_predefined = 'UNSET' + real(r8), public :: soil_layerstruct_userdefined(99) = rundef + integer, public :: soil_layerstruct_userdefined_nlevsoi = iundef !---------------------------------------------------------- ! plant hydraulic stress switch diff --git a/src/main/clm_varpar.F90 b/src/main/clm_varpar.F90 index c19caf5bd5..6d0e7b6bd6 100644 --- a/src/main/clm_varpar.F90 +++ b/src/main/clm_varpar.F90 @@ -11,7 +11,10 @@ module clm_varpar use clm_varctl , only: use_extralakelayers, use_vertsoilc use clm_varctl , only: use_century_decomp, use_c13, use_c14 use clm_varctl , only: iulog, use_crop, create_crop_landunit, irrigate - use clm_varctl , only: use_vichydro, soil_layerstruct + use clm_varctl , only: use_vichydro, rundef + use clm_varctl , only: soil_layerstruct_predefined + use clm_varctl , only: soil_layerstruct_userdefined + use clm_varctl , only: soil_layerstruct_userdefined_nlevsoi use clm_varctl , only: use_fates ! @@ -106,6 +109,7 @@ subroutine clm_varpar_init(actual_maxsoil_patches, actual_numcft) ! ! !LOCAL VARIABLES: ! + integer :: j ! loop index character(len=32) :: subname = 'clm_varpar_init' ! subroutine name !------------------------------------------------------------------------------ @@ -141,28 +145,56 @@ subroutine clm_varpar_init(actual_maxsoil_patches, actual_numcft) nlevsoifl = 10 nlevurb = 5 - if ( masterproc ) write(iulog, *) 'soil_layerstruct varpar ',soil_layerstruct - if ( soil_layerstruct == '10SL_3.5m' ) then - nlevsoi = nlevsoifl - nlevgrnd = 15 - else if ( soil_layerstruct == '23SL_3.5m' ) then - nlevsoi = 8 + nlev_equalspace - nlevgrnd = 15 + nlev_equalspace - else if ( soil_layerstruct == '49SL_10m' ) then - nlevsoi = 49 ! 10x10 + 9x100 + 30x300 = 1e4mm = 10m -! nlevsoi = 29 ! 10x10 + 9x100 + 10x300 = 4e3mm = 4m - nlevgrnd = nlevsoi+5 - else if ( soil_layerstruct == '20SL_8.5m' ) then - nlevsoi = 20 - nlevgrnd = nlevsoi+5 - else if ( soil_layerstruct == '5SL_3m' ) then - nlevsoi = 5 - nlevgrnd = 5 - else - write(iulog,*) subname//' ERROR: Unrecognized soil layer structure: ', trim(soil_layerstruct) - call shr_sys_abort(subname//' ERROR: Unrecognized soil layer structure') + + if ( masterproc ) write(iulog, *) 'soil_layerstruct_predefined varpar ', soil_layerstruct_predefined + if ( masterproc ) write(iulog, *) 'soil_layerstruct_userdefined varpar ', soil_layerstruct_userdefined + + if (soil_layerstruct_userdefined(1) /= rundef) then ! user defined soil layers + if (soil_layerstruct_predefined /= 'UNSET') then + write(iulog,*) subname//' ERROR: Both soil_layerstruct_predefined and soil_layer_userdefined have values' + call shr_sys_abort(subname//' ERROR: Cannot decide how to set the soil layer structure') + else + nlevgrnd = size(soil_layerstruct_userdefined) + ! loops backwards until it hits the last valid user-defined value + do j = nlevgrnd,1,-1 + if (soil_layerstruct_userdefined(j) /= rundef) then + exit + else + nlevgrnd = nlevgrnd - 1 + end if + end do + nlevsoi = soil_layerstruct_userdefined_nlevsoi ! read in namelist + if (nlevsoi >= nlevgrnd) then + write(iulog,*) subname//' ERROR: nlevsoi >= nlevgrnd; did you enter soil_layerstruct_userdefined_nlevsoi correctly in user_nl_clm?' + call shr_sys_abort(subname//' ERROR: nlevsoi must be less than nlevgrnd') + end if + end if + else ! pre-defined soil structure options + if ( soil_layerstruct_predefined == '10SL_3.5m' ) then + nlevsoi = nlevsoifl + nlevgrnd = 15 + else if ( soil_layerstruct_predefined == '23SL_3.5m' ) then + nlevsoi = 8 + nlev_equalspace + nlevgrnd = 15 + nlev_equalspace + else if ( soil_layerstruct_predefined == '49SL_10m' ) then + nlevsoi = 49 ! 10x10 + 9x100 + 30x300 = 1e4mm = 10m +! nlevsoi = 29 ! 10x10 + 9x100 + 10x300 = 4e3mm = 4m + nlevgrnd = nlevsoi+5 + else if ( soil_layerstruct_predefined == '20SL_8.5m' ) then + nlevsoi = 20 + nlevgrnd = nlevsoi+5 + else if ( soil_layerstruct_predefined == '4SL_2m' ) then + nlevsoi = 4 + nlevgrnd = 5 + else if (soil_layerstruct_predefined == 'UNSET') then + write(iulog,*) subname//' ERROR: Both soil_layerstruct_predefined and soil_layer_userdefined currently undefined' + call shr_sys_abort(subname//' ERROR: Cannot set the soil layer structure') + else + write(iulog,*) subname//' ERROR: Unrecognized pre-defined soil layer structure: ', trim(soil_layerstruct_predefined) + call shr_sys_abort(subname//' ERROR: Unrecognized pre-defined soil layer structure') + end if endif - if ( masterproc ) write(iulog, *) 'soil_layerstruct varpar ',soil_layerstruct,nlevsoi,nlevgrnd + if ( masterproc ) write(iulog, *) 'nlevsoi, nlevgrnd varpar ', nlevsoi, nlevgrnd if (use_vichydro) then nlayert = nlayer + (nlevgrnd -nlevsoi) diff --git a/src/main/controlMod.F90 b/src/main/controlMod.F90 index 59c4273277..3caea85a35 100644 --- a/src/main/controlMod.F90 +++ b/src/main/controlMod.F90 @@ -203,7 +203,8 @@ subroutine control_init( ) namelist /clm_inparm/ & clump_pproc, wrtdia, & create_crop_landunit, nsegspc, co2_ppmv, override_nsrest, & - albice, soil_layerstruct, subgridflag, & + albice, soil_layerstruct_predefined, soil_layerstruct_userdefined, & + soil_layerstruct_userdefined_nlevsoi, subgridflag, & irrigate, run_zero_weight_urban, all_active, & crop_fsat_equals_zero @@ -812,7 +813,9 @@ subroutine control_spmd() call mpi_bcast (scmlon, 1, MPI_REAL8,0, mpicom, ier) call mpi_bcast (co2_ppmv, 1, MPI_REAL8,0, mpicom, ier) call mpi_bcast (albice, 2, MPI_REAL8,0, mpicom, ier) - call mpi_bcast (soil_layerstruct,len(soil_layerstruct), MPI_CHARACTER, 0, mpicom, ier) + call mpi_bcast (soil_layerstruct_predefined,len(soil_layerstruct_predefined), MPI_CHARACTER, 0, mpicom, ier) + call mpi_bcast (soil_layerstruct_userdefined,size(soil_layerstruct_userdefined), MPI_REAL8, 0, mpicom, ier) + call mpi_bcast (soil_layerstruct_userdefined_nlevsoi, 1, MPI_INTEGER, 0, mpicom, ier) ! snow pack variables call mpi_bcast (nlevsno, 1, MPI_INTEGER, 0, mpicom, ier) @@ -1028,7 +1031,9 @@ subroutine control_print () end if write(iulog,*) ' land-ice albedos (unitless 0-1) = ', albice - write(iulog,*) ' soil layer structure = ', soil_layerstruct + write(iulog,*) ' pre-defined soil layer structure = ', soil_layerstruct_predefined + write(iulog,*) ' user-defined soil layer structure = ', soil_layerstruct_userdefined + write(iulog,*) ' user-defined number of soil layers = ', soil_layerstruct_userdefined_nlevsoi write(iulog,*) ' plant hydraulic stress = ', use_hydrstress write(iulog,*) ' dynamic roots = ', use_dynroot if (nsrest == nsrContinue) then diff --git a/src/main/initVerticalMod.F90 b/src/main/initVerticalMod.F90 index f4fc4cea4b..eef64a6aa5 100644 --- a/src/main/initVerticalMod.F90 +++ b/src/main/initVerticalMod.F90 @@ -17,7 +17,8 @@ module initVerticalMod use clm_varpar , only : nlevsoi, nlevsoifl, nlevurb use clm_varctl , only : fsurdat, iulog use clm_varctl , only : use_vancouver, use_mexicocity, use_vertsoilc, use_extralakelayers - use clm_varctl , only : use_bedrock, soil_layerstruct + use clm_varctl , only : use_bedrock, rundef + use clm_varctl , only : soil_layerstruct_predefined, soil_layerstruct_userdefined use clm_varctl , only : use_fates use clm_varcon , only : zlak, dzlak, zsoi, dzsoi, zisoi, dzsoi_decomp, spval, ispval, grlnd use column_varcon , only : icol_roof, icol_sunwall, icol_shadewall, is_hydrologically_active @@ -157,6 +158,7 @@ subroutine initVertical(bounds, glc_behavior, snow_depth, thick_wall, thick_roof integer :: ier ! error status real(r8) :: scalez = 0.025_r8 ! Soil layer thickness discretization (m) real(r8) :: thick_equal = 0.2 + character(len=20) :: calc_method ! soil layer calculation method real(r8) ,pointer :: zbedrock_in(:) ! read in - z_bedrock real(r8) ,pointer :: lakedepth_in(:) ! read in - lakedepth real(r8), allocatable :: zurb_wall(:,:) ! wall (layer node depth) @@ -213,41 +215,52 @@ subroutine initVertical(bounds, glc_behavior, snow_depth, thick_wall, thick_roof ! Soil layers and interfaces (assumed same for all non-lake patches) ! "0" refers to soil surface and "nlevsoi" refers to the bottom of model soil - if ( soil_layerstruct == '10SL_3.5m' ) then - do j = 1, nlevgrnd - zsoi(j) = scalez*(exp(0.5_r8*(j-0.5_r8))-1._r8) !node depths - enddo - - dzsoi(1) = 0.5_r8*(zsoi(1)+zsoi(2)) !thickness b/n two interfaces - do j = 2,nlevgrnd-1 - dzsoi(j)= 0.5_r8*(zsoi(j+1)-zsoi(j-1)) - enddo - dzsoi(nlevgrnd) = zsoi(nlevgrnd)-zsoi(nlevgrnd-1) - - zisoi(0) = 0._r8 - do j = 1, nlevgrnd-1 - zisoi(j) = 0.5_r8*(zsoi(j)+zsoi(j+1)) !interface depths - enddo - zisoi(nlevgrnd) = zsoi(nlevgrnd) + 0.5_r8*dzsoi(nlevgrnd) + if (soil_layerstruct_predefined == '10SL_3.5m' .or. soil_layerstruct_predefined == '23SL_3.5m') then + calc_method = 'node-based' ! node-based followed by error check + if (soil_layerstruct_userdefined(1) /= rundef) then + write(iulog,*) subname//' ERROR: Both soil_layerstruct_predefined and soil_layer_userdefined have values' + call shr_sys_abort(subname//' ERROR: Cannot decide how to set the soil layer structure') + end if + ! thickness-based (part 1) and error check + else if (soil_layerstruct_predefined == '49SL_10m' .or. & + soil_layerstruct_predefined == '20SL_8.5m' .or. & + soil_layerstruct_predefined == '4SL_2m') then + calc_method = 'thickness-based' + if (soil_layerstruct_userdefined(1) /= rundef) then + write(iulog,*) subname//' ERROR: Both soil_layerstruct_predefined and soil_layer_userdefined have values' + call shr_sys_abort(subname//' ERROR: Cannot decide how to set the soil layer structure') + end if + ! thickness-based (part 2) and error check + else if (soil_layerstruct_userdefined(1) /= rundef) then + calc_method = 'thickness-based' + if (soil_layerstruct_predefined /= 'UNSET') then + write(iulog,*) subname//' ERROR: Both soil_layerstruct_predefined and soil_layer_userdefined have values' + call shr_sys_abort(subname//' ERROR: Cannot decide how to set the soil layer structure') + end if + else ! error check + write(iulog,*) subname//' ERROR: Unrecognized pre-defined and user-defined soil layer structures: ', trim(soil_layerstruct_predefined), soil_layerstruct_userdefined + call endrun(subname//' ERROR: Unrecognized soil layer structure') + end if - else if ( soil_layerstruct == '23SL_3.5m' )then - ! Soil layer structure that starts with standard exponential - ! and then has several evenly spaced layers, then finishes off exponential. - ! this allows the upper soil to behave as standard, but then continues - ! with higher resolution to a deeper depth, so that, for example, permafrost - ! dynamics are not lost due to an inability to resolve temperature, moisture, - ! and biogeochemical dynamics at the base of the active layer - do j = 1, toplev_equalspace + if (calc_method == 'node-based') then + do j = 1, nlevgrnd zsoi(j) = scalez*(exp(0.5_r8*(j-0.5_r8))-1._r8) !node depths enddo - do j = toplev_equalspace+1,toplev_equalspace + nlev_equalspace - zsoi(j) = zsoi(j-1) + thick_equal - enddo - - do j = toplev_equalspace + nlev_equalspace +1, nlevgrnd - zsoi(j) = scalez*(exp(0.5_r8*((j - nlev_equalspace)-0.5_r8))-1._r8) + nlev_equalspace * thick_equal - enddo + if (soil_layerstruct_predefined == '23SL_3.5m') then + ! Soil layer structure that starts with standard exponential, + ! then has several evenly spaced layers and finishes off exponential. + ! This allows the upper soil to behave as standard, but then continues + ! with higher resolution to a deeper depth, so that, e.g., permafrost + ! dynamics are not lost due to an inability to resolve temperature, + ! moisture, and biogeochemical dynamics at the base of the active layer + do j = toplev_equalspace + 1, toplev_equalspace + nlev_equalspace + zsoi(j) = zsoi(j-1) + thick_equal + enddo + do j = toplev_equalspace + nlev_equalspace + 1, nlevgrnd + zsoi(j) = scalez * (exp(0.5_r8 * (j - nlev_equalspace - 0.5_r8)) - 1._r8) + nlev_equalspace * thick_equal + enddo + end if ! soil_layerstruct_predefined == '23SL_3.5m' dzsoi(1) = 0.5_r8*(zsoi(1)+zsoi(2)) !thickness b/n two interfaces do j = 2,nlevgrnd-1 @@ -257,24 +270,50 @@ subroutine initVertical(bounds, glc_behavior, snow_depth, thick_wall, thick_roof zisoi(0) = 0._r8 do j = 1, nlevgrnd-1 - zisoi(j) = 0.5_r8*(zsoi(j)+zsoi(j+1)) !interface depths + zisoi(j) = 0.5_r8*(zsoi(j)+zsoi(j+1)) !interface depths enddo zisoi(nlevgrnd) = zsoi(nlevgrnd) + 0.5_r8*dzsoi(nlevgrnd) - else if ( soil_layerstruct == '49SL_10m' ) then - !scs: 10 meter soil column, nlevsoi set to 49 in clm_varpar - do j = 1,10 - dzsoi(j)= 1.e-2_r8 !10mm layers - enddo - do j = 11,19 - dzsoi(j)= 1.e-1_r8 !100 mm layers - enddo - do j = 20,nlevsoi+1 !300 mm layers - dzsoi(j)= 3.e-1_r8 - enddo - do j = nlevsoi+2,nlevgrnd !10 meter bedrock layers - dzsoi(j)= 10._r8 - enddo + else if (calc_method == 'thickness-based') then + if (soil_layerstruct_userdefined(1) /= rundef) then + do j = 1, nlevgrnd + ! read dzsoi from user-entered namelist vector + dzsoi(j) = soil_layerstruct_userdefined(j) + end do + else if (soil_layerstruct_predefined == '49SL_10m') then + !scs: 10 meter soil column, nlevsoi set to 49 in clm_varpar + do j = 1, 10 + dzsoi(j) = 1.e-2_r8 ! 10-mm layers + enddo + do j = 11, 19 + dzsoi(j) = 1.e-1_r8 ! 100-mm layers + enddo + do j = 20, nlevsoi+1 ! 300-mm layers + dzsoi(j) = 3.e-1_r8 + enddo + do j = nlevsoi+2,nlevgrnd ! 10-m bedrock layers + dzsoi(j) = 10._r8 + enddo + else if (soil_layerstruct_predefined == '20SL_8.5m') then + do j = 1, 4 ! linear increase in layer thickness of... + dzsoi(j) = j * 0.02_r8 ! ...2 cm each layer + enddo + do j = 5, 13 + dzsoi(j) = dzsoi(4) + (j - 4) * 0.04_r8 ! ...4 cm each layer + enddo + do j = 14, nlevsoi + dzsoi(j) = dzsoi(13) + (j - 13) * 0.10_r8 ! ...10 cm each layer + enddo + do j = nlevsoi + 1, nlevgrnd ! bedrock layers + dzsoi(j) = dzsoi(nlevsoi) + (((j - nlevsoi) * 25._r8)**1.5_r8) / 100._r8 + enddo + else if (soil_layerstruct_predefined == '4SL_2m') then + dzsoi(1) = 0.1_r8 + dzsoi(2) = 0.3_r8 + dzsoi(3) = 0.6_r8 + dzsoi(4) = 1.0_r8 + dzsoi(5) = 1.0_r8 + end if ! thickness-based options zisoi(0) = 0._r8 do j = 1,nlevgrnd @@ -285,54 +324,17 @@ subroutine initVertical(bounds, glc_behavior, snow_depth, thick_wall, thick_roof zsoi(j) = 0.5*(zisoi(j-1) + zisoi(j)) enddo - else if ( soil_layerstruct == '20SL_8.5m' ) then - do j = 1,4 - dzsoi(j)= j*0.02_r8 ! linear increase in layer thickness of 2cm each layer - enddo - do j = 5,13 - dzsoi(j)= dzsoi(4)+(j-4)*0.04_r8 ! linear increase in layer thickness of 2cm each layer - enddo - do j = 14,nlevsoi - dzsoi(j)= dzsoi(13)+(j-13)*0.10_r8 ! linear increase in layer thickness of 2cm each layer - enddo - do j = nlevsoi+1,nlevgrnd !bedrock layers - dzsoi(j)= dzsoi(nlevsoi)+(((j-nlevsoi)*25._r8)**1.5_r8)/100._r8 ! bedrock layers - enddo - - zisoi(0) = 0._r8 - do j = 1,nlevgrnd - zisoi(j)= sum(dzsoi(1:j)) - enddo - - do j = 1, nlevgrnd - zsoi(j) = 0.5*(zisoi(j-1) + zisoi(j)) - enddo - else if ( soil_layerstruct == '5SL_3m' ) then - dzsoi(1)= 0.1_r8 - dzsoi(2)= 0.3_r8 - dzsoi(3)= 0.6_r8 - dzsoi(4)= 1.0_r8 - dzsoi(5)= 1.0_r8 - - zisoi(0) = 0._r8 - do j = 1,nlevgrnd - zisoi(j)= sum(dzsoi(1:j)) - enddo - - do j = 1, nlevgrnd - zsoi(j) = 0.5*(zisoi(j-1) + zisoi(j)) - enddo - else - write(iulog,*) subname//' ERROR: Unrecognized soil layer structure: ', trim(soil_layerstruct) - call endrun(subname//' ERROR: Unrecognized soil layer structure') - end if + else ! error check + write(iulog,*) subname//' ERROR: Unrecognized calc_method: ', trim(calc_method) + call endrun(subname//' ERROR: Unrecognized calc_method') + end if ! calc_method is node-based or thickness-based ! define a vertical grid spacing such that it is the normal dzsoi if ! nlevdecomp =nlevgrnd, or else 1 meter if (use_vertsoilc) then dzsoi_decomp = dzsoi !thickness b/n two interfaces else - dzsoi_decomp(1) = 1. + dzsoi_decomp(1) = 1._r8 end if if (masterproc) then