Skip to content

Commit c74d3f7

Browse files
slevis-lmwgslevis-lmwg
authored
Merge pull request #3893 from swensosc/wetland_update
branch_tags/wetlands.n01.ctsm5.4.029: Wetland update Merging Sean Swenson's wetland_update branch to the wetlands branch and tagging. The wetlands branch - will get periodic updates from master as deemed necessary - will also include new paramfile(s) with methane-related updates and updates from Linnia - is intended for Will's simulations for WIE-MIP and such - will get merged to master when deemed appropriate
2 parents ec0f752 + 3af9c70 commit c74d3f7

File tree

10 files changed

+1283
-491
lines changed

10 files changed

+1283
-491
lines changed

‎bld/CLMBuildNamelist.pm‎

Lines changed: 8 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -3497,11 +3497,14 @@ sub setup_logic_methane {
34973497
add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'finundation_method',
34983498
'use_cn'=>$nl_flags->{'use_cn'}, 'use_fates'=>$nl_flags->{'use_fates'} );
34993499
my $finundation_method = remove_leading_and_trailing_quotes($nl->get_value('finundation_method' ));
3500-
add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'stream_fldfilename_ch4finundated',
3501-
'finundation_method'=>$finundation_method);
3502-
if ($opts->{'driver'} eq "nuopc" ) {
3503-
add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'stream_meshfile_ch4finundated',
3504-
'finundation_method'=>$finundation_method);
3500+
# prognostic inundation does not require an input stream; other methods do
3501+
if($finundation_method ne 'h2osfc') {
3502+
add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'stream_fldfilename_ch4finundated',
3503+
'finundation_method'=>$finundation_method);
3504+
if ($opts->{'driver'} eq "nuopc" ) {
3505+
add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'stream_meshfile_ch4finundated',
3506+
'finundation_method'=>$finundation_method);
3507+
}
35053508
}
35063509
add_default($opts, $nl_flags->{'inputdata_rootdir'}, $definition, $defaults, $nl, 'use_aereoxid_prog',
35073510
'use_cn'=>$nl_flags->{'use_cn'}, 'use_fates'=>$nl_flags->{'use_fates'} );

‎src/biogeochem/ch4Mod.F90‎

Lines changed: 843 additions & 235 deletions
Large diffs are not rendered by default.

‎src/biogeophys/HillslopeHydrologyMod.F90‎

Lines changed: 82 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -53,6 +53,8 @@ module HillslopeHydrologyMod
5353
integer, private, parameter :: soil_profile_set_lowland_upland = 2
5454
integer, private, parameter :: soil_profile_linear = 3
5555

56+
real(r8), private, parameter :: floodplain_slope = 1e-3_r8
57+
5658
!-----------------------------------------------------------------------
5759

5860
contains
@@ -962,17 +964,32 @@ subroutine HillslopeStreamOutflow(bounds, &
962964

963965
integer :: c, l, g, i, j
964966
integer :: nstep
965-
real(r8) :: dtime ! land model time step (sec)
967+
real(r8) :: dtime ! land model time step (sec)
968+
real(r8) :: cross_sectional_area_channel ! cross sectional area of stream channel (m2)
966969
real(r8) :: cross_sectional_area ! cross sectional area of stream water (m2)
967970
real(r8) :: stream_depth ! depth of stream water (m)
968971
real(r8) :: hydraulic_radius ! cross sectional area divided by wetted perimeter (m)
969972
real(r8) :: flow_velocity ! flow velocity (m/s)
970973
real(r8) :: overbank_area ! area of water above bankfull (m2)
974+
real(r8) :: overbank_depth ! depth of water above bankfull (m)
975+
real(r8) :: dynamic_viscosity
976+
real(r8) :: bankfull_flow_velocity
977+
real(r8) :: kinematic_viscosity
978+
real(r8) :: length_scale
979+
real(r8) :: reynolds
971980
real(r8), parameter :: manning_roughness = 0.03_r8 ! manning roughness
972981
real(r8), parameter :: manning_exponent = 0.667_r8 ! manning exponent
973982

974983
integer, parameter :: overbank_method = 1 ! method to treat overbank stream storage; 1 = increase dynamic slope, 2 = increase flow area cross section, 3 = remove instantaneously
975984
logical :: active_stream
985+
real(r8), parameter :: Acoef = 0.02939 !mPa·s
986+
real(r8), parameter :: Bcoef = 507.88 !K
987+
real(r8), parameter :: Ccoef = 149.30 !K
988+
real(r8), parameter :: Tavg = 293.15 !20C + 273.15
989+
real(r8), parameter :: rhow = 1e3
990+
real(r8), parameter :: reynolds_thresh = 1000
991+
real(r8), parameter :: stream_conductivity = 1e-3_r8
992+
real(r8), parameter :: water_surface_slope = 5e-5_r8
976993
character(len=*), parameter :: subname = 'HillslopeStreamOutflow'
977994

978995
!-----------------------------------------------------------------------
@@ -998,25 +1015,71 @@ subroutine HillslopeStreamOutflow(bounds, &
9981015
if (lun%active(l) .and. active_stream) then
9991016
! Streamflow calculated from Manning equation
10001017
if (streamflow_method == streamflow_manning) then
1018+
cross_sectional_area_channel = lun%stream_channel_width(l)*lun%stream_channel_depth(l)
1019+
10011020
cross_sectional_area = stream_water_volume(l) &
10021021
/lun%stream_channel_length(l)
1003-
stream_depth = cross_sectional_area &
1004-
/lun%stream_channel_width(l)
1005-
hydraulic_radius = cross_sectional_area &
1006-
/(lun%stream_channel_width(l) + 2*stream_depth)
1022+
1023+
! calculate hydraulic radius
1024+
if(cross_sectional_area <= cross_sectional_area_channel) then
1025+
stream_depth = cross_sectional_area &
1026+
/lun%stream_channel_width(l)
1027+
overbank_depth = 0._r8
1028+
hydraulic_radius = cross_sectional_area &
1029+
/(lun%stream_channel_width(l) + 2*stream_depth)
1030+
else ! overbank conditions exist
1031+
stream_depth = lun%stream_channel_depth(l)
1032+
overbank_area = cross_sectional_area &
1033+
-cross_sectional_area_channel
1034+
! use small positive slope for floodplain
1035+
overbank_depth = sqrt(floodplain_slope*overbank_area)
1036+
1037+
hydraulic_radius = cross_sectional_area &
1038+
/(lun%stream_channel_width(l) &
1039+
+ 2*stream_depth &
1040+
+2*overbank_depth/floodplain_slope) ! sin ~ tan for small values
1041+
endif
10071042

10081043
if (hydraulic_radius <= 0._r8) then
10091044
volumetric_streamflow(l) = 0._r8
10101045
else
1011-
flow_velocity = (hydraulic_radius)**manning_exponent &
1046+
! check reynolds number using bankfull channel properties (static)
1047+
bankfull_flow_velocity = (lun%stream_channel_depth(l))**manning_exponent &
10121048
* sqrt(lun%stream_channel_slope(l)) &
10131049
/ manning_roughness
1050+
1051+
! check reynolds number dynamically
1052+
flow_velocity = (hydraulic_radius)**manning_exponent &
1053+
* sqrt(lun%stream_channel_slope(l)) &
1054+
/ manning_roughness
1055+
1056+
! Vogel-Fulcher-Tammann equation
1057+
dynamic_viscosity = Acoef * exp((Bcoef/(Tavg-Ccoef)))
1058+
kinematic_viscosity = dynamic_viscosity/rhow
1059+
1060+
! let length scale equal bankfull depth
1061+
!length_scale = lun%stream_channel_depth(l)
1062+
!reynolds = bankfull_flow_velocity*length_scale/kinematic_viscosity
1063+
! let length scale equal dynamic depth
1064+
length_scale = stream_depth + overbank_depth
1065+
reynolds = flow_velocity*length_scale/kinematic_viscosity
1066+
1067+
if( reynolds < reynolds_thresh) then
1068+
! use exponent = 2 instead of 0.66, and specified water surface slope
1069+
flow_velocity = (hydraulic_radius)**2 &
1070+
* water_surface_slope &
1071+
/ (manning_roughness)
1072+
else
1073+
flow_velocity = (hydraulic_radius)**manning_exponent &
1074+
* sqrt(lun%stream_channel_slope(l)) &
1075+
/ manning_roughness
1076+
endif
10141077
! overbank flow
1015-
if (stream_depth > lun%stream_channel_depth(l)) then
1078+
if (overbank_depth > 0._r8) then
10161079
if (overbank_method == 1) then
1017-
! try increasing dynamic slope
1018-
volumetric_streamflow(l) = cross_sectional_area * flow_velocity &
1019-
*(stream_depth/lun%stream_channel_depth(l))
1080+
! flow velocity already accounts for overbank conditions
1081+
1082+
volumetric_streamflow(l) = cross_sectional_area * flow_velocity
10201083
else if (overbank_method == 2) then
10211084
! try increasing flow area cross section
10221085
overbank_area = (stream_depth -lun%stream_channel_depth(l)) * 30._r8 * lun%stream_channel_width(l)
@@ -1080,6 +1143,7 @@ subroutine HillslopeUpdateStreamWater(bounds, waterstatebulk_inst, &
10801143
real(r8) :: qflx_drain_vol ! volumetric saturated drainage (m3/s)
10811144
real(r8) :: dtime ! land model time step (sec)
10821145
logical :: active_stream
1146+
real(r8) :: overbank_area, overbank_depth
10831147

10841148
character(len=*), parameter :: subname = 'HillslopeUpdateStreamWater'
10851149

@@ -1138,6 +1202,14 @@ subroutine HillslopeUpdateStreamWater(bounds, waterstatebulk_inst, &
11381202
/lun%stream_channel_length(l) &
11391203
/lun%stream_channel_width(l)
11401204

1205+
! recalculate for floodplain approximation
1206+
if (stream_water_depth(l)>lun%stream_channel_depth(l)) then
1207+
overbank_area = (stream_water_volume(l) &
1208+
/lun%stream_channel_length(l)) &
1209+
-lun%stream_channel_width(l)*lun%stream_channel_depth(l)
1210+
overbank_depth = sqrt(floodplain_slope*overbank_area)
1211+
stream_water_depth(l) = lun%stream_channel_depth(l) + overbank_depth
1212+
endif
11411213
end if
11421214
enddo
11431215

‎src/biogeophys/HydrologyNoDrainageMod.F90‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -322,7 +322,7 @@ subroutine HydrologyNoDrainage(bounds, &
322322

323323
call UpdateH2osfc(bounds, num_hydrologyc, filter_hydrologyc, &
324324
infiltration_excess_runoff_inst, &
325-
energyflux_inst, soilhydrology_inst, &
325+
energyflux_inst, soilhydrology_inst, soilstate_inst, &
326326
b_waterflux_inst, b_waterstate_inst, b_waterdiagnostic_inst)
327327

328328
call Infiltration(bounds, num_hydrologyc, filter_hydrologyc, &

‎src/biogeophys/SaturatedExcessRunoffMod.F90‎

Lines changed: 6 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -267,14 +267,17 @@ subroutine SaturatedExcessRunoff (this, bounds, num_hydrologyc, filter_hydrology
267267
endif
268268

269269
! ------------------------------------------------------------------------
270-
! Set fsat to zero for upland hillslope columns
270+
! Set fsat to zero for hillslope columns
271+
!
272+
! hillslope explicitly represents lateral water movement and
273+
! convergence, so the implicit topmodel approach is unnecessary
274+
!
271275
! ------------------------------------------------------------------------
272276
if (hillslope_fsat_equals_zero) then
273277
do fc = 1, num_hydrologyc
274278
c = filter_hydrologyc(fc)
275279
if(col%is_hillslope_column(c) .and. col%active(c)) then
276-
! Set fsat to zero for upland columns
277-
if (col%cold(c) /= ispval) fsat(c) = 0._r8
280+
fsat(c) = 0._r8
278281
endif
279282
end do
280283
endif

‎src/biogeophys/SoilHydrologyInitTimeConstMod.F90‎

Lines changed: 11 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -77,6 +77,7 @@ subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst, soilstate_inst
7777
use column_varcon , only : icol_shadewall, icol_road_perv, icol_road_imperv, icol_roof, icol_sunwall
7878
use fileutils , only : getfil
7979
use ncdio_pio , only : file_desc_t, ncd_io, ncd_pio_openfile, ncd_pio_closefile
80+
use SurfaceWaterMod , only : pc_hillslope
8081
!
8182
! !ARGUMENTS:
8283
type(bounds_type) , intent(in) :: bounds
@@ -99,6 +100,7 @@ subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst, soilstate_inst
99100
real(r8), pointer :: sandcol (:,:) ! column level sand fraction for calculating VIC parameters
100101
real(r8), pointer :: claycol (:,:) ! column level clay fraction for calculating VIC parameters
101102
real(r8), pointer :: om_fraccol (:,:) ! column level organic matter fraction for calculating VIC parameters
103+
real(r8) :: pc ! column-level pc parameter
102104
!-----------------------------------------------------------------------
103105
! Initialize VIC variables
104106

@@ -215,12 +217,19 @@ subroutine SoilHydrologyInitTimeConst(bounds, soilhydrology_inst, soilstate_inst
215217

216218
! determine h2osfc threshold ("fill & spill" concept)
217219
! set to zero for no h2osfc (w/frac_infclust =large)
218-
220+
221+
if(col%is_hillslope_column(c)) then
222+
pc = max(pc_hillslope(col%hill_slope(c)), 0.01_r8)
223+
else
224+
! use standard value
225+
pc = params_inst%pc
226+
endif
227+
219228
soilhydrology_inst%h2osfc_thresh_col(c) = 0._r8
220229
if (micro_sigma(c) > 1.e-6_r8 .and. (soilhydrology_inst%h2osfcflag /= 0)) then
221230
d = 0.0_r8
222231
do p = 1,4
223-
fd = 0.5_r8*(1.0_r8+shr_spfn_erf(d/(micro_sigma(c)*sqrt(2.0_r8)))) - params_inst%pc
232+
fd = 0.5_r8*(1.0_r8+shr_spfn_erf(d/(micro_sigma(c)*sqrt(2.0_r8)))) - pc
224233
dfdd = exp(-d**2/(2.0_r8*micro_sigma(c)**2))/(micro_sigma(c)*sqrt(2.0_r8*shr_const_pi))
225234
d = d - fd/dfdd
226235
enddo

0 commit comments

Comments
 (0)