From 947a32b1e18dc0546a7b63137a8a3fae99e42c09 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 15 Oct 2024 10:05:26 -0700 Subject: [PATCH 01/39] added debugging print fields to snow balance check --- components/elm/src/biogeophys/BalanceCheckMod.F90 | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index 7737f2adbe3d..490f8c99886e 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -553,6 +553,13 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & write(iulog,*)'qflx_snwcp_ice = ',qflx_snwcp_ice(indexc)*dtime write(iulog,*)'qflx_snwcp_liq = ',qflx_snwcp_liq(indexc)*dtime write(iulog,*)'qflx_sl_top_soil = ',qflx_sl_top_soil(indexc)*dtime + write(iulog,*)'qflx_snow_melt = ',qflx_snow_melt(indexc)*dtime + write(iulog,*)'frac_sno_eff = ',frac_sno_eff(indexc)*dtime + write(iulog,*)'qflx_rain_grnd_col= ',qflx_rain_grnd_col(indexc)*dtime + write(iulog,*)'qflx_snow_grnd_col= ',qflx_snow_grnd_col(indexc)*dtime + write(iulog,*)'qflx_h2osfc_to_ice= ',qflx_h2osfc_to_ice(indexc)*dtime + write(iulog,*)'qflx_snow_h2osfc= ',qflx_snow_h2osfc(indexc)*dtime + if (create_glacier_mec_landunit) then write(iulog,*)'qflx_glcice_frz = ',qflx_glcice_frz(indexc)*dtime end if From c1dfdb47b835ee6aa54a8a74ee15568c093230d3 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 22 Oct 2024 13:59:36 -0700 Subject: [PATCH 02/39] debugging snow balance error, fix snowcompaction equ --- .../elm/src/biogeophys/BalanceCheckMod.F90 | 69 +++++++++++++++---- .../elm/src/biogeophys/SnowHydrologyMod.F90 | 2 +- 2 files changed, 55 insertions(+), 16 deletions(-) diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index 490f8c99886e..0fd26bec1e64 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -534,6 +534,40 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & ' col_pp%itype= ',col_pp%itype(indexc), & ' lun_pp%itype= ',lun_pp%itype(col_pp%landunit(indexc)), & ' errh2osno= ',errh2osno(indexc) + write(iulog,*)'nstep = ',nstep + g = col_pp%gridcell(indexc) + l = col_pp%landunit(indexc) + write(iulog,*)'errh2osno = ',errh2osno(indexc) + write(iulog,*)'snl = ',col_pp%snl(indexc) + write(iulog,*)'h2osno = ',h2osno(indexc) + write(iulog,*)'h2osno_old = ',h2osno_old(indexc) + write(iulog,*)'snow_sources = ',snow_sources(indexc) + write(iulog,*)'snow_sinks = ',snow_sinks(indexc) + write(iulog,*)'qflx_prec_grnd = ',qflx_prec_grnd(indexc) + write(iulog,*)'qflx_sub_snow = ',qflx_sub_snow(indexc) + write(iulog,*)'qflx_evap_grnd = ',qflx_evap_grnd(indexc) + write(iulog,*)'qflx_top_soil = ',qflx_top_soil(indexc) + write(iulog,*)'qflx_dew_snow = ',qflx_dew_snow(indexc) + write(iulog,*)'qflx_dew_grnd = ',qflx_dew_grnd(indexc) + write(iulog,*)'qflx_snwcp_ice = ',qflx_snwcp_ice(indexc) + write(iulog,*)'qflx_snwcp_liq = ',qflx_snwcp_liq(indexc) + write(iulog,*)'qflx_sl_top_soil = ',qflx_sl_top_soil(indexc) + write(iulog,*)'qflx_snow_melt = ',qflx_snow_melt(indexc) + write(iulog,*)'frac_sno_eff = ',frac_sno_eff(indexc) + write(iulog,*)'qflx_rain_grnd_col= ',qflx_rain_grnd_col(indexc) + write(iulog,*)'qflx_snow_grnd_col= ',qflx_snow_grnd_col(indexc) + write(iulog,*)'qflx_h2osfc_to_ice= ',qflx_h2osfc_to_ice(indexc) + write(iulog,*)'qflx_snow_h2osfc= ',qflx_snow_h2osfc(indexc)*dtime + write(iulog,*)'(glc_dyn_runoff_routing(g))',(glc_dyn_runoff_routing(g)) + write(iulog,*)'lun_pp%itype(l))',lun_pp%itype(l) + write(iulog,*)'lun_pp%itype(l) == istsoil)',lun_pp%itype(l) == istsoil + write(iulog,*)'lun_pp%itype(l) == istcrop)',lun_pp%itype(l) == istcrop + write(iulog,*)'lun_pp%itype(l) == istwet)',lun_pp%itype(l) == istwet + + if (create_glacier_mec_landunit) then + write(iulog,*)'qflx_glcice_frz = ',qflx_glcice_frz(indexc) + end if + if (abs(errh2osno(indexc)) > 1.e-4_r8 .and. (nstep > 2) ) then write(iulog,*)'elm model is stopping - error is greater than 1e-4 (mm)' @@ -544,24 +578,29 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & write(iulog,*)'h2osno_old = ',h2osno_old(indexc) write(iulog,*)'snow_sources = ',snow_sources(indexc) write(iulog,*)'snow_sinks = ',snow_sinks(indexc) - write(iulog,*)'qflx_prec_grnd = ',qflx_prec_grnd(indexc)*dtime - write(iulog,*)'qflx_sub_snow = ',qflx_sub_snow(indexc)*dtime - write(iulog,*)'qflx_evap_grnd = ',qflx_evap_grnd(indexc)*dtime - write(iulog,*)'qflx_top_soil = ',qflx_top_soil(indexc)*dtime - write(iulog,*)'qflx_dew_snow = ',qflx_dew_snow(indexc)*dtime - write(iulog,*)'qflx_dew_grnd = ',qflx_dew_grnd(indexc)*dtime - write(iulog,*)'qflx_snwcp_ice = ',qflx_snwcp_ice(indexc)*dtime - write(iulog,*)'qflx_snwcp_liq = ',qflx_snwcp_liq(indexc)*dtime - write(iulog,*)'qflx_sl_top_soil = ',qflx_sl_top_soil(indexc)*dtime - write(iulog,*)'qflx_snow_melt = ',qflx_snow_melt(indexc)*dtime - write(iulog,*)'frac_sno_eff = ',frac_sno_eff(indexc)*dtime - write(iulog,*)'qflx_rain_grnd_col= ',qflx_rain_grnd_col(indexc)*dtime - write(iulog,*)'qflx_snow_grnd_col= ',qflx_snow_grnd_col(indexc)*dtime - write(iulog,*)'qflx_h2osfc_to_ice= ',qflx_h2osfc_to_ice(indexc)*dtime + write(iulog,*)'qflx_prec_grnd = ',qflx_prec_grnd(indexc) + write(iulog,*)'qflx_sub_snow = ',qflx_sub_snow(indexc) + write(iulog,*)'qflx_evap_grnd = ',qflx_evap_grnd(indexc) + write(iulog,*)'qflx_top_soil = ',qflx_top_soil(indexc) + write(iulog,*)'qflx_dew_snow = ',qflx_dew_snow(indexc) + write(iulog,*)'qflx_dew_grnd = ',qflx_dew_grnd(indexc) + write(iulog,*)'qflx_snwcp_ice = ',qflx_snwcp_ice(indexc) + write(iulog,*)'qflx_snwcp_liq = ',qflx_snwcp_liq(indexc) + write(iulog,*)'qflx_sl_top_soil = ',qflx_sl_top_soil(indexc) + write(iulog,*)'qflx_snow_melt = ',qflx_snow_melt(indexc) + write(iulog,*)'frac_sno_eff = ',frac_sno_eff(indexc) + write(iulog,*)'qflx_rain_grnd_col= ',qflx_rain_grnd_col(indexc) + write(iulog,*)'qflx_snow_grnd_col= ',qflx_snow_grnd_col(indexc) + write(iulog,*)'qflx_h2osfc_to_ice= ',qflx_h2osfc_to_ice(indexc) write(iulog,*)'qflx_snow_h2osfc= ',qflx_snow_h2osfc(indexc)*dtime + write(iulog,*)'(glc_dyn_runoff_routing(g))',(glc_dyn_runoff_routing(g)) + write(iulog,*)'lun_pp%itype(l))',lun_pp%itype(l) + write(iulog,*)'lun_pp%itype(l) == istsoil)',lun_pp%itype(l) == istsoil + write(iulog,*)'lun_pp%itype(l) == istcrop)',lun_pp%itype(l) == istcrop + write(iulog,*)'lun_pp%itype(l) == istwet)',lun_pp%itype(l) == istwet if (create_glacier_mec_landunit) then - write(iulog,*)'qflx_glcice_frz = ',qflx_glcice_frz(indexc)*dtime + write(iulog,*)'qflx_glcice_frz = ',qflx_glcice_frz(indexc) end if write(iulog,*)'elm model is stopping' call endrun(decomp_index=indexc, elmlevel=namec, msg=errmsg(__FILE__, __LINE__)) diff --git a/components/elm/src/biogeophys/SnowHydrologyMod.F90 b/components/elm/src/biogeophys/SnowHydrologyMod.F90 index f9270289d05d..49503da8ad8c 100644 --- a/components/elm/src/biogeophys/SnowHydrologyMod.F90 +++ b/components/elm/src/biogeophys/SnowHydrologyMod.F90 @@ -670,7 +670,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & if (bi > dm) ddz1 = ddz1*exp(-46.0e-3_r8*(bi-dm)) else ddz1_fresh = (-grav * (burden(c) + wx/2._r8)) / & - (0.007_r8 * bi**(4.75_r8 + td/40._r8)) + (0.007_r8 * min(max(bi,dm),denice)**(4.75_r8 + min(td,0._r8)/40._r8)) snw_ssa = 3.e6_r8 / (denice * snw_rds(c,j)) if (snw_ssa < 50._r8) then ddz1_fresh = ddz1_fresh * exp(-46.e-2_r8 * (50._r8 - snw_ssa)) From f4af4bd6906288473e25be8d5b4296bb53f87da2 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 29 Oct 2024 10:49:32 -0700 Subject: [PATCH 03/39] snow balance check error bugfix --- .../elm/src/biogeophys/BalanceCheckMod.F90 | 21 +++++++++++-------- .../src/biogeophys/HydrologyDrainageMod.F90 | 5 +++-- .../elm/src/biogeophys/SnowHydrologyMod.F90 | 4 +++- .../elm/src/biogeophys/SoilHydrologyMod.F90 | 10 ++++++++- .../elm/src/data_types/ColumnDataType.F90 | 2 ++ components/elm/src/main/lnd2atmMod.F90 | 4 +++- 6 files changed, 32 insertions(+), 14 deletions(-) diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index 0fd26bec1e64..5a1ae2441cec 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -164,7 +164,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & ! error = abs(precipitation - change of water storage - evaporation - runoff) ! ! !USES: - use elm_varcon , only : spval + use elm_varcon , only : spval, h2osno_max use column_varcon , only : icol_roof, icol_sunwall, icol_shadewall use column_varcon , only : icol_road_perv, icol_road_imperv use landunit_varcon , only : istice_mec, istdlak, istsoil,istcrop,istwet @@ -221,6 +221,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & qflx_surf_irrig_col => col_wf%qflx_surf_irrig , & ! Input: [real(r8) (:) ] real surface irrigation flux (mm H2O /s) qflx_over_supply_col => col_wf%qflx_over_supply , & ! Input: [real(r8) (:) ] over supply irrigation flux (mm H2O /s) qflx_snwcp_ice => col_wf%qflx_snwcp_ice , & ! Input: [real(r8) (:) ] excess snowfall due to snow capping (mm H2O /s) [+]` + qflx_ice_runoff_xs => col_wf%qflx_ice_runoff_xs , & ! Input: [real(r8) (:) ] qflx_evap_tot => col_wf%qflx_evap_tot , & ! Input: [real(r8) (:) ] qflx_evap_soi + qflx_evap_can + qflx_tran_veg qflx_dew_snow => col_wf%qflx_dew_snow , & ! Input: [real(r8) (:) ] surface dew added to snow pack (mm H2O /s) [+] qflx_sub_snow => col_wf%qflx_sub_snow , & ! Input: [real(r8) (:) ] sublimation rate from snow pack (mm H2O /s) [+] @@ -319,7 +320,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & errh2o(c) = endwb(c) - begwb(c) & - (forc_rain_col(c) + forc_snow_col(c) + qflx_floodc(c) + qflx_surf_irrig_col(c) + qflx_over_supply_col(c) & - qflx_evap_tot(c) - qflx_surf(c) - qflx_h2osfc_surf(c) & - - qflx_qrgwl(c) - qflx_drain(c) - qflx_drain_perched(c) - qflx_snwcp_ice(c) & + - qflx_qrgwl(c) - qflx_drain(c) - qflx_drain_perched(c) - qflx_snwcp_ice(c) - qflx_ice_runoff_xs(c) & - qflx_lateral(c) + qflx_h2orof_drain(c)) * dtime dwb(c) = (endwb(c)-begwb(c))/dtime @@ -391,6 +392,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & write(iulog,*)'qflx_lateral = ',qflx_lateral(indexc) write(iulog,*)'total_plant_stored_h2o_col = ',total_plant_stored_h2o_col(indexc) write(iulog,*)'qflx_h2orof_drain = ',qflx_h2orof_drain(indexc) + write(iulog,*)'qflx_ice_runoff_xs = ',qflx_ice_runoff_xs(indexc) write(iulog,*)'elm model is stopping' call endrun(decomp_index=indexc, elmlevel=namec, msg=errmsg(__FILE__, __LINE__)) @@ -421,6 +423,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & write(iulog,*)'qflx_lateral = ',qflx_lateral(indexc) write(iulog,*)'total_plant_stored_h2o_col = ',total_plant_stored_h2o_col(indexc) write(iulog,*)'qflx_h2orof_drain = ',qflx_h2orof_drain(indexc) + write(iulog,*)'qflx_ice_runoff_xs = ',qflx_ice_runoff_xs(indexc) write(iulog,*)'elm model is stopping' call endrun(decomp_index=indexc, elmlevel=namec, msg=errmsg(__FILE__, __LINE__)) end if @@ -490,16 +493,15 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & + qflx_snow_melt(c) + qflx_sl_top_soil(c) endif else ! firn model - snow_sources(c) = (qflx_snow_grnd_col(c) - qflx_snow_h2osfc(c) ) & - + frac_sno_eff(c) * (qflx_rain_grnd_col(c) & - + qflx_dew_snow(c) + qflx_dew_grnd(c) ) + qflx_h2osfc_to_ice(c) + snow_sources(c) = (qflx_snow_grnd_col(c) - qflx_snow_h2osfc(c) ) & + + frac_sno_eff(c) * (qflx_rain_grnd_col(c) & + + qflx_dew_snow(c) + qflx_dew_grnd(c) ) + qflx_h2osfc_to_ice(c) - snow_sinks(c) = frac_sno_eff(c) * (qflx_sub_snow(c) + qflx_evap_grnd(c)) & - + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & - + qflx_snow_melt(c) + qflx_sl_top_soil(c) + snow_sinks(c) = frac_sno_eff(c) * (qflx_sub_snow(c) + qflx_evap_grnd(c)) & + + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & + + qflx_snow_melt(c) + qflx_sl_top_soil(c) endif endif - if (glc_dyn_runoff_routing(g)) then ! Need to add qflx_glcice_frz to snow_sinks for the same reason as it is ! added to errh2o above - see the comment above for details. @@ -537,6 +539,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & write(iulog,*)'nstep = ',nstep g = col_pp%gridcell(indexc) l = col_pp%landunit(indexc) + write(iulog,*)'h2osno_max = ',h2osno_max write(iulog,*)'errh2osno = ',errh2osno(indexc) write(iulog,*)'snl = ',col_pp%snl(indexc) write(iulog,*)'h2osno = ',h2osno(indexc) diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index 850b7eedac52..da4006e9771c 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -7,7 +7,7 @@ module HydrologyDrainageMod use shr_kind_mod , only : r8 => shr_kind_r8 use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type - use elm_varctl , only : iulog, use_vichydro + use elm_varctl , only : iulog, use_vichydro, use_firn_percolation_and_compaction use elm_varcon , only : e_ice, denh2o, denice, rpi, spval use atm2lndType , only : atm2lnd_type use glc2lndMod , only : glc2lnd_type @@ -259,8 +259,9 @@ subroutine HydrologyDrainage(bounds, & ! glc_dyn_runoff_routing = true: in this case, melting ice runs off, and excess ! snow is sent to CISM, where it is converted to ice. These corrections are ! done here: - if (glc_dyn_runoff_routing(g) .and. lun_pp%itype(l)==istice_mec) then + ! this allows GLC melt to runoff to qflx_qrgwl! + ! If glc_dyn_runoff_routing=T, add meltwater from istice_mec ice columns to the runoff. ! Note: The meltwater contribution is computed in PhaseChanges (part of Biogeophysics2) qflx_qrgwl(c) = qflx_qrgwl(c) + qflx_glcice_melt(c) diff --git a/components/elm/src/biogeophys/SnowHydrologyMod.F90 b/components/elm/src/biogeophys/SnowHydrologyMod.F90 index 49503da8ad8c..74b8b20b92f4 100644 --- a/components/elm/src/biogeophys/SnowHydrologyMod.F90 +++ b/components/elm/src/biogeophys/SnowHydrologyMod.F90 @@ -2273,6 +2273,7 @@ subroutine SnowCapping(bounds, num_nolakec, filter_initc, num_snowc, filter_snow associate( & qflx_snwcp_ice => col_wf%qflx_snwcp_ice , & ! Output: [real(r8) (:) ] excess solid h2o due to snow capping (outgoing) (mm H2O /s) [+] qflx_snwcp_liq => col_wf%qflx_snwcp_liq , & ! Output: [real(r8) (:) ] excess liquid h2o due to snow capping (outgoing) (mm H2O /s) [+] + qflx_ice_runoff_xs => col_wf%qflx_ice_runoff_xs , & ! Input: [real(r8) (:) ] h2osoi_ice => col_ws%h2osoi_ice , & ! In/Out: [real(r8) (:,:) ] ice lens (kg/m2) h2osoi_liq => col_ws%h2osoi_liq , & ! In/Out: [real(r8) (:,:) ] liquid water (kg/m2) h2osno => col_ws%h2osno , & ! Input: [real(r8) (:) ] snow water (mm H2O) @@ -2294,6 +2295,7 @@ subroutine SnowCapping(bounds, num_nolakec, filter_initc, num_snowc, filter_snow c = filter_initc(fc) qflx_snwcp_ice(c) = 0.0_r8 qflx_snwcp_liq(c) = 0.0_r8 + qflx_ice_runoff_xs(c) = 0.0_r8 end do loop_columns: do fc = 1, num_snowc @@ -2342,7 +2344,7 @@ subroutine SnowCapping(bounds, num_nolakec, filter_initc, num_snowc, filter_snow mss_dst2(c,0) = mss_dst2(c,0) * frac_adjust mss_dst3(c,0) = mss_dst3(c,0) * frac_adjust mss_dst4(c,0) = mss_dst4(c,0) * frac_adjust - end if + end if end do loop_columns diff --git a/components/elm/src/biogeophys/SoilHydrologyMod.F90 b/components/elm/src/biogeophys/SoilHydrologyMod.F90 index 959fb2b69687..b6de90ed140a 100644 --- a/components/elm/src/biogeophys/SoilHydrologyMod.F90 +++ b/components/elm/src/biogeophys/SoilHydrologyMod.F90 @@ -987,6 +987,8 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte use elm_varctl , only : use_vsfm, use_var_soil_thick use SoilWaterMovementMod, only : zengdecker_2009_with_var_soil_thick use pftvarcon , only : rsub_top_globalmax + use LandunitType , only : lun_pp + use landunit_varcon , only : istice_mec, istice ! ! !ARGUMENTS: type(bounds_type) , intent(in) :: bounds @@ -1089,6 +1091,7 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte qflx_snwcp_liq => col_wf%qflx_snwcp_liq , & ! Output: [real(r8) (:) ] excess rainfall due to snow capping (mm H2O /s) [+] qflx_snwcp_ice => col_wf%qflx_snwcp_ice , & ! Output: [real(r8) (:) ] excess snowfall due to snow capping (mm H2O /s) [+] + qflx_ice_runoff_xs => col_wf%qflx_ice_runoff_xs , & ! Output: solid runoff from excess ice in soil (mm H2O /s) [+] !qflx_dew_grnd => col_wf%qflx_dew_grnd , & ! Output: [real(r8) (:) ] ground surface dew formation (mm H2O /s) [+] !qflx_dew_snow => col_wf%qflx_dew_snow , & ! Output: [real(r8) (:) ] surface dew added to snow pack (mm H2O /s) [+] !qflx_sub_snow => col_wf%qflx_sub_snow , & ! Output: [real(r8) (:) ] sublimation rate from snow pack (mm H2O /s) [+] @@ -1125,6 +1128,7 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte rsub_top(c) = 0._r8 fracice_rsub(c) = 0._r8 qflx_qrgwl(c) = 0._r8 + !qflx_ice_runoff_xs(c) = 0._r8 end do ! The layer index of the first unsaturated layer, i.e., the layer right above @@ -1527,7 +1531,11 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte ! add in ice check xs1(c) = max(max(h2osoi_ice(c,1),0._r8)-max(0._r8,(pondmx+watsat(c,1)*dzmm(c,1)-h2osoi_liq(c,1))),0._r8) h2osoi_ice(c,1) = min(max(0._r8,pondmx+watsat(c,1)*dzmm(c,1)-h2osoi_liq(c,1)), h2osoi_ice(c,1)) - qflx_snwcp_ice(c) = qflx_snwcp_ice(c) + xs1(c) / dtime + if (lun_pp%itype(col_pp%landunit(c)) == istice .or. lun_pp%itype(col_pp%landunit(c)) == istice_mec) then + qflx_snwcp_ice(c) = qflx_snwcp_ice(c) + xs1(c) / dtime + else + qflx_ice_runoff_xs(c) = qflx_ice_runoff_xs(c) + xs1(c) / dtime + endif end do ! Limit h2osoi_liq to be greater than or equal to watmin. diff --git a/components/elm/src/data_types/ColumnDataType.F90 b/components/elm/src/data_types/ColumnDataType.F90 index d0a4a10cd3d6..bafc3028eb9a 100644 --- a/components/elm/src/data_types/ColumnDataType.F90 +++ b/components/elm/src/data_types/ColumnDataType.F90 @@ -464,6 +464,7 @@ module ColumnDataType real(r8), pointer :: qflx_evap_grnd (:) => null() ! ground surface evaporation rate (mm H2O/s) [+] real(r8), pointer :: qflx_snwcp_liq (:) => null() ! excess rainfall due to snow capping (mm H2O /s) real(r8), pointer :: qflx_snwcp_ice (:) => null() ! excess snowfall due to snow capping (mm H2O /s) + real(r8), pointer :: qflx_ice_runoff_xs (:) => null() ! solid runoff from excess ice in soil (mm H2O /s) [+] real(r8), pointer :: qflx_tran_veg (:) => null() ! vegetation transpiration (mm H2O/s) (+ = to atm) real(r8), pointer :: qflx_dew_snow (:) => null() ! surface dew added to snow pack (mm H2O /s) [+] real(r8), pointer :: qflx_dew_grnd (:) => null() ! ground surface dew formation (mm H2O /s) [+] (+ = to atm); usually eflx_bot >= 0) @@ -5687,6 +5688,7 @@ subroutine col_wf_init(this, begc, endc) allocate(this%qflx_evap_grnd (begc:endc)) ; this%qflx_evap_grnd (:) = spval allocate(this%qflx_snwcp_liq (begc:endc)) ; this%qflx_snwcp_liq (:) = spval allocate(this%qflx_snwcp_ice (begc:endc)) ; this%qflx_snwcp_ice (:) = spval + allocate(this%qflx_ice_runoff_xs (begc:endc)) ; this%qflx_ice_runoff_xs (:) = 0._r8 allocate(this%qflx_tran_veg (begc:endc)) ; this%qflx_tran_veg (:) = spval allocate(this%qflx_dew_snow (begc:endc)) ; this%qflx_dew_snow (:) = spval allocate(this%qflx_dew_grnd (begc:endc)) ; this%qflx_dew_grnd (:) = spval diff --git a/components/elm/src/main/lnd2atmMod.F90 b/components/elm/src/main/lnd2atmMod.F90 index 2ebb57038051..978be5a1b6ff 100644 --- a/components/elm/src/main/lnd2atmMod.F90 +++ b/components/elm/src/main/lnd2atmMod.F90 @@ -224,6 +224,7 @@ subroutine lnd2atm(bounds, & qflx_qrgwl => col_wf%qflx_qrgwl , & qflx_rofliq_qgwl_grc => lnd2atm_vars%qflx_rofliq_qgwl_grc , & qflx_snwcp_ice => col_wf%qflx_snwcp_ice, & + qflx_ice_runoff_xs => col_wf%qflx_ice_runoff_xs, & qflx_rofice_grc => lnd2atm_vars%qflx_rofice_grc , & endwb => col_ws%endwb , & tws => grc_ws%tws , & @@ -418,7 +419,8 @@ subroutine lnd2atm(bounds, & qflx_qrgwl (bounds%begc:bounds%endc) , & qflx_rofliq_qgwl_grc(bounds%begg:bounds%endg) , & c2l_scale_type= urbanf, l2g_scale_type=unity ) - + + qflx_snwcp_ice(bounds%begc:bounds%endc) = qflx_snwcp_ice(bounds%begc:bounds%endc) + qflx_ice_runoff_xs(bounds%begc:bounds%endc) call c2g( bounds, & qflx_snwcp_ice (bounds%begc:bounds%endc) , & qflx_rofice_grc(bounds%begg:bounds%endg) , & From 629cea64492b73b1bcf3bd5e9e9e25857adb5bb1 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 29 Oct 2024 13:02:23 -0700 Subject: [PATCH 04/39] adding snow diagnostic fields --- .../elm/src/biogeophys/TotalWaterAndHeatMod.F90 | 3 +++ components/elm/src/data_types/ColumnDataType.F90 | 14 ++++++++++++++ 2 files changed, 17 insertions(+) diff --git a/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 b/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 index f222dca037e2..a43683aeadd8 100644 --- a/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 +++ b/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 @@ -418,6 +418,7 @@ subroutine ComputeHeatNonLake(bounds, num_nolakec, filter_nolakec, & real(r8) :: heat_dry_mass(bounds%begc:bounds%endc) ! sum of heat content: dry mass [J/m^2] real(r8) :: heat_ice(bounds%begc:bounds%endc) ! sum of heat content: ice [J/m^2] real(r8) :: latent_heat_liquid(bounds%begc:bounds%endc) ! sum of latent heat content of liquid water [J/m^2] + real(r8) :: sno_latent_heat_lcl(bounds%begc:bounds%endc,-nlevsno+1:0) ! sum of latent heat content of liquid water in snow [J/m^2] !character(len=*), parameter :: subname = 'ComputeHeatNonLake' !----------------------------------------------------------------------- @@ -438,6 +439,7 @@ subroutine ComputeHeatNonLake(bounds, num_nolakec, filter_nolakec, & h2osno => col_ws%h2osno, & ! snow water (mm H2O) h2osfc => col_ws%h2osfc, & ! surface water (mm H2O) h2ocan_patch => veg_ws%h2ocan, & ! canopy water (mm H2O) + sno_latent_heat_lcl => col_ws%sno_latent_heat, & ! snocan_patch => waterstate_inst%snocan_patch, & ! canopy snow water (mm H2O) total_plant_stored_h2o_col => col_ws%total_plant_stored_h2o, & ! Input: [real(r8) (:) ] water mass in plant tissues (kg m-2) wa => soilhydrology_inst%wa_col & ! water in the unconfined aquifer (mm) @@ -514,6 +516,7 @@ subroutine ComputeHeatNonLake(bounds, num_nolakec, filter_nolakec, & latent_heat_liquid = latent_heat_liquid(c)) heat_ice(c) = heat_ice(c) + & TempToHeat(t_soisno(c,j), (h2osoi_ice(c,j)*cpice)) + sno_latent_heat_lcl(c,j) = latent_heat_liquid(c) end do else if (h2osno(c) /= 0._r8) then ! No explicit snow layers, but there may still be some ice in h2osno (there is diff --git a/components/elm/src/data_types/ColumnDataType.F90 b/components/elm/src/data_types/ColumnDataType.F90 index bafc3028eb9a..7bf68e2630f6 100644 --- a/components/elm/src/data_types/ColumnDataType.F90 +++ b/components/elm/src/data_types/ColumnDataType.F90 @@ -124,6 +124,7 @@ module ColumnDataType real(r8), pointer :: swe_old (:,:) => null() ! initial snow water content (-nlevsno+1:0) (kg/m2) real(r8), pointer :: snw_rds (:,:) => null() ! col snow grain radius (-nlevsno+1:0) (m^-6, or microns) real(r8), pointer :: air_vol (:,:) => null() ! air filled porosity (m3/m3) + real(r8), pointer :: sno_latent_heat (:,:) => null() ! latent heat in snow [J/m^2] ! Derived water, ice, and snow variables, column aggregate real(r8), pointer :: qg_snow (:) => null() ! specific humidity over snow (kg H2O/kg moist air) real(r8), pointer :: qg_soil (:) => null() ! specific humidity over soil (kg H2O/kg moist air) @@ -1395,6 +1396,7 @@ subroutine col_ws_init(this, begc, endc, h2osno_input, snow_depth_input, watsat_ allocate(this%h2osoi_ice_old (begc:endc,-nlevsno+1:nlevgrnd)) ; this%h2osoi_ice_old (:,:) = spval allocate(this%bw (begc:endc,-nlevsno+1:0)) ; this%bw (:,:) = spval allocate(this%smp_l (begc:endc,-nlevsno+1:nlevgrnd)) ; this%smp_l (:,:) = spval + allocate(this%sno_latent_heat (begc:endc,-nlevsno+1:0)) ; this%sno_latent_heat (:,:) = spval allocate(this%soilp (begc:endc,1:nlevgrnd)) ; this%soilp (:,:) = 0._r8 allocate(this%swe_old (begc:endc,-nlevsno+1:0)) ; this%swe_old (:,:) = spval allocate(this%snw_rds (begc:endc,-nlevsno+1:0)) ; this%snw_rds (:,:) = spval @@ -1493,6 +1495,18 @@ subroutine col_ws_init(this, begc, endc, h2osno_input, snow_depth_input, watsat_ avgflag='A', long_name='Partial density of water in the snow pack (ice + liquid)', & ptr_col=data2dptr, no_snow_behavior=no_snow_normal, default='inactive') + this%qflx_snofrz_lyr(begc:endc,-nlevsno+1:0) = spval + data2dptr => this%qflx_snofrz_lyr(begc:endc,-nlevsno+1:0) + call hist_addfld2d (fname='QSNOFRZ_LYR', units='kg/m2/s', type2d='levsno',& + avgflag='I', long_name='layer snow freezing rate', & + ptr_col=data2dptr, no_snow_behavior=no_snow_normal, default='inactive') + + this%sno_latent_heat(begc:endc,-nlevsno+1:0) = spval + data2dptr => this%sno_latent_heat(:,-nlevsno+1:0) + call hist_addfld2d (fname='SNO_LTNT_HT', units='J/m^2', type2d='levsno', & + avgflag='I', long_name='latent heat in snow', & + ptr_col=data2dptr, no_snow_behavior=no_snow_normal, default='inactive') + this%snw_rds(begc:endc,-nlevsno+1:0) = spval data2dptr => this%snw_rds(:,-nlevsno+1:0) call hist_addfld2d (fname='SNO_GS', units='Microns', type2d='levsno', & From ac89b00b76675c6ab2642aceab4f0a1270be9229 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 29 Oct 2024 17:08:33 -0700 Subject: [PATCH 05/39] adding snow diag fields --- .../elm/src/biogeophys/TotalWaterAndHeatMod.F90 | 2 +- components/elm/src/data_types/ColumnDataType.F90 | 11 +++++------ 2 files changed, 6 insertions(+), 7 deletions(-) diff --git a/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 b/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 index a43683aeadd8..620626e2ddab 100644 --- a/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 +++ b/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 @@ -10,7 +10,7 @@ module TotalWaterAndHeatMod use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type use elm_varcon , only : cpice, cpliq, denh2o, tfrz, hfus, aquifer_water_baseline - use elm_varpar , only : nlevgrnd, nlevsoi, nlevurb, nlevlak + use elm_varpar , only : nlevgrnd, nlevsoi, nlevurb, nlevlak, nlevsno use subgridAveMod , only : p2c use SoilHydrologyType , only : soilhydrology_type use UrbanParamsType , only : urbanparams_type diff --git a/components/elm/src/data_types/ColumnDataType.F90 b/components/elm/src/data_types/ColumnDataType.F90 index 7bf68e2630f6..a3dbed405659 100644 --- a/components/elm/src/data_types/ColumnDataType.F90 +++ b/components/elm/src/data_types/ColumnDataType.F90 @@ -1495,12 +1495,6 @@ subroutine col_ws_init(this, begc, endc, h2osno_input, snow_depth_input, watsat_ avgflag='A', long_name='Partial density of water in the snow pack (ice + liquid)', & ptr_col=data2dptr, no_snow_behavior=no_snow_normal, default='inactive') - this%qflx_snofrz_lyr(begc:endc,-nlevsno+1:0) = spval - data2dptr => this%qflx_snofrz_lyr(begc:endc,-nlevsno+1:0) - call hist_addfld2d (fname='QSNOFRZ_LYR', units='kg/m2/s', type2d='levsno',& - avgflag='I', long_name='layer snow freezing rate', & - ptr_col=data2dptr, no_snow_behavior=no_snow_normal, default='inactive') - this%sno_latent_heat(begc:endc,-nlevsno+1:0) = spval data2dptr => this%sno_latent_heat(:,-nlevsno+1:0) call hist_addfld2d (fname='SNO_LTNT_HT', units='J/m^2', type2d='levsno', & @@ -5778,6 +5772,11 @@ subroutine col_wf_init(this, begc, endc) !----------------------------------------------------------------------- ! initialize history fields for select members of col_wf !----------------------------------------------------------------------- + this%qflx_snofrz_lyr(begc:endc,-nlevsno+1:0) = spval + call hist_addfld2d (fname='QSNOFRZ_LYR', units='kg/m2/s', type2d='levsno',& + avgflag='I', long_name='layer snow freezing rate', & + ptr_col=this%qflx_snofrz_lyr, no_snow_behavior=no_snow_normal, default='inactive') + this%qflx_snwcp_ice(begc:endc) = spval call hist_addfld1d (fname='QSNWCPICE_NODYNLNDUSE', units='mm H2O/s', & avgflag='A', long_name='excess snowfall due to snow capping not including correction for land use change', & From a0af6eb387d22bd4d62ce78ab6e7f91b76a954bc Mon Sep 17 00:00:00 2001 From: Chloe Date: Wed, 13 Nov 2024 05:43:38 -0800 Subject: [PATCH 06/39] adding snow diagnostic fields --- components/elm/src/biogeophys/LakeTemperatureMod.F90 | 11 +++++++---- components/elm/src/biogeophys/SoilTemperatureMod.F90 | 4 ++++ components/elm/src/data_types/ColumnDataType.F90 | 9 ++++++++- 3 files changed, 19 insertions(+), 5 deletions(-) diff --git a/components/elm/src/biogeophys/LakeTemperatureMod.F90 b/components/elm/src/biogeophys/LakeTemperatureMod.F90 index 13ed09d1df71..f17a4a32a18e 100644 --- a/components/elm/src/biogeophys/LakeTemperatureMod.F90 +++ b/components/elm/src/biogeophys/LakeTemperatureMod.F90 @@ -1291,6 +1291,7 @@ subroutine PhaseChange_Lake (bounds, num_lakec, filter_lakec, cv, cv_lake, lhabs qflx_snofrz_lyr => col_wf%qflx_snofrz_lyr , & ! Input: [real(r8) (:,:) ] snow freezing rate (positive definite) (col,lyr) [kg m-2 s-1] qflx_snow_melt => col_wf%qflx_snow_melt , & ! Output: [real(r8) (:) ] net snow melt qflx_snomelt => col_wf%qflx_snomelt , & ! Output: [real(r8) (:) ] snow melt (mm H2O /s) + qflx_snomelt_lyr=> col_wf%qflx_snomelt_lyr, & ! Output: [real(r8) (:) ] snow melt per layer (mm H2O /s) qflx_snofrz_col => col_wf%qflx_snofrz , & ! Output: [real(r8) (:) ] column-integrated snow freezing rate (kg m-2 s-1) [+] t_soisno => col_es%t_soisno , & ! Input: [real(r8) (:,:) ] soil temperature (Kelvin) @@ -1306,10 +1307,11 @@ subroutine PhaseChange_Lake (bounds, num_lakec, filter_lakec, cv, cv_lake, lhabs do fc = 1,num_lakec c = filter_lakec(fc) - qflx_snomelt(c) = 0._r8 - eflx_snomelt(c) = 0._r8 - lhabs(c) = 0._r8 - qflx_snow_melt(c) = 0._r8 + qflx_snomelt(c) = 0._r8 + qflx_snomelt_lyr(c,:)= 0._r8 + eflx_snomelt(c) = 0._r8 + lhabs(c) = 0._r8 + qflx_snow_melt(c) = 0._r8 end do do j = -nlevsno+1,0 @@ -1399,6 +1401,7 @@ subroutine PhaseChange_Lake (bounds, num_lakec, filter_lakec, cv, cv_lake, lhabs if (j <= 0) then !snow imelt(c,j) = 1 qflx_snomelt(c) = qflx_snomelt(c) + melt/dtime + qflx_snomelt_lyr(c,j) = melt/dtime end if else if (t_soisno(c,j) < tfrz .and. h2osoi_liq(c,j) > 0._r8) then !freezing dophasechangeflag = .true. diff --git a/components/elm/src/biogeophys/SoilTemperatureMod.F90 b/components/elm/src/biogeophys/SoilTemperatureMod.F90 index d4a1074bf4a7..7420e857bce7 100644 --- a/components/elm/src/biogeophys/SoilTemperatureMod.F90 +++ b/components/elm/src/biogeophys/SoilTemperatureMod.F90 @@ -1370,6 +1370,7 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & qflx_glcice => col_wf%qflx_glcice , & ! Output: [real(r8) (:) ] flux of new glacier ice (mm H2O/s) [+ = ice grows] qflx_glcice_melt => col_wf%qflx_glcice_melt , & ! Output: [real(r8) (:) ] ice melt (positive definite) (mm H2O/s) qflx_snomelt => col_wf%qflx_snomelt , & ! Output: [real(r8) (:) ] snow melt (mm H2O /s) + qflx_snomelt_lyr => col_wf%qflx_snomelt_lyr , & ! Output: [real(r8) (:) ] snow melt (mm H2O /s) eflx_snomelt => col_ef%eflx_snomelt , & ! Output: [real(r8) (:) ] snow melt heat flux (W/m**2) eflx_snomelt_r => col_ef%eflx_snomelt_r , & ! Output: [real(r8) (:) ] rural snow melt heat flux (W/m**2) @@ -1389,6 +1390,7 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & l = col_pp%landunit(c) qflx_snomelt(c) = 0._r8 + qflx_snomelt_lyr(c,-nlevsno+1:0) = 0._r8 xmf(c) = 0._r8 qflx_snofrz_lyr(c,-nlevsno+1:0) = 0._r8 qflx_snofrz_col(c) = 0._r8 @@ -1565,6 +1567,7 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & qflx_snomelt(c) = max(0._r8,(temp1-h2osno(c)))/dtime ! kg/(m2 s) xmf(c) = hfus*qflx_snomelt(c) qflx_snow_melt(c) = qflx_snomelt(c) + !qflx_snomelt_lyr(c,j) = qflx_snomelt(c) endif endif @@ -1623,6 +1626,7 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & if (imelt(c,j) == 1 .AND. j < 1) then qflx_snomelt(c) = qflx_snomelt(c) + max(0._r8,(wice0(c,j)-h2osoi_ice(c,j)))/dtime + qflx_snomelt_lyr(c,j) = qflx_snomelt(c) endif diff --git a/components/elm/src/data_types/ColumnDataType.F90 b/components/elm/src/data_types/ColumnDataType.F90 index a3dbed405659..fa66a5b77cfa 100644 --- a/components/elm/src/data_types/ColumnDataType.F90 +++ b/components/elm/src/data_types/ColumnDataType.F90 @@ -494,6 +494,7 @@ module ColumnDataType real(r8), pointer :: qflx_sl_top_soil (:) => null() ! liquid water + ice from layer above soil to top soil layer or sent to qflx_qrgwl (mm H2O/s) real(r8), pointer :: qflx_snomelt (:) => null() ! snow melt (mm H2O /s) real(r8), pointer :: qflx_snow_melt (:) => null() ! snow melt (net) + real(r8), pointer :: qflx_snomelt_lyr (:,:) => null() ! snow melt (net) real(r8), pointer :: qflx_qrgwl (:) => null() ! qflx_surf at glaciers, wetlands, lakes real(r8), pointer :: qflx_runoff (:) => null() ! total runoff (qflx_drain+qflx_surf+qflx_qrgwl) (mm H2O /s) real(r8), pointer :: qflx_runoff_r (:) => null() ! Rural total runoff (qflx_drain+qflx_surf+qflx_qrgwl) (mm H2O /s) @@ -5724,6 +5725,7 @@ subroutine col_wf_init(this, begc, endc) allocate(this%qflx_floodc (begc:endc)) ; this%qflx_floodc (:) = spval allocate(this%qflx_sl_top_soil (begc:endc)) ; this%qflx_sl_top_soil (:) = spval allocate(this%qflx_snomelt (begc:endc)) ; this%qflx_snomelt (:) = spval + allocate(this%qflx_snomelt_lyr (begc:endc,-nlevsno+1:0)) ; this%qflx_snomelt_lyr (:,:) = spval allocate(this%qflx_snow_melt (begc:endc)) ; this%qflx_snow_melt (:) = spval allocate(this%qflx_qrgwl (begc:endc)) ; this%qflx_qrgwl (:) = spval allocate(this%qflx_runoff (begc:endc)) ; this%qflx_runoff (:) = spval @@ -5827,6 +5829,11 @@ subroutine col_wf_init(this, begc, endc) avgflag='A', long_name='snow melt', & ptr_col=this%qflx_snow_melt, c2l_scale_type='urbanf') + this%qflx_snomelt_lyr(begc:endc,-nlevsno+1:0) = spval + call hist_addfld2d (fname='QSNOMELT_LYR', units='mm/s',type2d='levsno',& + avgflag='A', long_name='snow melt per snow layer', & + ptr_col=this%qflx_snomelt_lyr,no_snow_behavior=no_snow_normal, c2l_scale_type='urbanf') + this%qflx_qrgwl(begc:endc) = spval call hist_addfld1d (fname='QRGWL', units='mm/s', & avgflag='A', long_name='surface runoff at glaciers (liquid only), wetlands, lakes', & @@ -5856,7 +5863,7 @@ subroutine col_wf_init(this, begc, endc) this%qflx_snofrz(begc:endc) = spval call hist_addfld1d (fname='QSNOFRZ', units='kg/m2/s', & avgflag='A', long_name='column-integrated snow freezing rate', & - ptr_col=this%qflx_snofrz, set_lake=spval, c2l_scale_type='urbanf', default='inactive') + ptr_col=this%qflx_snofrz, set_lake=spval, c2l_scale_type='urbanf') if (create_glacier_mec_landunit) then this%qflx_glcice(begc:endc) = spval From 36f65ab3de041e809e84a59a6fee172082ce9cea Mon Sep 17 00:00:00 2001 From: Chloe Date: Mon, 25 Nov 2024 15:41:16 -0800 Subject: [PATCH 07/39] cleaning up changes for PR submission --- .../elm/src/biogeophys/BalanceCheckMod.F90 | 55 +++---------------- .../elm/src/biogeophys/SnowHydrologyMod.F90 | 2 +- .../elm/src/biogeophys/SoilHydrologyMod.F90 | 1 - .../src/biogeophys/TotalWaterAndHeatMod.F90 | 5 +- .../elm/src/data_types/ColumnDataType.F90 | 8 --- 5 files changed, 10 insertions(+), 61 deletions(-) diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index 5a1ae2441cec..ec8b02f90895 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -493,13 +493,13 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & + qflx_snow_melt(c) + qflx_sl_top_soil(c) endif else ! firn model - snow_sources(c) = (qflx_snow_grnd_col(c) - qflx_snow_h2osfc(c) ) & - + frac_sno_eff(c) * (qflx_rain_grnd_col(c) & - + qflx_dew_snow(c) + qflx_dew_grnd(c) ) + qflx_h2osfc_to_ice(c) + snow_sources(c) = (qflx_snow_grnd_col(c) - qflx_snow_h2osfc(c) ) & + + frac_sno_eff(c) * (qflx_rain_grnd_col(c) & + + qflx_dew_snow(c) + qflx_dew_grnd(c) ) + qflx_h2osfc_to_ice(c) - snow_sinks(c) = frac_sno_eff(c) * (qflx_sub_snow(c) + qflx_evap_grnd(c)) & - + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & - + qflx_snow_melt(c) + qflx_sl_top_soil(c) + snow_sinks(c) = frac_sno_eff(c) * (qflx_sub_snow(c) + qflx_evap_grnd(c)) & + + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & + + qflx_snow_melt(c) + qflx_sl_top_soil(c) endif endif if (glc_dyn_runoff_routing(g)) then @@ -536,42 +536,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & ' col_pp%itype= ',col_pp%itype(indexc), & ' lun_pp%itype= ',lun_pp%itype(col_pp%landunit(indexc)), & ' errh2osno= ',errh2osno(indexc) - write(iulog,*)'nstep = ',nstep - g = col_pp%gridcell(indexc) - l = col_pp%landunit(indexc) - write(iulog,*)'h2osno_max = ',h2osno_max - write(iulog,*)'errh2osno = ',errh2osno(indexc) - write(iulog,*)'snl = ',col_pp%snl(indexc) - write(iulog,*)'h2osno = ',h2osno(indexc) - write(iulog,*)'h2osno_old = ',h2osno_old(indexc) - write(iulog,*)'snow_sources = ',snow_sources(indexc) - write(iulog,*)'snow_sinks = ',snow_sinks(indexc) - write(iulog,*)'qflx_prec_grnd = ',qflx_prec_grnd(indexc) - write(iulog,*)'qflx_sub_snow = ',qflx_sub_snow(indexc) - write(iulog,*)'qflx_evap_grnd = ',qflx_evap_grnd(indexc) - write(iulog,*)'qflx_top_soil = ',qflx_top_soil(indexc) - write(iulog,*)'qflx_dew_snow = ',qflx_dew_snow(indexc) - write(iulog,*)'qflx_dew_grnd = ',qflx_dew_grnd(indexc) - write(iulog,*)'qflx_snwcp_ice = ',qflx_snwcp_ice(indexc) - write(iulog,*)'qflx_snwcp_liq = ',qflx_snwcp_liq(indexc) - write(iulog,*)'qflx_sl_top_soil = ',qflx_sl_top_soil(indexc) - write(iulog,*)'qflx_snow_melt = ',qflx_snow_melt(indexc) - write(iulog,*)'frac_sno_eff = ',frac_sno_eff(indexc) - write(iulog,*)'qflx_rain_grnd_col= ',qflx_rain_grnd_col(indexc) - write(iulog,*)'qflx_snow_grnd_col= ',qflx_snow_grnd_col(indexc) - write(iulog,*)'qflx_h2osfc_to_ice= ',qflx_h2osfc_to_ice(indexc) - write(iulog,*)'qflx_snow_h2osfc= ',qflx_snow_h2osfc(indexc)*dtime - write(iulog,*)'(glc_dyn_runoff_routing(g))',(glc_dyn_runoff_routing(g)) - write(iulog,*)'lun_pp%itype(l))',lun_pp%itype(l) - write(iulog,*)'lun_pp%itype(l) == istsoil)',lun_pp%itype(l) == istsoil - write(iulog,*)'lun_pp%itype(l) == istcrop)',lun_pp%itype(l) == istcrop - write(iulog,*)'lun_pp%itype(l) == istwet)',lun_pp%itype(l) == istwet - - if (create_glacier_mec_landunit) then - write(iulog,*)'qflx_glcice_frz = ',qflx_glcice_frz(indexc) - end if - - + if (abs(errh2osno(indexc)) > 1.e-4_r8 .and. (nstep > 2) ) then write(iulog,*)'elm model is stopping - error is greater than 1e-4 (mm)' write(iulog,*)'nstep = ',nstep @@ -597,13 +562,9 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & write(iulog,*)'qflx_h2osfc_to_ice= ',qflx_h2osfc_to_ice(indexc) write(iulog,*)'qflx_snow_h2osfc= ',qflx_snow_h2osfc(indexc)*dtime write(iulog,*)'(glc_dyn_runoff_routing(g))',(glc_dyn_runoff_routing(g)) - write(iulog,*)'lun_pp%itype(l))',lun_pp%itype(l) - write(iulog,*)'lun_pp%itype(l) == istsoil)',lun_pp%itype(l) == istsoil - write(iulog,*)'lun_pp%itype(l) == istcrop)',lun_pp%itype(l) == istcrop - write(iulog,*)'lun_pp%itype(l) == istwet)',lun_pp%itype(l) == istwet if (create_glacier_mec_landunit) then - write(iulog,*)'qflx_glcice_frz = ',qflx_glcice_frz(indexc) + write(iulog,*)'qflx_glcice_frz = ',qflx_glcice_frz(indexc)*dtime` end if write(iulog,*)'elm model is stopping' call endrun(decomp_index=indexc, elmlevel=namec, msg=errmsg(__FILE__, __LINE__)) diff --git a/components/elm/src/biogeophys/SnowHydrologyMod.F90 b/components/elm/src/biogeophys/SnowHydrologyMod.F90 index 74b8b20b92f4..57e486d0ce90 100644 --- a/components/elm/src/biogeophys/SnowHydrologyMod.F90 +++ b/components/elm/src/biogeophys/SnowHydrologyMod.F90 @@ -2344,7 +2344,7 @@ subroutine SnowCapping(bounds, num_nolakec, filter_initc, num_snowc, filter_snow mss_dst2(c,0) = mss_dst2(c,0) * frac_adjust mss_dst3(c,0) = mss_dst3(c,0) * frac_adjust mss_dst4(c,0) = mss_dst4(c,0) * frac_adjust - end if + end if end do loop_columns diff --git a/components/elm/src/biogeophys/SoilHydrologyMod.F90 b/components/elm/src/biogeophys/SoilHydrologyMod.F90 index b6de90ed140a..01a13bde3ba4 100644 --- a/components/elm/src/biogeophys/SoilHydrologyMod.F90 +++ b/components/elm/src/biogeophys/SoilHydrologyMod.F90 @@ -1128,7 +1128,6 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte rsub_top(c) = 0._r8 fracice_rsub(c) = 0._r8 qflx_qrgwl(c) = 0._r8 - !qflx_ice_runoff_xs(c) = 0._r8 end do ! The layer index of the first unsaturated layer, i.e., the layer right above diff --git a/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 b/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 index 620626e2ddab..f222dca037e2 100644 --- a/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 +++ b/components/elm/src/biogeophys/TotalWaterAndHeatMod.F90 @@ -10,7 +10,7 @@ module TotalWaterAndHeatMod use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type use elm_varcon , only : cpice, cpliq, denh2o, tfrz, hfus, aquifer_water_baseline - use elm_varpar , only : nlevgrnd, nlevsoi, nlevurb, nlevlak, nlevsno + use elm_varpar , only : nlevgrnd, nlevsoi, nlevurb, nlevlak use subgridAveMod , only : p2c use SoilHydrologyType , only : soilhydrology_type use UrbanParamsType , only : urbanparams_type @@ -418,7 +418,6 @@ subroutine ComputeHeatNonLake(bounds, num_nolakec, filter_nolakec, & real(r8) :: heat_dry_mass(bounds%begc:bounds%endc) ! sum of heat content: dry mass [J/m^2] real(r8) :: heat_ice(bounds%begc:bounds%endc) ! sum of heat content: ice [J/m^2] real(r8) :: latent_heat_liquid(bounds%begc:bounds%endc) ! sum of latent heat content of liquid water [J/m^2] - real(r8) :: sno_latent_heat_lcl(bounds%begc:bounds%endc,-nlevsno+1:0) ! sum of latent heat content of liquid water in snow [J/m^2] !character(len=*), parameter :: subname = 'ComputeHeatNonLake' !----------------------------------------------------------------------- @@ -439,7 +438,6 @@ subroutine ComputeHeatNonLake(bounds, num_nolakec, filter_nolakec, & h2osno => col_ws%h2osno, & ! snow water (mm H2O) h2osfc => col_ws%h2osfc, & ! surface water (mm H2O) h2ocan_patch => veg_ws%h2ocan, & ! canopy water (mm H2O) - sno_latent_heat_lcl => col_ws%sno_latent_heat, & ! snocan_patch => waterstate_inst%snocan_patch, & ! canopy snow water (mm H2O) total_plant_stored_h2o_col => col_ws%total_plant_stored_h2o, & ! Input: [real(r8) (:) ] water mass in plant tissues (kg m-2) wa => soilhydrology_inst%wa_col & ! water in the unconfined aquifer (mm) @@ -516,7 +514,6 @@ subroutine ComputeHeatNonLake(bounds, num_nolakec, filter_nolakec, & latent_heat_liquid = latent_heat_liquid(c)) heat_ice(c) = heat_ice(c) + & TempToHeat(t_soisno(c,j), (h2osoi_ice(c,j)*cpice)) - sno_latent_heat_lcl(c,j) = latent_heat_liquid(c) end do else if (h2osno(c) /= 0._r8) then ! No explicit snow layers, but there may still be some ice in h2osno (there is diff --git a/components/elm/src/data_types/ColumnDataType.F90 b/components/elm/src/data_types/ColumnDataType.F90 index fa66a5b77cfa..0329430fbbd2 100644 --- a/components/elm/src/data_types/ColumnDataType.F90 +++ b/components/elm/src/data_types/ColumnDataType.F90 @@ -124,7 +124,6 @@ module ColumnDataType real(r8), pointer :: swe_old (:,:) => null() ! initial snow water content (-nlevsno+1:0) (kg/m2) real(r8), pointer :: snw_rds (:,:) => null() ! col snow grain radius (-nlevsno+1:0) (m^-6, or microns) real(r8), pointer :: air_vol (:,:) => null() ! air filled porosity (m3/m3) - real(r8), pointer :: sno_latent_heat (:,:) => null() ! latent heat in snow [J/m^2] ! Derived water, ice, and snow variables, column aggregate real(r8), pointer :: qg_snow (:) => null() ! specific humidity over snow (kg H2O/kg moist air) real(r8), pointer :: qg_soil (:) => null() ! specific humidity over soil (kg H2O/kg moist air) @@ -1397,7 +1396,6 @@ subroutine col_ws_init(this, begc, endc, h2osno_input, snow_depth_input, watsat_ allocate(this%h2osoi_ice_old (begc:endc,-nlevsno+1:nlevgrnd)) ; this%h2osoi_ice_old (:,:) = spval allocate(this%bw (begc:endc,-nlevsno+1:0)) ; this%bw (:,:) = spval allocate(this%smp_l (begc:endc,-nlevsno+1:nlevgrnd)) ; this%smp_l (:,:) = spval - allocate(this%sno_latent_heat (begc:endc,-nlevsno+1:0)) ; this%sno_latent_heat (:,:) = spval allocate(this%soilp (begc:endc,1:nlevgrnd)) ; this%soilp (:,:) = 0._r8 allocate(this%swe_old (begc:endc,-nlevsno+1:0)) ; this%swe_old (:,:) = spval allocate(this%snw_rds (begc:endc,-nlevsno+1:0)) ; this%snw_rds (:,:) = spval @@ -1496,12 +1494,6 @@ subroutine col_ws_init(this, begc, endc, h2osno_input, snow_depth_input, watsat_ avgflag='A', long_name='Partial density of water in the snow pack (ice + liquid)', & ptr_col=data2dptr, no_snow_behavior=no_snow_normal, default='inactive') - this%sno_latent_heat(begc:endc,-nlevsno+1:0) = spval - data2dptr => this%sno_latent_heat(:,-nlevsno+1:0) - call hist_addfld2d (fname='SNO_LTNT_HT', units='J/m^2', type2d='levsno', & - avgflag='I', long_name='latent heat in snow', & - ptr_col=data2dptr, no_snow_behavior=no_snow_normal, default='inactive') - this%snw_rds(begc:endc,-nlevsno+1:0) = spval data2dptr => this%snw_rds(:,-nlevsno+1:0) call hist_addfld2d (fname='SNO_GS', units='Microns', type2d='levsno', & From 194444d944b36da0c65828e90a725b1e7e43f2e4 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 3 Dec 2024 15:05:57 -0800 Subject: [PATCH 08/39] syntax error -minor change --- components/elm/src/biogeophys/BalanceCheckMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index ec8b02f90895..3ce501475dc5 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -564,7 +564,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & write(iulog,*)'(glc_dyn_runoff_routing(g))',(glc_dyn_runoff_routing(g)) if (create_glacier_mec_landunit) then - write(iulog,*)'qflx_glcice_frz = ',qflx_glcice_frz(indexc)*dtime` + write(iulog,*)'qflx_glcice_frz = ',qflx_glcice_frz(indexc)*dtime end if write(iulog,*)'elm model is stopping' call endrun(decomp_index=indexc, elmlevel=namec, msg=errmsg(__FILE__, __LINE__)) From c230f3a03fa942cc280f8fc832280e63a85db9d9 Mon Sep 17 00:00:00 2001 From: Chloe Date: Mon, 16 Dec 2024 16:40:48 -0800 Subject: [PATCH 09/39] remove use_firn_percolation_and_compaction for use_extrasnowlayers --- .../bld/namelist_files/namelist_definition.xml | 7 ------- components/elm/src/biogeophys/AerosolMod.F90 | 4 ++-- .../elm/src/biogeophys/BalanceCheckMod.F90 | 6 +++--- .../elm/src/biogeophys/CanopyHydrologyMod.F90 | 12 ++++++------ .../src/biogeophys/HydrologyDrainageMod.F90 | 2 +- .../src/biogeophys/HydrologyNoDrainageMod.F90 | 4 ++-- .../elm/src/biogeophys/LakeFluxesMod.F90 | 4 ++-- .../elm/src/biogeophys/LakeHydrologyMod.F90 | 18 +++++++++--------- .../elm/src/biogeophys/SnowHydrologyMod.F90 | 16 ++++++++-------- .../elm/src/biogeophys/SnowSnicarMod.F90 | 6 +++--- .../elm/src/biogeophys/SoilFluxesMod.F90 | 4 ++-- components/elm/src/main/controlMod.F90 | 4 +--- components/elm/src/main/elm_driver.F90 | 6 +++--- components/elm/src/main/elm_varctl.F90 | 1 - 14 files changed, 42 insertions(+), 52 deletions(-) diff --git a/components/elm/bld/namelist_files/namelist_definition.xml b/components/elm/bld/namelist_files/namelist_definition.xml index c76271b7baf7..f6c45e05859b 100644 --- a/components/elm/bld/namelist_files/namelist_definition.xml +++ b/components/elm/bld/namelist_files/namelist_definition.xml @@ -766,13 +766,6 @@ snowpack conditions on glaciers and ice sheets, including a semi-empirical firn densification model. - -This option enables more realistic -snowpack conditions on glaciers and ice sheets, including a semi-empirical -firn densification model. uses the same physics as use_extrasnowlayers but -uses the original 5 snow layer scheme - Toggle to use new snow thermal conductivity that relies on snow temperature and density. diff --git a/components/elm/src/biogeophys/AerosolMod.F90 b/components/elm/src/biogeophys/AerosolMod.F90 index 81fc8d62abc5..d97c48dfd29b 100644 --- a/components/elm/src/biogeophys/AerosolMod.F90 +++ b/components/elm/src/biogeophys/AerosolMod.F90 @@ -5,7 +5,7 @@ module AerosolMod use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type use elm_varpar , only : nlevsno - use elm_varctl , only : use_firn_percolation_and_compaction + use elm_varctl , only : use_extrasnowlayers use elm_time_manager , only : get_step_size use atm2lndType , only : atm2lnd_type use AerosolType , only : aerosol_type @@ -107,7 +107,7 @@ subroutine AerosolMasses(bounds, num_on, filter_on, num_off, filter_off, aerosol ! layer mass of snow: snowmass = h2osoi_ice(c,j) + h2osoi_liq(c,j) - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then ! Correct the top layer aerosol mass to account for snow capping. ! This approach conserves the aerosol mass concentration ! (but not the aerosol amss) when snow-capping is invoked diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index 3ce501475dc5..ae26ef0a96a2 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -9,7 +9,7 @@ module BalanceCheckMod use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type use abortutils , only : endrun - use elm_varctl , only : iulog, use_var_soil_thick, use_firn_percolation_and_compaction + use elm_varctl , only : iulog, use_var_soil_thick, use_extrasnowlayers use elm_varcon , only : namep, namec use GetGlobalValuesMod , only : GetGlobalIndex use atm2lndType , only : atm2lnd_type @@ -448,7 +448,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) + qflx_sl_top_soil(c) if (lun_pp%itype(l) == istdlak) then - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then if ( do_capsnow(c)) then snow_sources(c) = qflx_snow_grnd_col(c) & + frac_sno_eff(c) * (qflx_dew_snow(c) + qflx_dew_grnd(c) ) @@ -476,7 +476,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & endif if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop .or. lun_pp%itype(l) == istwet ) then - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then if (do_capsnow(c)) then snow_sources(c) = frac_sno_eff(c) * (qflx_dew_snow(c) + qflx_dew_grnd(c) ) & + qflx_h2osfc_to_ice(c) + qflx_prec_grnd(c) diff --git a/components/elm/src/biogeophys/CanopyHydrologyMod.F90 b/components/elm/src/biogeophys/CanopyHydrologyMod.F90 index 5381b1358f8e..17debec4bcba 100755 --- a/components/elm/src/biogeophys/CanopyHydrologyMod.F90 +++ b/components/elm/src/biogeophys/CanopyHydrologyMod.F90 @@ -16,7 +16,7 @@ module CanopyHydrologyMod use shr_sys_mod , only : shr_sys_flush use decompMod , only : bounds_type use abortutils , only : endrun - use elm_varctl , only : iulog, tw_irr, extra_gw_irr, irrigate, use_firn_percolation_and_compaction + use elm_varctl , only : iulog, tw_irr, extra_gw_irr, irrigate use LandunitType , only : lun_pp use atm2lndType , only : atm2lnd_type use AerosolType , only : aerosol_type @@ -433,7 +433,7 @@ subroutine CanopyHydrology(bounds, & qflx_prec_grnd(p) = qflx_prec_grnd_snow(p) + qflx_prec_grnd_rain(p) - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then if (do_capsnow(c)) then qflx_snwcp_liq(p) = qflx_prec_grnd_rain(p) qflx_snwcp_ice(p) = qflx_prec_grnd_snow(p) @@ -491,7 +491,7 @@ subroutine CanopyHydrology(bounds, & ! Determine snow height and snow water - if (use_firn_percolation_and_compaction) then + if (use_extrasnowlayers) then call NewSnowBulkDensity(bounds, num_nolakec, filter_nolakec, & top_as, bifall(bounds%begc:bounds%endc)) end if @@ -517,13 +517,13 @@ subroutine CanopyHydrology(bounds, & swe_old(c,j)=h2osoi_liq(c,j)+h2osoi_ice(c,j) enddo - if (do_capsnow(c) .and. .not. use_firn_percolation_and_compaction) then + if (do_capsnow(c) .and. .not. use_extrasnowlayers) then dz_snowf = 0._r8 newsnow(c) = qflx_snow_grnd_col(c) * dtime frac_sno(c)=1._r8 int_snow(c) = 5.e2_r8 else - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then if (forc_t(t) > tfrz + 2._r8) then bifall(c)=50._r8 + 1.7_r8*(17.0_r8)**1.5_r8 else if (forc_t(t) > tfrz - 15._r8) then @@ -662,7 +662,7 @@ subroutine CanopyHydrology(bounds, & ! as the surface air temperature newnode = 0 ! flag for when snow node will be initialized - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then if (snl(c) == 0 .and. qflx_snow_grnd_col(c) > 0.0_r8 .and. frac_sno(c)*snow_depth(c) >= 0.01_r8) then newnode = 1 snl(c) = -1 diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index da4006e9771c..d44f20c17eda 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -7,7 +7,7 @@ module HydrologyDrainageMod use shr_kind_mod , only : r8 => shr_kind_r8 use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type - use elm_varctl , only : iulog, use_vichydro, use_firn_percolation_and_compaction + use elm_varctl , only : iulog, use_vichydro use elm_varcon , only : e_ice, denh2o, denice, rpi, spval use atm2lndType , only : atm2lnd_type use glc2lndMod , only : glc2lnd_type diff --git a/components/elm/src/biogeophys/HydrologyNoDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyNoDrainageMod.F90 index fa6fe25c4ae7..d6bf3b34d136 100644 --- a/components/elm/src/biogeophys/HydrologyNoDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyNoDrainageMod.F90 @@ -7,7 +7,7 @@ Module HydrologyNoDrainageMod use shr_kind_mod , only : r8 => shr_kind_r8 use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type - use elm_varctl , only : iulog, use_vichydro, use_extrasnowlayers, use_firn_percolation_and_compaction + use elm_varctl , only : iulog, use_vichydro, use_extrasnowlayers use elm_varcon , only : e_ice, denh2o, denice, rpi, spval use atm2lndType , only : atm2lnd_type use lnd2atmType , only : lnd2atm_type @@ -300,7 +300,7 @@ subroutine HydrologyNoDrainage(bounds, & endif #endif - if (use_firn_percolation_and_compaction) then + if (use_extrasnowlayers) then call SnowCapping(bounds, num_nolakec, filter_nolakec, num_snowc, filter_snowc, & aerosol_vars) end if diff --git a/components/elm/src/biogeophys/LakeFluxesMod.F90 b/components/elm/src/biogeophys/LakeFluxesMod.F90 index 92733834c775..30e05a57a9be 100644 --- a/components/elm/src/biogeophys/LakeFluxesMod.F90 +++ b/components/elm/src/biogeophys/LakeFluxesMod.F90 @@ -52,7 +52,7 @@ subroutine LakeFluxes(bounds, num_lakec, filter_lakec, num_lakep, filter_lakep, use elm_varpar , only : nlevlak use elm_varcon , only : hvap, hsub, hfus, cpair, cpliq, tkwat, tkice, tkair use elm_varcon , only : sb, vkc, grav, denh2o, tfrz, spval, zsno - use elm_varctl , only : iulog, use_lch4, use_firn_percolation_and_compaction + use elm_varctl , only : iulog, use_lch4, use_extrasnowlayers use LakeCon , only : betavis, z0frzlake, tdmax, emg_lake use LakeCon , only : lake_use_old_fcrit_minz0 use LakeCon , only : minz0lake, cur0, cus, curm, fcrit @@ -695,7 +695,7 @@ subroutine LakeFluxes(bounds, num_lakec, filter_lakec, num_lakep, filter_lakep, qflx_dirct_rain(p) = 0._r8 qflx_leafdrip(p) = 0._r8 - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then ! Because they will be used in pft2col initialize here. ! This will be overwritten in LakeHydrology qflx_snwcp_ice(p) = 0._r8 diff --git a/components/elm/src/biogeophys/LakeHydrologyMod.F90 b/components/elm/src/biogeophys/LakeHydrologyMod.F90 index f0bf1d373fcb..7391a2d80fb0 100644 --- a/components/elm/src/biogeophys/LakeHydrologyMod.F90 +++ b/components/elm/src/biogeophys/LakeHydrologyMod.F90 @@ -74,7 +74,7 @@ subroutine LakeHydrology(bounds, & !$acc routine seq use elm_varcon , only : denh2o, denice, spval, hfus, tfrz, cpliq, cpice use elm_varpar , only : nlevsno, nlevgrnd, nlevsoi - use elm_varctl , only : iulog, use_extrasnowlayers, use_lake_wat_storage, use_firn_percolation_and_compaction + use elm_varctl , only : iulog, use_extrasnowlayers, use_lake_wat_storage use elm_time_manager, only : get_step_size use SnowHydrologyMod, only : SnowCompaction, CombineSnowLayers, SnowWater, BuildSnowFilter use SnowHydrologyMod, only : DivideSnowLayers, DivideExtraSnowLayers, SnowCapping @@ -256,7 +256,7 @@ subroutine LakeHydrology(bounds, & qflx_dirct_rain(p) = 0._r8 qflx_leafdrip(p) = 0._r8 - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then if (do_capsnow(c)) then qflx_snwcp_ice(p) = qflx_prec_grnd_snow(p) qflx_snwcp_liq(p) = qflx_prec_grnd_rain(p) @@ -288,7 +288,7 @@ subroutine LakeHydrology(bounds, & ! U.S.Department of Agriculture Forest Service, Project F, ! Progress Rep. 1, Alta Avalanche Study Center:Snow Layer Densification. - if (do_capsnow(c) .and. .not. use_firn_percolation_and_compaction) then + if (do_capsnow(c) .and. .not. use_extrasnowlayers) then dz_snowf = 0._r8 else if (forc_t(t) > tfrz + 2._r8) then @@ -400,7 +400,7 @@ subroutine LakeHydrology(bounds, & ! Update the pft-level qflx_snowcap ! This was moved in from Hydrology2 to keep all pft-level ! calculations out of Hydrology2 - if (do_capsnow(c) .and. .not. use_firn_percolation_and_compaction) then + if (do_capsnow(c) .and. .not. use_extrasnowlayers) then qflx_snwcp_ice(p) = qflx_snwcp_ice(p) + qflx_dew_snow(p) qflx_snwcp_liq(p) = qflx_snwcp_liq(p) + qflx_dew_grnd(p) end if @@ -422,7 +422,7 @@ subroutine LakeHydrology(bounds, & ! Update snow pack for dew & sub. h2osno_temp = h2osno(c) - if (do_capsnow(c) .and. .not. use_firn_percolation_and_compaction) then + if (do_capsnow(c) .and. .not. use_extrasnowlayers) then h2osno(c) = h2osno(c) - qflx_sub_snow(p)*dtime qflx_snwcp_ice(p) = qflx_snwcp_ice(p) + qflx_dew_snow(p) qflx_snwcp_liq(p) = qflx_snwcp_liq(p) + qflx_dew_grnd(p) @@ -437,7 +437,7 @@ subroutine LakeHydrology(bounds, & end if end if - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then qflx_snwcp_ice_col(c) = qflx_snwcp_ice(p) qflx_snwcp_liq_col(c) = qflx_snwcp_liq(p) end if @@ -485,7 +485,7 @@ subroutine LakeHydrology(bounds, & num_shlakesnowc, filter_shlakesnowc, num_shlakenosnowc, filter_shlakenosnowc, & atm2lnd_vars, aerosol_vars) - if (use_firn_percolation_and_compaction) then + if (use_extrasnowlayers) then call SnowCapping(bounds, num_lakec, filter_lakec, num_shlakesnowc, filter_shlakesnowc, & aerosol_vars) end if @@ -719,7 +719,7 @@ subroutine LakeHydrology(bounds, & if (use_lake_wat_storage) then qflx_qrgwl(c) = 0._r8 if (wslake(c) >= 5000._r8) then - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then qflx_snwcp = qflx_snwcp_ice(p) else qflx_snwcp = qflx_snwcp_ice_col(c) @@ -732,7 +732,7 @@ subroutine LakeHydrology(bounds, & (endwb(c) - begwb(c)) endwb(c) = endwb(c) + wslake(c) else - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then qflx_snwcp = qflx_snwcp_ice(p) else qflx_snwcp = qflx_snwcp_ice_col(c) diff --git a/components/elm/src/biogeophys/SnowHydrologyMod.F90 b/components/elm/src/biogeophys/SnowHydrologyMod.F90 index 57e486d0ce90..46e1d7089b2c 100644 --- a/components/elm/src/biogeophys/SnowHydrologyMod.F90 +++ b/components/elm/src/biogeophys/SnowHydrologyMod.F90 @@ -18,7 +18,7 @@ module SnowHydrologyMod use decompMod , only : bounds_type use abortutils , only : endrun use elm_varpar , only : nlevsno - use elm_varctl , only : iulog, use_extrasnowlayers, use_firn_percolation_and_compaction + use elm_varctl , only : iulog, use_extrasnowlayers use elm_varcon , only : namec, h2osno_max use atm2lndType , only : atm2lnd_type use AerosolType , only : aerosol_type @@ -233,7 +233,7 @@ subroutine SnowWater(bounds, & c = filter_snowc(fc) l=col_pp%landunit(c) - if (do_capsnow(c) .and. .not. use_firn_percolation_and_compaction) then + if (do_capsnow(c) .and. .not. use_extrasnowlayers) then wgdif = h2osoi_ice(c,snl(c)+1) - frac_sno_eff(c)*qflx_sub_snow(c)*dtime h2osoi_ice(c,snl(c)+1) = wgdif if (wgdif < 0._r8) then @@ -623,7 +623,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & ! Begin calculation - note that the following column loops are only invoked if snl(c) < 0 - if (use_firn_percolation_and_compaction) then + if (use_extrasnowlayers) then do fc = 1, num_snowc c = filter_snowc(fc) burden(c) = 0._r8 @@ -637,7 +637,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & do j = -nlevsno+1, 0 do fc = 1, num_snowc c = filter_snowc(fc) - if (use_firn_percolation_and_compaction) then + if (use_extrasnowlayers) then t = col_pp%topounit(c) end if if (j >= snl(c)+1) then @@ -647,7 +647,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & ! If void is negative, then increase dz such that void = 0. ! This should be done for any landunit, but for now is done only for glacier_mec 1andunits. - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then ! I don't think the next 5 lines are necessary (removed in CLMv5) l = col_pp%landunit(c) if (ltype(l)==istice_mec .and. void < 0._r8) then @@ -665,7 +665,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & dexpf = exp(-c4*td) ! Settling as a result of destructive metamorphism - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then ddz1 = -c3*dexpf if (bi > dm) ddz1 = ddz1*exp(-46.0e-3_r8*(bi-dm)) else @@ -685,7 +685,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & if (h2osoi_liq(c,j) > 0.01_r8*dz(c,j)*frac_sno(c)) ddz1=ddz1*c5 ! Compaction due to overburden - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then ddz2 = -(burden(c)+wx/2._r8)*exp(-0.08_r8*td - c2*bi)/eta0 else p_gls = max(denice / bi, 1._r8) * grav * (burden(c) + wx/2._r8) @@ -724,7 +724,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & ddz3 = 0._r8 end if - if (use_firn_percolation_and_compaction) then + if (use_extrasnowlayers) then ! Compaction occurring due to wind drift call WindDriftCompaction( & bi = bi, & diff --git a/components/elm/src/biogeophys/SnowSnicarMod.F90 b/components/elm/src/biogeophys/SnowSnicarMod.F90 index de161c1d4c9d..d6d6de75ffaa 100644 --- a/components/elm/src/biogeophys/SnowSnicarMod.F90 +++ b/components/elm/src/biogeophys/SnowSnicarMod.F90 @@ -11,7 +11,7 @@ module SnowSnicarMod use shr_kind_mod , only : r8 => shr_kind_r8 use shr_sys_mod , only : shr_sys_flush use shr_log_mod , only : errMsg => shr_log_errMsg - use elm_varctl , only : iulog, use_firn_percolation_and_compaction + use elm_varctl , only : iulog, use_extrasnowlayers use elm_varcon , only : namec use shr_const_mod , only : SHR_CONST_RHOICE use abortutils , only : endrun @@ -1435,13 +1435,13 @@ subroutine SnowAge_grain(bounds, & ! RE-FREEZING ! ! new snowfall [kg/m2] - if (do_capsnow(c_idx) .and. .not. use_firn_percolation_and_compaction) then + if (do_capsnow(c_idx) .and. .not. use_extrasnowlayers) then newsnow = max(0._r8, (qflx_snwcp_ice(c_idx)*dtime)) else newsnow = max(0._r8, (qflx_snow_grnd_col(c_idx)*dtime)) endif - if (use_firn_percolation_and_compaction) then + if (use_extrasnowlayers) then snw_rds_refrz = 1500._r8 else snw_rds_refrz = 1000._r8 diff --git a/components/elm/src/biogeophys/SoilFluxesMod.F90 b/components/elm/src/biogeophys/SoilFluxesMod.F90 index 1b221a617a83..080e6aceb394 100644 --- a/components/elm/src/biogeophys/SoilFluxesMod.F90 +++ b/components/elm/src/biogeophys/SoilFluxesMod.F90 @@ -8,7 +8,7 @@ module SoilFluxesMod use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type use abortutils , only : endrun - use elm_varctl , only : iulog, use_firn_percolation_and_compaction + use elm_varctl , only : iulog use perfMod_GPU use elm_varpar , only : nlevsno, nlevgrnd, nlevurb, max_patch_per_col use atm2lndType , only : atm2lnd_type @@ -349,7 +349,7 @@ subroutine SoilFluxes (bounds, num_urbanl, filter_urbanl, & ! This was moved in from Hydrology2 to keep all pft-level ! calculations out of Hydrology2 - if (col_pp%snl(c) < 0 .and. do_capsnow(c) .and. .not. use_firn_percolation_and_compaction) then + if (col_pp%snl(c) < 0 .and. do_capsnow(c) .and. .not. use_extrasnowlayers) then qflx_snwcp_liq(p) = qflx_snwcp_liq(p)+frac_sno_eff(c)*qflx_dew_grnd(p) qflx_snwcp_ice(p) = qflx_snwcp_ice(p)+frac_sno_eff(c)*qflx_dew_snow(p) end if diff --git a/components/elm/src/main/controlMod.F90 b/components/elm/src/main/controlMod.F90 index dc7fa975dbb4..9db909a83fbb 100755 --- a/components/elm/src/main/controlMod.F90 +++ b/components/elm/src/main/controlMod.F90 @@ -289,7 +289,7 @@ subroutine control_init( ) namelist /elm_inparm/ & use_nofire, use_lch4, use_vertsoilc, use_extralakelayers, & use_vichydro, use_century_decomp, use_cn, use_crop, use_snicar_frc, & - use_snicar_ad, use_firn_percolation_and_compaction, use_extrasnowlayers,& + use_snicar_ad, use_extrasnowlayers,& use_T_rho_dependent_snowthk, use_vancouver, use_mexicocity, use_noio ! cpl_bypass variables @@ -731,7 +731,6 @@ subroutine control_spmd() call mpi_bcast (use_vertsoilc, 1, MPI_LOGICAL, 0, mpicom, ier) call mpi_bcast (use_extralakelayers, 1, MPI_LOGICAL, 0, mpicom, ier) call mpi_bcast (use_extrasnowlayers, 1, MPI_LOGICAL, 0, mpicom, ier) - call mpi_bcast (use_firn_percolation_and_compaction, 1, MPI_LOGICAL, 0, mpicom, ier) call mpi_bcast (use_T_rho_dependent_snowthk, 1, MPI_LOGICAL, 0, mpicom, ier) call mpi_bcast (use_vichydro, 1, MPI_LOGICAL, 0, mpicom, ier) call mpi_bcast (use_century_decomp, 1, MPI_LOGICAL, 0, mpicom, ier) @@ -1037,7 +1036,6 @@ subroutine control_print () write(iulog,*) ' use_lake_wat_storage = ', use_lake_wat_storage write(iulog,*) ' use_extralakelayers = ', use_extralakelayers write(iulog,*) ' use_extrasnowlayers = ', use_extrasnowlayers - write(iulog,*) ' use_firn_percolation_and_compaction = ', use_firn_percolation_and_compaction write(iulog,*) ' use_T_rho_dependent_snowthk = ', use_T_rho_dependent_snowthk write(iulog,*) ' use_vichydro = ', use_vichydro write(iulog,*) ' use_century_decomp = ', use_century_decomp diff --git a/components/elm/src/main/elm_driver.F90 b/components/elm/src/main/elm_driver.F90 index 51f08bf235c3..322dd3775b61 100644 --- a/components/elm/src/main/elm_driver.F90 +++ b/components/elm/src/main/elm_driver.F90 @@ -12,7 +12,7 @@ module elm_driver use shr_sys_mod , only : shr_sys_flush use shr_log_mod , only : errMsg => shr_log_errMsg use elm_varpar , only : nlevtrc_soil, nlevsoi - use elm_varctl , only : wrtdia, iulog, create_glacier_mec_landunit, use_fates, use_betr, use_firn_percolation_and_compaction + use elm_varctl , only : wrtdia, iulog, create_glacier_mec_landunit, use_fates, use_betr use elm_varctl , only : use_cn, use_lch4, use_voc, use_noio, use_c13, use_c14 use elm_varctl , only : use_erosion, use_fates_sp, use_fan use elm_varctl , only : mpi_sync_nstep_freq @@ -1622,7 +1622,7 @@ subroutine elm_drv_init(bounds, & ! Save snow mass at previous time step h2osno_old(c) = h2osno(c) - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then ! Decide whether to cap snow if (h2osno(c) > h2osno_max) then do_capsnow(c) = .true. @@ -1790,7 +1790,7 @@ subroutine elm_drv_patch2col (bounds, num_nolakec, filter_nolakec, & qflx_snow_grnd_patch(bounds%begp:bounds%endp), & qflx_snow_grnd_col (bounds%begc:bounds%endc)) - if (.not. use_firn_percolation_and_compaction) then + if (.not. use_extrasnowlayers) then call p2c (bounds, num_allc, filter_allc, & veg_wf%qflx_snwcp_liq(bounds%begp:bounds%endp), & col_wf%qflx_snwcp_liq(bounds%begc:bounds%endc)) diff --git a/components/elm/src/main/elm_varctl.F90 b/components/elm/src/main/elm_varctl.F90 index fb9ca0dbd24e..f1a1733391fd 100644 --- a/components/elm/src/main/elm_varctl.F90 +++ b/components/elm/src/main/elm_varctl.F90 @@ -375,7 +375,6 @@ module elm_varctl logical, public :: use_snicar_frc = .false. logical, public :: use_snicar_ad = .false. logical, public :: use_extrasnowlayers = .false. - logical, public :: use_firn_percolation_and_compaction = .false. logical, public :: use_vancouver = .false. logical, public :: use_mexicocity = .false. logical, public :: use_noio = .false. From ee58bb0e9d690e016e3906f91a90da68145c11ee Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 17 Dec 2024 08:53:20 -0800 Subject: [PATCH 10/39] minor mods to remove use_firn... for use_extrasnowlayers --- components/elm/src/biogeophys/CanopyHydrologyMod.F90 | 2 +- components/elm/src/biogeophys/SoilFluxesMod.F90 | 2 +- components/elm/src/main/elm_driver.F90 | 2 +- 3 files changed, 3 insertions(+), 3 deletions(-) diff --git a/components/elm/src/biogeophys/CanopyHydrologyMod.F90 b/components/elm/src/biogeophys/CanopyHydrologyMod.F90 index 17debec4bcba..e24e9b486d88 100755 --- a/components/elm/src/biogeophys/CanopyHydrologyMod.F90 +++ b/components/elm/src/biogeophys/CanopyHydrologyMod.F90 @@ -116,7 +116,7 @@ subroutine CanopyHydrology(bounds, & use elm_varcon , only : hfus, denice, zlnd, rpi, spval, tfrz use column_varcon , only : icol_roof, icol_sunwall, icol_shadewall use landunit_varcon , only : istcrop, istice, istwet, istsoil, istice_mec, istdlak - use elm_varctl , only : subgridflag + use elm_varctl , only : subgridflag, use_extrasnowlayers use elm_varpar , only : nlevsoi,nlevsno use elm_varsur , only : wt_lunit, f_grd, f_surf use atm2lndType , only : atm2lnd_type diff --git a/components/elm/src/biogeophys/SoilFluxesMod.F90 b/components/elm/src/biogeophys/SoilFluxesMod.F90 index 080e6aceb394..a83ae505ff7c 100644 --- a/components/elm/src/biogeophys/SoilFluxesMod.F90 +++ b/components/elm/src/biogeophys/SoilFluxesMod.F90 @@ -8,7 +8,7 @@ module SoilFluxesMod use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type use abortutils , only : endrun - use elm_varctl , only : iulog + use elm_varctl , only : iulog, use_extrasnowlayers use perfMod_GPU use elm_varpar , only : nlevsno, nlevgrnd, nlevurb, max_patch_per_col use atm2lndType , only : atm2lnd_type diff --git a/components/elm/src/main/elm_driver.F90 b/components/elm/src/main/elm_driver.F90 index 322dd3775b61..06668f6cad07 100644 --- a/components/elm/src/main/elm_driver.F90 +++ b/components/elm/src/main/elm_driver.F90 @@ -14,7 +14,7 @@ module elm_driver use elm_varpar , only : nlevtrc_soil, nlevsoi use elm_varctl , only : wrtdia, iulog, create_glacier_mec_landunit, use_fates, use_betr use elm_varctl , only : use_cn, use_lch4, use_voc, use_noio, use_c13, use_c14 - use elm_varctl , only : use_erosion, use_fates_sp, use_fan + use elm_varctl , only : use_erosion, use_fates_sp, use_fan,use_extrasnowlayers use elm_varctl , only : mpi_sync_nstep_freq use elm_time_manager , only : get_step_size, get_curr_date, get_ref_date, get_nstep, is_beg_curr_day, get_curr_time_string use elm_time_manager , only : get_curr_calday, get_days_per_year From 09824c82e4356260ea97b1fb2540e7d2b2b9ac59 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 7 Jan 2025 12:37:58 -0800 Subject: [PATCH 11/39] Revert "remove use_firn_percolation_and_compaction for use_extrasnowlayers" This reverts commit c230f3a03fa942cc280f8fc832280e63a85db9d9. --- .../bld/namelist_files/namelist_definition.xml | 7 +++++++ components/elm/src/biogeophys/AerosolMod.F90 | 4 ++-- .../elm/src/biogeophys/BalanceCheckMod.F90 | 6 +++--- .../elm/src/biogeophys/CanopyHydrologyMod.F90 | 12 ++++++------ .../src/biogeophys/HydrologyDrainageMod.F90 | 2 +- .../src/biogeophys/HydrologyNoDrainageMod.F90 | 4 ++-- .../elm/src/biogeophys/LakeFluxesMod.F90 | 4 ++-- .../elm/src/biogeophys/LakeHydrologyMod.F90 | 18 +++++++++--------- .../elm/src/biogeophys/SnowHydrologyMod.F90 | 16 ++++++++-------- .../elm/src/biogeophys/SnowSnicarMod.F90 | 6 +++--- .../elm/src/biogeophys/SoilFluxesMod.F90 | 4 ++-- components/elm/src/main/controlMod.F90 | 4 +++- components/elm/src/main/elm_driver.F90 | 6 +++--- components/elm/src/main/elm_varctl.F90 | 1 + 14 files changed, 52 insertions(+), 42 deletions(-) diff --git a/components/elm/bld/namelist_files/namelist_definition.xml b/components/elm/bld/namelist_files/namelist_definition.xml index f6c45e05859b..c76271b7baf7 100644 --- a/components/elm/bld/namelist_files/namelist_definition.xml +++ b/components/elm/bld/namelist_files/namelist_definition.xml @@ -766,6 +766,13 @@ snowpack conditions on glaciers and ice sheets, including a semi-empirical firn densification model. + +This option enables more realistic +snowpack conditions on glaciers and ice sheets, including a semi-empirical +firn densification model. uses the same physics as use_extrasnowlayers but +uses the original 5 snow layer scheme + Toggle to use new snow thermal conductivity that relies on snow temperature and density. diff --git a/components/elm/src/biogeophys/AerosolMod.F90 b/components/elm/src/biogeophys/AerosolMod.F90 index d97c48dfd29b..81fc8d62abc5 100644 --- a/components/elm/src/biogeophys/AerosolMod.F90 +++ b/components/elm/src/biogeophys/AerosolMod.F90 @@ -5,7 +5,7 @@ module AerosolMod use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type use elm_varpar , only : nlevsno - use elm_varctl , only : use_extrasnowlayers + use elm_varctl , only : use_firn_percolation_and_compaction use elm_time_manager , only : get_step_size use atm2lndType , only : atm2lnd_type use AerosolType , only : aerosol_type @@ -107,7 +107,7 @@ subroutine AerosolMasses(bounds, num_on, filter_on, num_off, filter_off, aerosol ! layer mass of snow: snowmass = h2osoi_ice(c,j) + h2osoi_liq(c,j) - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then ! Correct the top layer aerosol mass to account for snow capping. ! This approach conserves the aerosol mass concentration ! (but not the aerosol amss) when snow-capping is invoked diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index ae26ef0a96a2..3ce501475dc5 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -9,7 +9,7 @@ module BalanceCheckMod use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type use abortutils , only : endrun - use elm_varctl , only : iulog, use_var_soil_thick, use_extrasnowlayers + use elm_varctl , only : iulog, use_var_soil_thick, use_firn_percolation_and_compaction use elm_varcon , only : namep, namec use GetGlobalValuesMod , only : GetGlobalIndex use atm2lndType , only : atm2lnd_type @@ -448,7 +448,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) + qflx_sl_top_soil(c) if (lun_pp%itype(l) == istdlak) then - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then if ( do_capsnow(c)) then snow_sources(c) = qflx_snow_grnd_col(c) & + frac_sno_eff(c) * (qflx_dew_snow(c) + qflx_dew_grnd(c) ) @@ -476,7 +476,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & endif if (lun_pp%itype(l) == istsoil .or. lun_pp%itype(l) == istcrop .or. lun_pp%itype(l) == istwet ) then - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then if (do_capsnow(c)) then snow_sources(c) = frac_sno_eff(c) * (qflx_dew_snow(c) + qflx_dew_grnd(c) ) & + qflx_h2osfc_to_ice(c) + qflx_prec_grnd(c) diff --git a/components/elm/src/biogeophys/CanopyHydrologyMod.F90 b/components/elm/src/biogeophys/CanopyHydrologyMod.F90 index e24e9b486d88..5838fd448002 100755 --- a/components/elm/src/biogeophys/CanopyHydrologyMod.F90 +++ b/components/elm/src/biogeophys/CanopyHydrologyMod.F90 @@ -16,7 +16,7 @@ module CanopyHydrologyMod use shr_sys_mod , only : shr_sys_flush use decompMod , only : bounds_type use abortutils , only : endrun - use elm_varctl , only : iulog, tw_irr, extra_gw_irr, irrigate + use elm_varctl , only : iulog, tw_irr, extra_gw_irr, irrigate, use_firn_percolation_and_compaction use LandunitType , only : lun_pp use atm2lndType , only : atm2lnd_type use AerosolType , only : aerosol_type @@ -433,7 +433,7 @@ subroutine CanopyHydrology(bounds, & qflx_prec_grnd(p) = qflx_prec_grnd_snow(p) + qflx_prec_grnd_rain(p) - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then if (do_capsnow(c)) then qflx_snwcp_liq(p) = qflx_prec_grnd_rain(p) qflx_snwcp_ice(p) = qflx_prec_grnd_snow(p) @@ -491,7 +491,7 @@ subroutine CanopyHydrology(bounds, & ! Determine snow height and snow water - if (use_extrasnowlayers) then + if (use_firn_percolation_and_compaction) then call NewSnowBulkDensity(bounds, num_nolakec, filter_nolakec, & top_as, bifall(bounds%begc:bounds%endc)) end if @@ -517,13 +517,13 @@ subroutine CanopyHydrology(bounds, & swe_old(c,j)=h2osoi_liq(c,j)+h2osoi_ice(c,j) enddo - if (do_capsnow(c) .and. .not. use_extrasnowlayers) then + if (do_capsnow(c) .and. .not. use_firn_percolation_and_compaction) then dz_snowf = 0._r8 newsnow(c) = qflx_snow_grnd_col(c) * dtime frac_sno(c)=1._r8 int_snow(c) = 5.e2_r8 else - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then if (forc_t(t) > tfrz + 2._r8) then bifall(c)=50._r8 + 1.7_r8*(17.0_r8)**1.5_r8 else if (forc_t(t) > tfrz - 15._r8) then @@ -662,7 +662,7 @@ subroutine CanopyHydrology(bounds, & ! as the surface air temperature newnode = 0 ! flag for when snow node will be initialized - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then if (snl(c) == 0 .and. qflx_snow_grnd_col(c) > 0.0_r8 .and. frac_sno(c)*snow_depth(c) >= 0.01_r8) then newnode = 1 snl(c) = -1 diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index d44f20c17eda..da4006e9771c 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -7,7 +7,7 @@ module HydrologyDrainageMod use shr_kind_mod , only : r8 => shr_kind_r8 use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type - use elm_varctl , only : iulog, use_vichydro + use elm_varctl , only : iulog, use_vichydro, use_firn_percolation_and_compaction use elm_varcon , only : e_ice, denh2o, denice, rpi, spval use atm2lndType , only : atm2lnd_type use glc2lndMod , only : glc2lnd_type diff --git a/components/elm/src/biogeophys/HydrologyNoDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyNoDrainageMod.F90 index d6bf3b34d136..fa6fe25c4ae7 100644 --- a/components/elm/src/biogeophys/HydrologyNoDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyNoDrainageMod.F90 @@ -7,7 +7,7 @@ Module HydrologyNoDrainageMod use shr_kind_mod , only : r8 => shr_kind_r8 use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type - use elm_varctl , only : iulog, use_vichydro, use_extrasnowlayers + use elm_varctl , only : iulog, use_vichydro, use_extrasnowlayers, use_firn_percolation_and_compaction use elm_varcon , only : e_ice, denh2o, denice, rpi, spval use atm2lndType , only : atm2lnd_type use lnd2atmType , only : lnd2atm_type @@ -300,7 +300,7 @@ subroutine HydrologyNoDrainage(bounds, & endif #endif - if (use_extrasnowlayers) then + if (use_firn_percolation_and_compaction) then call SnowCapping(bounds, num_nolakec, filter_nolakec, num_snowc, filter_snowc, & aerosol_vars) end if diff --git a/components/elm/src/biogeophys/LakeFluxesMod.F90 b/components/elm/src/biogeophys/LakeFluxesMod.F90 index 30e05a57a9be..92733834c775 100644 --- a/components/elm/src/biogeophys/LakeFluxesMod.F90 +++ b/components/elm/src/biogeophys/LakeFluxesMod.F90 @@ -52,7 +52,7 @@ subroutine LakeFluxes(bounds, num_lakec, filter_lakec, num_lakep, filter_lakep, use elm_varpar , only : nlevlak use elm_varcon , only : hvap, hsub, hfus, cpair, cpliq, tkwat, tkice, tkair use elm_varcon , only : sb, vkc, grav, denh2o, tfrz, spval, zsno - use elm_varctl , only : iulog, use_lch4, use_extrasnowlayers + use elm_varctl , only : iulog, use_lch4, use_firn_percolation_and_compaction use LakeCon , only : betavis, z0frzlake, tdmax, emg_lake use LakeCon , only : lake_use_old_fcrit_minz0 use LakeCon , only : minz0lake, cur0, cus, curm, fcrit @@ -695,7 +695,7 @@ subroutine LakeFluxes(bounds, num_lakec, filter_lakec, num_lakep, filter_lakep, qflx_dirct_rain(p) = 0._r8 qflx_leafdrip(p) = 0._r8 - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then ! Because they will be used in pft2col initialize here. ! This will be overwritten in LakeHydrology qflx_snwcp_ice(p) = 0._r8 diff --git a/components/elm/src/biogeophys/LakeHydrologyMod.F90 b/components/elm/src/biogeophys/LakeHydrologyMod.F90 index 7391a2d80fb0..f0bf1d373fcb 100644 --- a/components/elm/src/biogeophys/LakeHydrologyMod.F90 +++ b/components/elm/src/biogeophys/LakeHydrologyMod.F90 @@ -74,7 +74,7 @@ subroutine LakeHydrology(bounds, & !$acc routine seq use elm_varcon , only : denh2o, denice, spval, hfus, tfrz, cpliq, cpice use elm_varpar , only : nlevsno, nlevgrnd, nlevsoi - use elm_varctl , only : iulog, use_extrasnowlayers, use_lake_wat_storage + use elm_varctl , only : iulog, use_extrasnowlayers, use_lake_wat_storage, use_firn_percolation_and_compaction use elm_time_manager, only : get_step_size use SnowHydrologyMod, only : SnowCompaction, CombineSnowLayers, SnowWater, BuildSnowFilter use SnowHydrologyMod, only : DivideSnowLayers, DivideExtraSnowLayers, SnowCapping @@ -256,7 +256,7 @@ subroutine LakeHydrology(bounds, & qflx_dirct_rain(p) = 0._r8 qflx_leafdrip(p) = 0._r8 - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then if (do_capsnow(c)) then qflx_snwcp_ice(p) = qflx_prec_grnd_snow(p) qflx_snwcp_liq(p) = qflx_prec_grnd_rain(p) @@ -288,7 +288,7 @@ subroutine LakeHydrology(bounds, & ! U.S.Department of Agriculture Forest Service, Project F, ! Progress Rep. 1, Alta Avalanche Study Center:Snow Layer Densification. - if (do_capsnow(c) .and. .not. use_extrasnowlayers) then + if (do_capsnow(c) .and. .not. use_firn_percolation_and_compaction) then dz_snowf = 0._r8 else if (forc_t(t) > tfrz + 2._r8) then @@ -400,7 +400,7 @@ subroutine LakeHydrology(bounds, & ! Update the pft-level qflx_snowcap ! This was moved in from Hydrology2 to keep all pft-level ! calculations out of Hydrology2 - if (do_capsnow(c) .and. .not. use_extrasnowlayers) then + if (do_capsnow(c) .and. .not. use_firn_percolation_and_compaction) then qflx_snwcp_ice(p) = qflx_snwcp_ice(p) + qflx_dew_snow(p) qflx_snwcp_liq(p) = qflx_snwcp_liq(p) + qflx_dew_grnd(p) end if @@ -422,7 +422,7 @@ subroutine LakeHydrology(bounds, & ! Update snow pack for dew & sub. h2osno_temp = h2osno(c) - if (do_capsnow(c) .and. .not. use_extrasnowlayers) then + if (do_capsnow(c) .and. .not. use_firn_percolation_and_compaction) then h2osno(c) = h2osno(c) - qflx_sub_snow(p)*dtime qflx_snwcp_ice(p) = qflx_snwcp_ice(p) + qflx_dew_snow(p) qflx_snwcp_liq(p) = qflx_snwcp_liq(p) + qflx_dew_grnd(p) @@ -437,7 +437,7 @@ subroutine LakeHydrology(bounds, & end if end if - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then qflx_snwcp_ice_col(c) = qflx_snwcp_ice(p) qflx_snwcp_liq_col(c) = qflx_snwcp_liq(p) end if @@ -485,7 +485,7 @@ subroutine LakeHydrology(bounds, & num_shlakesnowc, filter_shlakesnowc, num_shlakenosnowc, filter_shlakenosnowc, & atm2lnd_vars, aerosol_vars) - if (use_extrasnowlayers) then + if (use_firn_percolation_and_compaction) then call SnowCapping(bounds, num_lakec, filter_lakec, num_shlakesnowc, filter_shlakesnowc, & aerosol_vars) end if @@ -719,7 +719,7 @@ subroutine LakeHydrology(bounds, & if (use_lake_wat_storage) then qflx_qrgwl(c) = 0._r8 if (wslake(c) >= 5000._r8) then - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then qflx_snwcp = qflx_snwcp_ice(p) else qflx_snwcp = qflx_snwcp_ice_col(c) @@ -732,7 +732,7 @@ subroutine LakeHydrology(bounds, & (endwb(c) - begwb(c)) endwb(c) = endwb(c) + wslake(c) else - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then qflx_snwcp = qflx_snwcp_ice(p) else qflx_snwcp = qflx_snwcp_ice_col(c) diff --git a/components/elm/src/biogeophys/SnowHydrologyMod.F90 b/components/elm/src/biogeophys/SnowHydrologyMod.F90 index 46e1d7089b2c..57e486d0ce90 100644 --- a/components/elm/src/biogeophys/SnowHydrologyMod.F90 +++ b/components/elm/src/biogeophys/SnowHydrologyMod.F90 @@ -18,7 +18,7 @@ module SnowHydrologyMod use decompMod , only : bounds_type use abortutils , only : endrun use elm_varpar , only : nlevsno - use elm_varctl , only : iulog, use_extrasnowlayers + use elm_varctl , only : iulog, use_extrasnowlayers, use_firn_percolation_and_compaction use elm_varcon , only : namec, h2osno_max use atm2lndType , only : atm2lnd_type use AerosolType , only : aerosol_type @@ -233,7 +233,7 @@ subroutine SnowWater(bounds, & c = filter_snowc(fc) l=col_pp%landunit(c) - if (do_capsnow(c) .and. .not. use_extrasnowlayers) then + if (do_capsnow(c) .and. .not. use_firn_percolation_and_compaction) then wgdif = h2osoi_ice(c,snl(c)+1) - frac_sno_eff(c)*qflx_sub_snow(c)*dtime h2osoi_ice(c,snl(c)+1) = wgdif if (wgdif < 0._r8) then @@ -623,7 +623,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & ! Begin calculation - note that the following column loops are only invoked if snl(c) < 0 - if (use_extrasnowlayers) then + if (use_firn_percolation_and_compaction) then do fc = 1, num_snowc c = filter_snowc(fc) burden(c) = 0._r8 @@ -637,7 +637,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & do j = -nlevsno+1, 0 do fc = 1, num_snowc c = filter_snowc(fc) - if (use_extrasnowlayers) then + if (use_firn_percolation_and_compaction) then t = col_pp%topounit(c) end if if (j >= snl(c)+1) then @@ -647,7 +647,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & ! If void is negative, then increase dz such that void = 0. ! This should be done for any landunit, but for now is done only for glacier_mec 1andunits. - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then ! I don't think the next 5 lines are necessary (removed in CLMv5) l = col_pp%landunit(c) if (ltype(l)==istice_mec .and. void < 0._r8) then @@ -665,7 +665,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & dexpf = exp(-c4*td) ! Settling as a result of destructive metamorphism - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then ddz1 = -c3*dexpf if (bi > dm) ddz1 = ddz1*exp(-46.0e-3_r8*(bi-dm)) else @@ -685,7 +685,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & if (h2osoi_liq(c,j) > 0.01_r8*dz(c,j)*frac_sno(c)) ddz1=ddz1*c5 ! Compaction due to overburden - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then ddz2 = -(burden(c)+wx/2._r8)*exp(-0.08_r8*td - c2*bi)/eta0 else p_gls = max(denice / bi, 1._r8) * grav * (burden(c) + wx/2._r8) @@ -724,7 +724,7 @@ subroutine SnowCompaction(bounds, num_snowc, filter_snowc, & ddz3 = 0._r8 end if - if (use_extrasnowlayers) then + if (use_firn_percolation_and_compaction) then ! Compaction occurring due to wind drift call WindDriftCompaction( & bi = bi, & diff --git a/components/elm/src/biogeophys/SnowSnicarMod.F90 b/components/elm/src/biogeophys/SnowSnicarMod.F90 index d6d6de75ffaa..de161c1d4c9d 100644 --- a/components/elm/src/biogeophys/SnowSnicarMod.F90 +++ b/components/elm/src/biogeophys/SnowSnicarMod.F90 @@ -11,7 +11,7 @@ module SnowSnicarMod use shr_kind_mod , only : r8 => shr_kind_r8 use shr_sys_mod , only : shr_sys_flush use shr_log_mod , only : errMsg => shr_log_errMsg - use elm_varctl , only : iulog, use_extrasnowlayers + use elm_varctl , only : iulog, use_firn_percolation_and_compaction use elm_varcon , only : namec use shr_const_mod , only : SHR_CONST_RHOICE use abortutils , only : endrun @@ -1435,13 +1435,13 @@ subroutine SnowAge_grain(bounds, & ! RE-FREEZING ! ! new snowfall [kg/m2] - if (do_capsnow(c_idx) .and. .not. use_extrasnowlayers) then + if (do_capsnow(c_idx) .and. .not. use_firn_percolation_and_compaction) then newsnow = max(0._r8, (qflx_snwcp_ice(c_idx)*dtime)) else newsnow = max(0._r8, (qflx_snow_grnd_col(c_idx)*dtime)) endif - if (use_extrasnowlayers) then + if (use_firn_percolation_and_compaction) then snw_rds_refrz = 1500._r8 else snw_rds_refrz = 1000._r8 diff --git a/components/elm/src/biogeophys/SoilFluxesMod.F90 b/components/elm/src/biogeophys/SoilFluxesMod.F90 index a83ae505ff7c..1b221a617a83 100644 --- a/components/elm/src/biogeophys/SoilFluxesMod.F90 +++ b/components/elm/src/biogeophys/SoilFluxesMod.F90 @@ -8,7 +8,7 @@ module SoilFluxesMod use shr_log_mod , only : errMsg => shr_log_errMsg use decompMod , only : bounds_type use abortutils , only : endrun - use elm_varctl , only : iulog, use_extrasnowlayers + use elm_varctl , only : iulog, use_firn_percolation_and_compaction use perfMod_GPU use elm_varpar , only : nlevsno, nlevgrnd, nlevurb, max_patch_per_col use atm2lndType , only : atm2lnd_type @@ -349,7 +349,7 @@ subroutine SoilFluxes (bounds, num_urbanl, filter_urbanl, & ! This was moved in from Hydrology2 to keep all pft-level ! calculations out of Hydrology2 - if (col_pp%snl(c) < 0 .and. do_capsnow(c) .and. .not. use_extrasnowlayers) then + if (col_pp%snl(c) < 0 .and. do_capsnow(c) .and. .not. use_firn_percolation_and_compaction) then qflx_snwcp_liq(p) = qflx_snwcp_liq(p)+frac_sno_eff(c)*qflx_dew_grnd(p) qflx_snwcp_ice(p) = qflx_snwcp_ice(p)+frac_sno_eff(c)*qflx_dew_snow(p) end if diff --git a/components/elm/src/main/controlMod.F90 b/components/elm/src/main/controlMod.F90 index 9db909a83fbb..dc7fa975dbb4 100755 --- a/components/elm/src/main/controlMod.F90 +++ b/components/elm/src/main/controlMod.F90 @@ -289,7 +289,7 @@ subroutine control_init( ) namelist /elm_inparm/ & use_nofire, use_lch4, use_vertsoilc, use_extralakelayers, & use_vichydro, use_century_decomp, use_cn, use_crop, use_snicar_frc, & - use_snicar_ad, use_extrasnowlayers,& + use_snicar_ad, use_firn_percolation_and_compaction, use_extrasnowlayers,& use_T_rho_dependent_snowthk, use_vancouver, use_mexicocity, use_noio ! cpl_bypass variables @@ -731,6 +731,7 @@ subroutine control_spmd() call mpi_bcast (use_vertsoilc, 1, MPI_LOGICAL, 0, mpicom, ier) call mpi_bcast (use_extralakelayers, 1, MPI_LOGICAL, 0, mpicom, ier) call mpi_bcast (use_extrasnowlayers, 1, MPI_LOGICAL, 0, mpicom, ier) + call mpi_bcast (use_firn_percolation_and_compaction, 1, MPI_LOGICAL, 0, mpicom, ier) call mpi_bcast (use_T_rho_dependent_snowthk, 1, MPI_LOGICAL, 0, mpicom, ier) call mpi_bcast (use_vichydro, 1, MPI_LOGICAL, 0, mpicom, ier) call mpi_bcast (use_century_decomp, 1, MPI_LOGICAL, 0, mpicom, ier) @@ -1036,6 +1037,7 @@ subroutine control_print () write(iulog,*) ' use_lake_wat_storage = ', use_lake_wat_storage write(iulog,*) ' use_extralakelayers = ', use_extralakelayers write(iulog,*) ' use_extrasnowlayers = ', use_extrasnowlayers + write(iulog,*) ' use_firn_percolation_and_compaction = ', use_firn_percolation_and_compaction write(iulog,*) ' use_T_rho_dependent_snowthk = ', use_T_rho_dependent_snowthk write(iulog,*) ' use_vichydro = ', use_vichydro write(iulog,*) ' use_century_decomp = ', use_century_decomp diff --git a/components/elm/src/main/elm_driver.F90 b/components/elm/src/main/elm_driver.F90 index 06668f6cad07..c4169c6b04ea 100644 --- a/components/elm/src/main/elm_driver.F90 +++ b/components/elm/src/main/elm_driver.F90 @@ -12,7 +12,7 @@ module elm_driver use shr_sys_mod , only : shr_sys_flush use shr_log_mod , only : errMsg => shr_log_errMsg use elm_varpar , only : nlevtrc_soil, nlevsoi - use elm_varctl , only : wrtdia, iulog, create_glacier_mec_landunit, use_fates, use_betr + use elm_varctl , only : wrtdia, iulog, create_glacier_mec_landunit, use_fates, use_betr, use_firn_percolation_and_compaction use elm_varctl , only : use_cn, use_lch4, use_voc, use_noio, use_c13, use_c14 use elm_varctl , only : use_erosion, use_fates_sp, use_fan,use_extrasnowlayers use elm_varctl , only : mpi_sync_nstep_freq @@ -1622,7 +1622,7 @@ subroutine elm_drv_init(bounds, & ! Save snow mass at previous time step h2osno_old(c) = h2osno(c) - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then ! Decide whether to cap snow if (h2osno(c) > h2osno_max) then do_capsnow(c) = .true. @@ -1790,7 +1790,7 @@ subroutine elm_drv_patch2col (bounds, num_nolakec, filter_nolakec, & qflx_snow_grnd_patch(bounds%begp:bounds%endp), & qflx_snow_grnd_col (bounds%begc:bounds%endc)) - if (.not. use_extrasnowlayers) then + if (.not. use_firn_percolation_and_compaction) then call p2c (bounds, num_allc, filter_allc, & veg_wf%qflx_snwcp_liq(bounds%begp:bounds%endp), & col_wf%qflx_snwcp_liq(bounds%begc:bounds%endc)) diff --git a/components/elm/src/main/elm_varctl.F90 b/components/elm/src/main/elm_varctl.F90 index f1a1733391fd..fb9ca0dbd24e 100644 --- a/components/elm/src/main/elm_varctl.F90 +++ b/components/elm/src/main/elm_varctl.F90 @@ -375,6 +375,7 @@ module elm_varctl logical, public :: use_snicar_frc = .false. logical, public :: use_snicar_ad = .false. logical, public :: use_extrasnowlayers = .false. + logical, public :: use_firn_percolation_and_compaction = .false. logical, public :: use_vancouver = .false. logical, public :: use_mexicocity = .false. logical, public :: use_noio = .false. From 8e0b8996803df94c872c7400e99d59c9dd920899 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 7 Jan 2025 12:40:04 -0800 Subject: [PATCH 12/39] Revert "minor mods to remove use_firn... for use_extrasnowlayers" This reverts commit ee58bb0e9d690e016e3906f91a90da68145c11ee. --- components/elm/src/biogeophys/CanopyHydrologyMod.F90 | 2 +- components/elm/src/main/elm_driver.F90 | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/components/elm/src/biogeophys/CanopyHydrologyMod.F90 b/components/elm/src/biogeophys/CanopyHydrologyMod.F90 index 5838fd448002..5381b1358f8e 100755 --- a/components/elm/src/biogeophys/CanopyHydrologyMod.F90 +++ b/components/elm/src/biogeophys/CanopyHydrologyMod.F90 @@ -116,7 +116,7 @@ subroutine CanopyHydrology(bounds, & use elm_varcon , only : hfus, denice, zlnd, rpi, spval, tfrz use column_varcon , only : icol_roof, icol_sunwall, icol_shadewall use landunit_varcon , only : istcrop, istice, istwet, istsoil, istice_mec, istdlak - use elm_varctl , only : subgridflag, use_extrasnowlayers + use elm_varctl , only : subgridflag use elm_varpar , only : nlevsoi,nlevsno use elm_varsur , only : wt_lunit, f_grd, f_surf use atm2lndType , only : atm2lnd_type diff --git a/components/elm/src/main/elm_driver.F90 b/components/elm/src/main/elm_driver.F90 index c4169c6b04ea..51f08bf235c3 100644 --- a/components/elm/src/main/elm_driver.F90 +++ b/components/elm/src/main/elm_driver.F90 @@ -14,7 +14,7 @@ module elm_driver use elm_varpar , only : nlevtrc_soil, nlevsoi use elm_varctl , only : wrtdia, iulog, create_glacier_mec_landunit, use_fates, use_betr, use_firn_percolation_and_compaction use elm_varctl , only : use_cn, use_lch4, use_voc, use_noio, use_c13, use_c14 - use elm_varctl , only : use_erosion, use_fates_sp, use_fan,use_extrasnowlayers + use elm_varctl , only : use_erosion, use_fates_sp, use_fan use elm_varctl , only : mpi_sync_nstep_freq use elm_time_manager , only : get_step_size, get_curr_date, get_ref_date, get_nstep, is_beg_curr_day, get_curr_time_string use elm_time_manager , only : get_curr_calday, get_days_per_year From 8ad01ed58270504146fea9f46ad0bd642f8bd182 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 7 Jan 2025 12:48:46 -0800 Subject: [PATCH 13/39] water balance error testing --- components/elm/src/biogeophys/BalanceCheckMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index 3ce501475dc5..c4e70d5073a1 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -320,7 +320,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & errh2o(c) = endwb(c) - begwb(c) & - (forc_rain_col(c) + forc_snow_col(c) + qflx_floodc(c) + qflx_surf_irrig_col(c) + qflx_over_supply_col(c) & - qflx_evap_tot(c) - qflx_surf(c) - qflx_h2osfc_surf(c) & - - qflx_qrgwl(c) - qflx_drain(c) - qflx_drain_perched(c) - qflx_snwcp_ice(c) - qflx_ice_runoff_xs(c) & + - qflx_qrgwl(c) - qflx_drain(c) - qflx_drain_perched(c) - qflx_snwcp_ice(c) & - qflx_lateral(c) + qflx_h2orof_drain(c)) * dtime dwb(c) = (endwb(c)-begwb(c))/dtime From c354b3edcd204cc12b27483d2d900fc6ef5679b1 Mon Sep 17 00:00:00 2001 From: Chloe Date: Wed, 8 Jan 2025 12:33:41 -0800 Subject: [PATCH 14/39] added Steve P's G +GLC 20thTR config --- cime_config/allactive/config_compsets.xml | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/cime_config/allactive/config_compsets.xml b/cime_config/allactive/config_compsets.xml index d5a52a8bac38..f474f4932ebc 100755 --- a/cime_config/allactive/config_compsets.xml +++ b/cime_config/allactive/config_compsets.xml @@ -459,6 +459,10 @@ 1850_EAM%CMIP6_ELM%SPBC_MPASSI_MPASO_MOSART_MALI_SWAV + + BGWCYCL20TR + 20TRSOI_EAM%CMIP6_ELM%CNPRDCTCBCTOP_MPASSI_MPASO_MOSART_MALI_SWAV + From 0dea30f8b3ba323140955e22a3fb3775ac3eb104 Mon Sep 17 00:00:00 2001 From: Chloe Date: Wed, 8 Jan 2025 13:31:32 -0800 Subject: [PATCH 15/39] Revert "water balance error testing" This reverts commit 8ad01ed58270504146fea9f46ad0bd642f8bd182. --- components/elm/src/biogeophys/BalanceCheckMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index c4e70d5073a1..3ce501475dc5 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -320,7 +320,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & errh2o(c) = endwb(c) - begwb(c) & - (forc_rain_col(c) + forc_snow_col(c) + qflx_floodc(c) + qflx_surf_irrig_col(c) + qflx_over_supply_col(c) & - qflx_evap_tot(c) - qflx_surf(c) - qflx_h2osfc_surf(c) & - - qflx_qrgwl(c) - qflx_drain(c) - qflx_drain_perched(c) - qflx_snwcp_ice(c) & + - qflx_qrgwl(c) - qflx_drain(c) - qflx_drain_perched(c) - qflx_snwcp_ice(c) - qflx_ice_runoff_xs(c) & - qflx_lateral(c) + qflx_h2orof_drain(c)) * dtime dwb(c) = (endwb(c)-begwb(c))/dtime From 028add438a3b1a65224dff25ae30fd4e08ee6387 Mon Sep 17 00:00:00 2001 From: Chloe Date: Wed, 8 Jan 2025 13:34:05 -0800 Subject: [PATCH 16/39] water balance error testing --- components/elm/src/biogeophys/SoilHydrologyMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/components/elm/src/biogeophys/SoilHydrologyMod.F90 b/components/elm/src/biogeophys/SoilHydrologyMod.F90 index 01a13bde3ba4..30d45d6f119f 100644 --- a/components/elm/src/biogeophys/SoilHydrologyMod.F90 +++ b/components/elm/src/biogeophys/SoilHydrologyMod.F90 @@ -1530,7 +1530,7 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte ! add in ice check xs1(c) = max(max(h2osoi_ice(c,1),0._r8)-max(0._r8,(pondmx+watsat(c,1)*dzmm(c,1)-h2osoi_liq(c,1))),0._r8) h2osoi_ice(c,1) = min(max(0._r8,pondmx+watsat(c,1)*dzmm(c,1)-h2osoi_liq(c,1)), h2osoi_ice(c,1)) - if (lun_pp%itype(col_pp%landunit(c)) == istice .or. lun_pp%itype(col_pp%landunit(c)) == istice_mec) then + if ( (lun_pp%itype(col_pp%landunit(c)) == istice .or. lun_pp%itype(col_pp%landunit(c)) == istice_mec) .and. (.not. use_extrasnowlayers)) then qflx_snwcp_ice(c) = qflx_snwcp_ice(c) + xs1(c) / dtime else qflx_ice_runoff_xs(c) = qflx_ice_runoff_xs(c) + xs1(c) / dtime From d89b1f7d34e1476c85345d9a95ee81e6186bf892 Mon Sep 17 00:00:00 2001 From: Chloe Date: Thu, 9 Jan 2025 08:30:10 -0800 Subject: [PATCH 17/39] water balance error testing --- components/elm/src/biogeophys/SoilHydrologyMod.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/components/elm/src/biogeophys/SoilHydrologyMod.F90 b/components/elm/src/biogeophys/SoilHydrologyMod.F90 index 30d45d6f119f..671333192ea7 100644 --- a/components/elm/src/biogeophys/SoilHydrologyMod.F90 +++ b/components/elm/src/biogeophys/SoilHydrologyMod.F90 @@ -984,7 +984,7 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte use elm_varpar , only : nlevsoi, nlevgrnd, nlayer, nlayert use elm_varcon , only : pondmx, tfrz, watmin,rpi, secspday, nlvic use column_varcon , only : icol_roof, icol_road_imperv, icol_road_perv - use elm_varctl , only : use_vsfm, use_var_soil_thick + use elm_varctl , only : use_vsfm, use_var_soil_thick, use_firn_percolation_and_compaction use SoilWaterMovementMod, only : zengdecker_2009_with_var_soil_thick use pftvarcon , only : rsub_top_globalmax use LandunitType , only : lun_pp @@ -1530,7 +1530,7 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte ! add in ice check xs1(c) = max(max(h2osoi_ice(c,1),0._r8)-max(0._r8,(pondmx+watsat(c,1)*dzmm(c,1)-h2osoi_liq(c,1))),0._r8) h2osoi_ice(c,1) = min(max(0._r8,pondmx+watsat(c,1)*dzmm(c,1)-h2osoi_liq(c,1)), h2osoi_ice(c,1)) - if ( (lun_pp%itype(col_pp%landunit(c)) == istice .or. lun_pp%itype(col_pp%landunit(c)) == istice_mec) .and. (.not. use_extrasnowlayers)) then + if ( (lun_pp%itype(col_pp%landunit(c)) == istice .or. lun_pp%itype(col_pp%landunit(c)) == istice_mec) .and. (.not. use_firn_percolation_and_compaction)) then qflx_snwcp_ice(c) = qflx_snwcp_ice(c) + xs1(c) / dtime else qflx_ice_runoff_xs(c) = qflx_ice_runoff_xs(c) + xs1(c) / dtime From 87521c4eb6e8ddb522085f982b203a6652dd1608 Mon Sep 17 00:00:00 2001 From: Chloe Date: Thu, 9 Jan 2025 11:32:53 -0800 Subject: [PATCH 18/39] water balance debugging --- components/elm/src/biogeophys/SoilHydrologyMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/components/elm/src/biogeophys/SoilHydrologyMod.F90 b/components/elm/src/biogeophys/SoilHydrologyMod.F90 index 671333192ea7..3064e8045d66 100644 --- a/components/elm/src/biogeophys/SoilHydrologyMod.F90 +++ b/components/elm/src/biogeophys/SoilHydrologyMod.F90 @@ -1530,7 +1530,7 @@ subroutine Drainage(bounds, num_hydrologyc, filter_hydrologyc, num_urbanc, filte ! add in ice check xs1(c) = max(max(h2osoi_ice(c,1),0._r8)-max(0._r8,(pondmx+watsat(c,1)*dzmm(c,1)-h2osoi_liq(c,1))),0._r8) h2osoi_ice(c,1) = min(max(0._r8,pondmx+watsat(c,1)*dzmm(c,1)-h2osoi_liq(c,1)), h2osoi_ice(c,1)) - if ( (lun_pp%itype(col_pp%landunit(c)) == istice .or. lun_pp%itype(col_pp%landunit(c)) == istice_mec) .and. (.not. use_firn_percolation_and_compaction)) then + if ( (lun_pp%itype(col_pp%landunit(c)) == istice .or. lun_pp%itype(col_pp%landunit(c)) == istice_mec) .or. (.not. use_firn_percolation_and_compaction)) then qflx_snwcp_ice(c) = qflx_snwcp_ice(c) + xs1(c) / dtime else qflx_ice_runoff_xs(c) = qflx_ice_runoff_xs(c) + xs1(c) / dtime From 52fdd1bbfa5ea7798408ab572287e28ba4b0a7b3 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 14 Jan 2025 11:02:38 -0800 Subject: [PATCH 19/39] minor clean up for PR --- components/elm/src/biogeophys/BalanceCheckMod.F90 | 12 ++++++------ components/elm/src/biogeophys/SoilTemperatureMod.F90 | 1 - 2 files changed, 6 insertions(+), 7 deletions(-) diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index 3ce501475dc5..209a9977dd01 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -465,13 +465,13 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & + qflx_snow_melt(c) + qflx_sl_top_soil(c) endif else ! firn model - snow_sources(c) = qflx_snow_grnd_col(c) & - + frac_sno_eff(c) * (qflx_rain_grnd_col(c) & - + qflx_dew_snow(c) + qflx_dew_grnd(c) ) + snow_sources(c) = qflx_snow_grnd_col(c) & + + frac_sno_eff(c) * (qflx_rain_grnd_col(c) & + + qflx_dew_snow(c) + qflx_dew_grnd(c) ) - snow_sinks(c) = frac_sno_eff(c) * (qflx_sub_snow(c) + qflx_evap_grnd(c) ) & - + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & - + qflx_snow_melt(c) + qflx_sl_top_soil(c) + snow_sinks(c) = frac_sno_eff(c) * (qflx_sub_snow(c) + qflx_evap_grnd(c) ) & + + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & + + qflx_snow_melt(c) + qflx_sl_top_soil(c) endif endif diff --git a/components/elm/src/biogeophys/SoilTemperatureMod.F90 b/components/elm/src/biogeophys/SoilTemperatureMod.F90 index 7420e857bce7..bef47d47ff69 100644 --- a/components/elm/src/biogeophys/SoilTemperatureMod.F90 +++ b/components/elm/src/biogeophys/SoilTemperatureMod.F90 @@ -1567,7 +1567,6 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & qflx_snomelt(c) = max(0._r8,(temp1-h2osno(c)))/dtime ! kg/(m2 s) xmf(c) = hfus*qflx_snomelt(c) qflx_snow_melt(c) = qflx_snomelt(c) - !qflx_snomelt_lyr(c,j) = qflx_snomelt(c) endif endif From d807b9a98a728dcfcae7c9ab927333f43ac567f4 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 14 Jan 2025 11:16:48 -0800 Subject: [PATCH 20/39] minor changes for PR --- .../elm/src/biogeophys/BalanceCheckMod.F90 | 24 +++++++++---------- 1 file changed, 12 insertions(+), 12 deletions(-) diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index 209a9977dd01..a56cc57fc45a 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -465,13 +465,13 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & + qflx_snow_melt(c) + qflx_sl_top_soil(c) endif else ! firn model - snow_sources(c) = qflx_snow_grnd_col(c) & - + frac_sno_eff(c) * (qflx_rain_grnd_col(c) & - + qflx_dew_snow(c) + qflx_dew_grnd(c) ) + snow_sources(c) = qflx_snow_grnd_col(c) & + + frac_sno_eff(c) * (qflx_rain_grnd_col(c) & + + qflx_dew_snow(c) + qflx_dew_grnd(c) ) - snow_sinks(c) = frac_sno_eff(c) * (qflx_sub_snow(c) + qflx_evap_grnd(c) ) & - + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & - + qflx_snow_melt(c) + qflx_sl_top_soil(c) + snow_sinks(c) = frac_sno_eff(c) * (qflx_sub_snow(c) + qflx_evap_grnd(c) ) & + + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & + + qflx_snow_melt(c) + qflx_sl_top_soil(c) endif endif @@ -493,13 +493,13 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & + qflx_snow_melt(c) + qflx_sl_top_soil(c) endif else ! firn model - snow_sources(c) = (qflx_snow_grnd_col(c) - qflx_snow_h2osfc(c) ) & - + frac_sno_eff(c) * (qflx_rain_grnd_col(c) & - + qflx_dew_snow(c) + qflx_dew_grnd(c) ) + qflx_h2osfc_to_ice(c) + snow_sources(c) = (qflx_snow_grnd_col(c) - qflx_snow_h2osfc(c) ) & + + frac_sno_eff(c) * (qflx_rain_grnd_col(c) & + + qflx_dew_snow(c) + qflx_dew_grnd(c) ) + qflx_h2osfc_to_ice(c) - snow_sinks(c) = frac_sno_eff(c) * (qflx_sub_snow(c) + qflx_evap_grnd(c)) & - + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & - + qflx_snow_melt(c) + qflx_sl_top_soil(c) + snow_sinks(c) = frac_sno_eff(c) * (qflx_sub_snow(c) + qflx_evap_grnd(c)) & + + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & + + qflx_snow_melt(c) + qflx_sl_top_soil(c) endif endif if (glc_dyn_runoff_routing(g)) then From 8bb5682263eeea8a1c3a97d63556b8ce2576273b Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 14 Jan 2025 11:19:16 -0800 Subject: [PATCH 21/39] minor changes for PR --- components/elm/src/biogeophys/BalanceCheckMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index a56cc57fc45a..861fd88e5b3f 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -536,7 +536,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & ' col_pp%itype= ',col_pp%itype(indexc), & ' lun_pp%itype= ',lun_pp%itype(col_pp%landunit(indexc)), & ' errh2osno= ',errh2osno(indexc) - + if (abs(errh2osno(indexc)) > 1.e-4_r8 .and. (nstep > 2) ) then write(iulog,*)'elm model is stopping - error is greater than 1e-4 (mm)' write(iulog,*)'nstep = ',nstep From 13693ef8e621949d744728a8acdde979512b6f96 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 14 Jan 2025 11:21:50 -0800 Subject: [PATCH 22/39] minor changes for PR --- components/elm/src/biogeophys/BalanceCheckMod.F90 | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index 861fd88e5b3f..e60cd0c7ed39 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -495,13 +495,14 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & else ! firn model snow_sources(c) = (qflx_snow_grnd_col(c) - qflx_snow_h2osfc(c) ) & + frac_sno_eff(c) * (qflx_rain_grnd_col(c) & - + qflx_dew_snow(c) + qflx_dew_grnd(c) ) + qflx_h2osfc_to_ice(c) + + qflx_dew_snow(c) + qflx_dew_grnd(c) ) + qflx_h2osfc_to_ice(c) snow_sinks(c) = frac_sno_eff(c) * (qflx_sub_snow(c) + qflx_evap_grnd(c)) & + qflx_snwcp_ice(c) + qflx_snwcp_liq(c) & + qflx_snow_melt(c) + qflx_sl_top_soil(c) endif endif + if (glc_dyn_runoff_routing(g)) then ! Need to add qflx_glcice_frz to snow_sinks for the same reason as it is ! added to errh2o above - see the comment above for details. From 6426ef69a84c93a2133dc22c7425f108cb04924b Mon Sep 17 00:00:00 2001 From: Chloe Date: Wed, 18 Sep 2024 10:03:19 -0700 Subject: [PATCH 23/39] remove conditional if (create_glacier_mec_landunit) then to include ELM/GLC SMB QICE fields - makes QICE fields in h0 by default --- .../elm/src/data_types/ColumnDataType.F90 | 30 +++++++++---------- 1 file changed, 14 insertions(+), 16 deletions(-) diff --git a/components/elm/src/data_types/ColumnDataType.F90 b/components/elm/src/data_types/ColumnDataType.F90 index 0329430fbbd2..1de17aafd583 100644 --- a/components/elm/src/data_types/ColumnDataType.F90 +++ b/components/elm/src/data_types/ColumnDataType.F90 @@ -5857,22 +5857,20 @@ subroutine col_wf_init(this, begc, endc) avgflag='A', long_name='column-integrated snow freezing rate', & ptr_col=this%qflx_snofrz, set_lake=spval, c2l_scale_type='urbanf') - if (create_glacier_mec_landunit) then - this%qflx_glcice(begc:endc) = spval - call hist_addfld1d (fname='QICE', units='mm/s', & - avgflag='A', long_name='ice growth/melt', & - ptr_col=this%qflx_glcice, l2g_scale_type='ice') - - this%qflx_glcice_frz(begc:endc) = spval - call hist_addfld1d (fname='QICE_FRZ', units='mm/s', & - avgflag='A', long_name='ice growth', & - ptr_col=this%qflx_glcice_frz, l2g_scale_type='ice') - - this%qflx_glcice_melt(begc:endc) = spval - call hist_addfld1d (fname='QICE_MELT', units='mm/s', & - avgflag='A', long_name='ice melt', & - ptr_col=this%qflx_glcice_melt, l2g_scale_type='ice') - endif + this%qflx_glcice(begc:endc) = spval + call hist_addfld1d (fname='QICE', units='mm/s', & + avgflag='A', long_name='ice growth/melt', & + ptr_col=this%qflx_glcice, l2g_scale_type='ice') + + this%qflx_glcice_frz(begc:endc) = spval + call hist_addfld1d (fname='QICE_FRZ', units='mm/s', & + avgflag='A', long_name='ice growth', & + ptr_col=this%qflx_glcice_frz, l2g_scale_type='ice') + + this%qflx_glcice_melt(begc:endc) = spval + call hist_addfld1d (fname='QICE_MELT', units='mm/s', & + avgflag='A', long_name='ice melt', & + ptr_col=this%qflx_glcice_melt, l2g_scale_type='ice') ! As defined here, snow_sources - snow_sinks will equal the change in h2osno at any ! given time step but only if there is at least one snow layer (for all landunits From 0694b85b4e93b2e39d257051fa2f9de0948fce4d Mon Sep 17 00:00:00 2001 From: Chloe Date: Wed, 2 Oct 2024 13:11:36 -0700 Subject: [PATCH 24/39] mods to QICE field calcs, added conditional to include land types isice --- components/elm/src/biogeophys/HydrologyDrainageMod.F90 | 2 +- components/elm/src/biogeophys/SoilTemperatureMod.F90 | 3 +-- 2 files changed, 2 insertions(+), 3 deletions(-) diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index da4006e9771c..1c083c5bfa12 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -222,7 +222,7 @@ subroutine HydrologyDrainage(bounds, & g = col_pp%gridcell(c) ! In the following, we convert glc_snow_persistence_max_days to r8 to avoid overflow if ( (snow_persistence(c) >= (real(glc_snow_persistence_max_days, r8) * secspday)) & - .or. lun_pp%itype(l) == istice_mec) then + .or. lun_pp%itype(l) == istice_mec .or. lun_pp%itype(l) == istice) then qflx_glcice_frz(c) = qflx_snwcp_ice(c) qflx_glcice(c) = qflx_glcice(c) + qflx_glcice_frz(c) if (glc_dyn_runoff_routing(g)) qflx_snwcp_ice(c) = 0._r8 diff --git a/components/elm/src/biogeophys/SoilTemperatureMod.F90 b/components/elm/src/biogeophys/SoilTemperatureMod.F90 index bef47d47ff69..d03aab75ceb1 100644 --- a/components/elm/src/biogeophys/SoilTemperatureMod.F90 +++ b/components/elm/src/biogeophys/SoilTemperatureMod.F90 @@ -1646,8 +1646,7 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & ! as computed in HydrologyDrainageMod.F90. l = col_pp%landunit(c) - if (lun_pp%itype(l)==istice_mec) then - + if ( (lun_pp%itype(l)==istice) .or. (lun_pp%itype(l)==istice_mec) ) then if (j>=1 .and. h2osoi_liq(c,j) > 0._r8) then ! ice layer with meltwater ! melting corresponds to a negative ice flux qflx_glcice_melt(c) = qflx_glcice_melt(c) + h2osoi_liq(c,j)/dtime From f18ae84a0c186ddb4c01782ea7e324d095ed7b78 Mon Sep 17 00:00:00 2001 From: Chloe Date: Wed, 2 Oct 2024 14:49:34 -0700 Subject: [PATCH 25/39] QICE h0 chnages -- not b4b --- components/elm/src/biogeophys/HydrologyDrainageMod.F90 | 2 +- components/elm/src/biogeophys/SoilTemperatureMod.F90 | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index 1c083c5bfa12..cc22bb29d53e 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -47,7 +47,7 @@ subroutine HydrologyDrainage(bounds, & ! ! !USES: !$acc routine seq - use landunit_varcon , only : istice, istwet, istsoil, istice_mec, istcrop + use landunit_varcon , only : istice, istwet, istsoil, istice_mec, istcrop, istice use column_varcon , only : icol_roof, icol_road_imperv, icol_road_perv, icol_sunwall, icol_shadewall use elm_varcon , only : denh2o, denice, secspday use elm_varctl , only : glc_snow_persistence_max_days, use_vichydro, use_betr diff --git a/components/elm/src/biogeophys/SoilTemperatureMod.F90 b/components/elm/src/biogeophys/SoilTemperatureMod.F90 index d03aab75ceb1..9118a9439c2b 100644 --- a/components/elm/src/biogeophys/SoilTemperatureMod.F90 +++ b/components/elm/src/biogeophys/SoilTemperatureMod.F90 @@ -1316,7 +1316,7 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & use elm_varctl , only : iulog use elm_varcon , only : tfrz, hfus, grav use column_varcon , only : icol_roof, icol_sunwall, icol_shadewall, icol_road_perv - use landunit_varcon , only : istsoil, istcrop, istice_mec + use landunit_varcon , only : istsoil, istcrop, istice_mec,istice ! ! !ARGUMENTS: type(bounds_type) , intent(in) :: bounds From 725b410926129b048779158f78c329302b606abe Mon Sep 17 00:00:00 2001 From: Chloe Date: Thu, 3 Oct 2024 13:15:40 -0700 Subject: [PATCH 26/39] added conditionals to calc SMB (QICE vars) even for non MEC runs - these changes should be b4b and just add additional SMB output of what would be passed to mali if we had an active GLC --- .../elm/src/biogeophys/HydrologyDrainageMod.F90 | 8 +++++++- components/elm/src/biogeophys/SoilTemperatureMod.F90 | 12 +++++++++++- 2 files changed, 18 insertions(+), 2 deletions(-) diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index cc22bb29d53e..ac748febaf47 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -222,11 +222,17 @@ subroutine HydrologyDrainage(bounds, & g = col_pp%gridcell(c) ! In the following, we convert glc_snow_persistence_max_days to r8 to avoid overflow if ( (snow_persistence(c) >= (real(glc_snow_persistence_max_days, r8) * secspday)) & - .or. lun_pp%itype(l) == istice_mec .or. lun_pp%itype(l) == istice) then + .or. lun_pp%itype(l) == istice_mec) then qflx_glcice_frz(c) = qflx_snwcp_ice(c) qflx_glcice(c) = qflx_glcice(c) + qflx_glcice_frz(c) if (glc_dyn_runoff_routing(g)) qflx_snwcp_ice(c) = 0._r8 end if + + if (lun_pp%itype(l) == istice) then + qflx_glcice_frz(c) = qflx_snwcp_ice(c) + qflx_glcice(c) = qflx_glcice(c) + qflx_glcice_frz(c) + end if + end do ! Determine wetland and land ice hydrology (must be placed here diff --git a/components/elm/src/biogeophys/SoilTemperatureMod.F90 b/components/elm/src/biogeophys/SoilTemperatureMod.F90 index 9118a9439c2b..be525add7230 100644 --- a/components/elm/src/biogeophys/SoilTemperatureMod.F90 +++ b/components/elm/src/biogeophys/SoilTemperatureMod.F90 @@ -1646,7 +1646,7 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & ! as computed in HydrologyDrainageMod.F90. l = col_pp%landunit(c) - if ( (lun_pp%itype(l)==istice) .or. (lun_pp%itype(l)==istice_mec) ) then + if ( lun_pp%itype(l)==istice_mec) then if (j>=1 .and. h2osoi_liq(c,j) > 0._r8) then ! ice layer with meltwater ! melting corresponds to a negative ice flux qflx_glcice_melt(c) = qflx_glcice_melt(c) + h2osoi_liq(c,j)/dtime @@ -1658,6 +1658,16 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & endif ! liquid water is present endif ! istice_mec + ! for diagnostic QICE SMB output only - + ! these are to calculate SMB even without MECs + if ( lun_pp%itype(l)==istice) then + if (j>=1 .and. h2osoi_liq(c,j) > 0._r8) then ! ice layer with meltwater + ! melting corresponds to a negative ice flux + qflx_glcice_melt(c) = qflx_glcice_melt(c) + h2osoi_liq(c,j)/dtime + qflx_glcice(c) = qflx_glcice(c) - h2osoi_liq(c,j)/dtime + endif ! liquid water is present + endif ! istice_mec + end do ! end of column-loop enddo ! end of level-loop From 68767736c698ad411d55af1f45a7455024104934 Mon Sep 17 00:00:00 2001 From: Chloe Date: Fri, 4 Oct 2024 08:02:57 -0700 Subject: [PATCH 27/39] modified conditional for outputting diagnostic QICE fields --- .../src/biogeophys/HydrologyDrainageMod.F90 | 19 ++++++++++--------- 1 file changed, 10 insertions(+), 9 deletions(-) diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index ac748febaf47..97233433b4d8 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -222,15 +222,16 @@ subroutine HydrologyDrainage(bounds, & g = col_pp%gridcell(c) ! In the following, we convert glc_snow_persistence_max_days to r8 to avoid overflow if ( (snow_persistence(c) >= (real(glc_snow_persistence_max_days, r8) * secspday)) & - .or. lun_pp%itype(l) == istice_mec) then - qflx_glcice_frz(c) = qflx_snwcp_ice(c) - qflx_glcice(c) = qflx_glcice(c) + qflx_glcice_frz(c) - if (glc_dyn_runoff_routing(g)) qflx_snwcp_ice(c) = 0._r8 - end if - - if (lun_pp%itype(l) == istice) then - qflx_glcice_frz(c) = qflx_snwcp_ice(c) - qflx_glcice(c) = qflx_glcice(c) + qflx_glcice_frz(c) + .or. lun_pp%itype(l) == istice_mec .or. lun_pp%itype(l) == istice) then + + if (lun_pp%itype(l) == istice) then + qflx_glcice_frz(c) = qflx_snwcp_ice(c) + qflx_glcice(c) = qflx_glcice(c) + qflx_glcice_frz(c) + else + qflx_glcice_frz(c) = qflx_snwcp_ice(c) + qflx_glcice(c) = qflx_glcice(c) + qflx_glcice_frz(c) + if (glc_dyn_runoff_routing(g)) qflx_snwcp_ice(c) = 0._r8 + end if ! lun_pp%itype(l) == istice end if end do From 51ece732359a2aaf54eb9bf049ec80371a976414 Mon Sep 17 00:00:00 2001 From: Chloe Date: Fri, 4 Oct 2024 09:07:08 -0700 Subject: [PATCH 28/39] mods to QICE conditions to get b4b --- components/elm/src/biogeophys/HydrologyDrainageMod.F90 | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index 97233433b4d8..2c11de226d7c 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -224,13 +224,14 @@ subroutine HydrologyDrainage(bounds, & if ( (snow_persistence(c) >= (real(glc_snow_persistence_max_days, r8) * secspday)) & .or. lun_pp%itype(l) == istice_mec .or. lun_pp%itype(l) == istice) then - if (lun_pp%itype(l) == istice) then + if ( (snow_persistence(c) >= (real(glc_snow_persistence_max_days, r8) * secspday)) & + .or. lun_pp%itype(l) == istice_mec) then qflx_glcice_frz(c) = qflx_snwcp_ice(c) qflx_glcice(c) = qflx_glcice(c) + qflx_glcice_frz(c) - else + if (glc_dyn_runoff_routing(g)) qflx_snwcp_ice(c) = 0._r8 + else qflx_glcice_frz(c) = qflx_snwcp_ice(c) qflx_glcice(c) = qflx_glcice(c) + qflx_glcice_frz(c) - if (glc_dyn_runoff_routing(g)) qflx_snwcp_ice(c) = 0._r8 end if ! lun_pp%itype(l) == istice end if From eb8b4139c04e9648aa45060d7d8a7d4bca43b864 Mon Sep 17 00:00:00 2001 From: Chloe Date: Fri, 4 Oct 2024 10:51:49 -0700 Subject: [PATCH 29/39] QICE calc mods --- components/elm/src/biogeophys/BalanceCheckMod.F90 | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index e60cd0c7ed39..8d59bd9e5fe6 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -506,7 +506,11 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & if (glc_dyn_runoff_routing(g)) then ! Need to add qflx_glcice_frz to snow_sinks for the same reason as it is ! added to errh2o above - see the comment above for details. - snow_sinks(c) = snow_sinks(c) + qflx_glcice_frz(c) + if (lun_pp%itype(l) == istice) then + snow_sinks(c) = snow_sinks(c) + else + snow_sinks(c) = snow_sinks(c) + qflx_glcice_frz(c) + end if end if errh2osno(c) = (h2osno(c) - h2osno_old(c)) - (snow_sources(c) - snow_sinks(c)) * dtime From f5195b8c4f620a258ef1b5845ce8cfa13eee307a Mon Sep 17 00:00:00 2001 From: Chloe Date: Fri, 4 Oct 2024 13:17:06 -0700 Subject: [PATCH 30/39] mods to QICE scheme, if not GLC/MEC output diagnostic QICE fields that do not change h2o mass balance --- .../elm/src/biogeophys/BalanceCheckMod.F90 | 8 +-- .../src/biogeophys/HydrologyDrainageMod.F90 | 25 +++++---- .../elm/src/biogeophys/SoilTemperatureMod.F90 | 6 ++- .../elm/src/data_types/ColumnDataType.F90 | 52 ++++++++++++++----- 4 files changed, 56 insertions(+), 35 deletions(-) diff --git a/components/elm/src/biogeophys/BalanceCheckMod.F90 b/components/elm/src/biogeophys/BalanceCheckMod.F90 index 8d59bd9e5fe6..0781321006ee 100755 --- a/components/elm/src/biogeophys/BalanceCheckMod.F90 +++ b/components/elm/src/biogeophys/BalanceCheckMod.F90 @@ -167,7 +167,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & use elm_varcon , only : spval, h2osno_max use column_varcon , only : icol_roof, icol_sunwall, icol_shadewall use column_varcon , only : icol_road_perv, icol_road_imperv - use landunit_varcon , only : istice_mec, istdlak, istsoil,istcrop,istwet + use landunit_varcon , only : istice_mec, istice, istdlak, istsoil,istcrop,istwet use elm_varctl , only : create_glacier_mec_landunit use elm_initializeMod , only : surfalb_vars use CanopyStateType , only : canopystate_type @@ -506,11 +506,7 @@ subroutine ColWaterBalanceCheck( bounds, num_do_smb_c, filter_do_smb_c, & if (glc_dyn_runoff_routing(g)) then ! Need to add qflx_glcice_frz to snow_sinks for the same reason as it is ! added to errh2o above - see the comment above for details. - if (lun_pp%itype(l) == istice) then - snow_sinks(c) = snow_sinks(c) - else - snow_sinks(c) = snow_sinks(c) + qflx_glcice_frz(c) - end if + snow_sinks(c) = snow_sinks(c) + qflx_glcice_frz(c) end if errh2osno(c) = (h2osno(c) - h2osno_old(c)) - (snow_sources(c) - snow_sinks(c)) * dtime diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index 2c11de226d7c..1c40981ad4d1 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -120,8 +120,9 @@ subroutine HydrologyDrainage(bounds, & qflx_runoff_r => col_wf%qflx_runoff_r , & ! Output: [real(r8) (:) ] Rural total runoff (qflx_drain+qflx_surf+qflx_qrgwl) (mm H2O /s) qflx_snwcp_ice => col_wf%qflx_snwcp_ice , & ! Output: [real(r8) (:) ] excess snowfall due to snow capping (mm H2O /s) [+]` qflx_glcice => col_wf%qflx_glcice , & ! Output: [real(r8) (:) ] flux of new glacier ice (mm H2O /s) - qflx_glcice_frz => col_wf%qflx_glcice_frz & ! Output: [real(r8) (:) ] ice growth (positive definite) (mm H2O/s) - ) + qflx_glcice_frz => col_wf%qflx_glcice_frz , & ! Output: [real(r8) (:) ] ice growth (positive definite) (mm H2O/s) + qflx_glcice_diag => col_wf%qflx_glcice_diag , & ! Output: [real(r8) (:) ] flux of new glacier ice (mm H2O/s) - diagnostic, no MECs or GLC + qflx_glcice_frz_diag => col_wf%qflx_glcice_frz_diag & ! Output: [real(r8) (:) ] ice growth (positive definite) (mm H2O/s)) - diagnostic, no MECs or GLC ! Determine time step and step size @@ -222,19 +223,17 @@ subroutine HydrologyDrainage(bounds, & g = col_pp%gridcell(c) ! In the following, we convert glc_snow_persistence_max_days to r8 to avoid overflow if ( (snow_persistence(c) >= (real(glc_snow_persistence_max_days, r8) * secspday)) & - .or. lun_pp%itype(l) == istice_mec .or. lun_pp%itype(l) == istice) then - - if ( (snow_persistence(c) >= (real(glc_snow_persistence_max_days, r8) * secspday)) & - .or. lun_pp%itype(l) == istice_mec) then - qflx_glcice_frz(c) = qflx_snwcp_ice(c) - qflx_glcice(c) = qflx_glcice(c) + qflx_glcice_frz(c) - if (glc_dyn_runoff_routing(g)) qflx_snwcp_ice(c) = 0._r8 - else - qflx_glcice_frz(c) = qflx_snwcp_ice(c) - qflx_glcice(c) = qflx_glcice(c) + qflx_glcice_frz(c) - end if ! lun_pp%itype(l) == istice + .or. lun_pp%itype(l) == istice_mec ) then + qflx_glcice_frz(c) = qflx_snwcp_ice(c) + qflx_glcice(c) = qflx_glcice(c) + qflx_glcice_frz(c) + if (glc_dyn_runoff_routing(g)) qflx_snwcp_ice(c) = 0._r8 end if + if (lun_pp%itype(l)==istice) then + qflx_glcice_frz_diags(c) = qflx_snwcp_ice(c) + qflx_glcice_diags(c) = qflx_glcice_diags(c) + qflx_glcice_frz_diags(c) + endif + end do ! Determine wetland and land ice hydrology (must be placed here diff --git a/components/elm/src/biogeophys/SoilTemperatureMod.F90 b/components/elm/src/biogeophys/SoilTemperatureMod.F90 index be525add7230..b03068fde4fc 100644 --- a/components/elm/src/biogeophys/SoilTemperatureMod.F90 +++ b/components/elm/src/biogeophys/SoilTemperatureMod.F90 @@ -1369,6 +1369,8 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & qflx_snofrz_col => col_wf%qflx_snofrz , & ! Output: [real(r8) (:) ] column-integrated snow freezing rate (positive definite) [kg m-2 s-1] qflx_glcice => col_wf%qflx_glcice , & ! Output: [real(r8) (:) ] flux of new glacier ice (mm H2O/s) [+ = ice grows] qflx_glcice_melt => col_wf%qflx_glcice_melt , & ! Output: [real(r8) (:) ] ice melt (positive definite) (mm H2O/s) + qflx_glcice_diag => col_wf%qflx_glcice_diag , & ! Output: [real(r8) (:) ] flux of new glacier ice (mm H2O/s) [+ = ice grows] + qflx_glcice_melt_diag => col_wf%qflx_glcice_melt_diag , & ! Output: [real(r8) (:) ] ice melt (positive definite) (mm H2O/s) qflx_snomelt => col_wf%qflx_snomelt , & ! Output: [real(r8) (:) ] snow melt (mm H2O /s) qflx_snomelt_lyr => col_wf%qflx_snomelt_lyr , & ! Output: [real(r8) (:) ] snow melt (mm H2O /s) @@ -1663,8 +1665,8 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & if ( lun_pp%itype(l)==istice) then if (j>=1 .and. h2osoi_liq(c,j) > 0._r8) then ! ice layer with meltwater ! melting corresponds to a negative ice flux - qflx_glcice_melt(c) = qflx_glcice_melt(c) + h2osoi_liq(c,j)/dtime - qflx_glcice(c) = qflx_glcice(c) - h2osoi_liq(c,j)/dtime + qflx_glcice_melt_diags(c) = qflx_glcice_melt_diags(c) + h2osoi_liq(c,j)/dtime + qflx_glcice_diags(c) = qflx_glcice_daigs(c) - h2osoi_liq(c,j)/dtime endif ! liquid water is present endif ! istice_mec diff --git a/components/elm/src/data_types/ColumnDataType.F90 b/components/elm/src/data_types/ColumnDataType.F90 index 1de17aafd583..c6a54b3314c5 100644 --- a/components/elm/src/data_types/ColumnDataType.F90 +++ b/components/elm/src/data_types/ColumnDataType.F90 @@ -504,6 +504,9 @@ module ColumnDataType real(r8), pointer :: qflx_glcice (:) => null() ! net flux of new glacial ice (growth - melt) (mm H2O/s), passed to GLC real(r8), pointer :: qflx_glcice_frz (:) => null() ! ice growth (positive definite) (mm H2O/s) real(r8), pointer :: qflx_glcice_melt (:) => null() ! ice melt (positive definite) (mm H2O/s) + real(r8), pointer :: qflx_glcice_diag (:) => null() ! net flux of new glacial ice (growth - melt) (mm H2O/s), passed to GLC + real(r8), pointer :: qflx_glcice_frz_diag (:) => null() ! ice growth (positive definite) (mm H2O/s) + real(r8), pointer :: qflx_glcice_melt_daig(:) => null() ! ice melt (positive definite) (mm H2O/s) real(r8), pointer :: qflx_drain_vr (:,:) => null() ! liquid water lost as drainage (m /time step) real(r8), pointer :: qflx_h2osfc2topsoi (:) => null() ! liquid water coming from surface standing water top soil (mm H2O/s) real(r8), pointer :: qflx_snow2topsoi (:) => null() ! liquid water coming from residual snow to topsoil (mm H2O/s) @@ -5729,6 +5732,9 @@ subroutine col_wf_init(this, begc, endc) allocate(this%qflx_glcice (begc:endc)) ; this%qflx_glcice (:) = spval allocate(this%qflx_glcice_frz (begc:endc)) ; this%qflx_glcice_frz (:) = spval allocate(this%qflx_glcice_melt (begc:endc)) ; this%qflx_glcice_melt (:) = spval + allocate(this%qflx_glcice_diag (begc:endc)) ; this%qflx_glcice_daig (:) = spval + allocate(this%qflx_glcice_frz_diag (begc:endc)) ; this%qflx_glcice_frz_diag (:) = spval + allocate(this%qflx_glcice_melt_diag (begc:endc)) ; this%qflx_glcice_melt_diag(:) = spval allocate(this%qflx_drain_vr (begc:endc,1:nlevgrnd)) ; this%qflx_drain_vr (:,:) = spval allocate(this%qflx_h2osfc2topsoi (begc:endc)) ; this%qflx_h2osfc2topsoi (:) = spval allocate(this%qflx_snow2topsoi (begc:endc)) ; this%qflx_snow2topsoi (:) = spval @@ -5856,21 +5862,39 @@ subroutine col_wf_init(this, begc, endc) call hist_addfld1d (fname='QSNOFRZ', units='kg/m2/s', & avgflag='A', long_name='column-integrated snow freezing rate', & ptr_col=this%qflx_snofrz, set_lake=spval, c2l_scale_type='urbanf') + + if (create_glacier_mec_landunit) then + this%qflx_glcice(begc:endc) = spval + call hist_addfld1d (fname='QICE', units='mm/s', & + avgflag='A', long_name='ice growth/melt', & + ptr_col=this%qflx_glcice, l2g_scale_type='ice') + + this%qflx_glcice_frz(begc:endc) = spval + call hist_addfld1d (fname='QICE_FRZ', units='mm/s', & + avgflag='A', long_name='ice growth', & + ptr_col=this%qflx_glcice_frz, l2g_scale_type='ice') + + this%qflx_glcice_melt(begc:endc) = spval + call hist_addfld1d (fname='QICE_MELT', units='mm/s', & + avgflag='A', long_name='ice melt', & + ptr_col=this%qflx_glcice_melt, l2g_scale_type='ice') + else + this%qflx_glcice_diag(begc:endc) = spval + call hist_addfld1d (fname='QICE', units='mm/s', & + avgflag='A', long_name='diagnostic ice growth/melt (no active GLC/MECs)', & + ptr_col=this%qflx_glcice_diag, l2g_scale_type='ice') + + this%qflx_glcice_frz_diag(begc:endc) = spval + call hist_addfld1d (fname='QICE_FRZ', units='mm/s', & + avgflag='A', long_name='diagnostic ice growth (no active GLC/MECs)', & + ptr_col=this%qflx_glcice_frz_diag, l2g_scale_type='ice') + + this%qflx_glcice_melt_diag(begc:endc) = spval + call hist_addfld1d (fname='QICE_MELT', units='mm/s', & + avgflag='A', long_name='diagnostic ice melt (no active GLC/MECs)', & + ptr_col=this%qflx_glcice_melt_diag, l2g_scale_type='ice') + end if - this%qflx_glcice(begc:endc) = spval - call hist_addfld1d (fname='QICE', units='mm/s', & - avgflag='A', long_name='ice growth/melt', & - ptr_col=this%qflx_glcice, l2g_scale_type='ice') - - this%qflx_glcice_frz(begc:endc) = spval - call hist_addfld1d (fname='QICE_FRZ', units='mm/s', & - avgflag='A', long_name='ice growth', & - ptr_col=this%qflx_glcice_frz, l2g_scale_type='ice') - - this%qflx_glcice_melt(begc:endc) = spval - call hist_addfld1d (fname='QICE_MELT', units='mm/s', & - avgflag='A', long_name='ice melt', & - ptr_col=this%qflx_glcice_melt, l2g_scale_type='ice') ! As defined here, snow_sources - snow_sinks will equal the change in h2osno at any ! given time step but only if there is at least one snow layer (for all landunits From c316747bb9f3e1c8d6adbc728581fc9b84690344 Mon Sep 17 00:00:00 2001 From: Chloe Date: Mon, 7 Oct 2024 13:48:09 -0700 Subject: [PATCH 31/39] generated new diagnostic QICE fields for non MEC/GLC sims --- components/elm/src/biogeophys/HydrologyDrainageMod.F90 | 5 +++-- components/elm/src/biogeophys/SoilTemperatureMod.F90 | 4 ++-- components/elm/src/data_types/ColumnDataType.F90 | 4 ++-- 3 files changed, 7 insertions(+), 6 deletions(-) diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index 1c40981ad4d1..7b58508f399f 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -123,6 +123,7 @@ subroutine HydrologyDrainage(bounds, & qflx_glcice_frz => col_wf%qflx_glcice_frz , & ! Output: [real(r8) (:) ] ice growth (positive definite) (mm H2O/s) qflx_glcice_diag => col_wf%qflx_glcice_diag , & ! Output: [real(r8) (:) ] flux of new glacier ice (mm H2O/s) - diagnostic, no MECs or GLC qflx_glcice_frz_diag => col_wf%qflx_glcice_frz_diag & ! Output: [real(r8) (:) ] ice growth (positive definite) (mm H2O/s)) - diagnostic, no MECs or GLC + ) ! Determine time step and step size @@ -230,8 +231,8 @@ subroutine HydrologyDrainage(bounds, & end if if (lun_pp%itype(l)==istice) then - qflx_glcice_frz_diags(c) = qflx_snwcp_ice(c) - qflx_glcice_diags(c) = qflx_glcice_diags(c) + qflx_glcice_frz_diags(c) + qflx_glcice_frz_diag(c) = qflx_snwcp_ice(c) + qflx_glcice_diag(c) = qflx_glcice_diag(c) + qflx_glcice_frz_diag(c) endif end do diff --git a/components/elm/src/biogeophys/SoilTemperatureMod.F90 b/components/elm/src/biogeophys/SoilTemperatureMod.F90 index b03068fde4fc..bab7cfc61aee 100644 --- a/components/elm/src/biogeophys/SoilTemperatureMod.F90 +++ b/components/elm/src/biogeophys/SoilTemperatureMod.F90 @@ -1665,8 +1665,8 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & if ( lun_pp%itype(l)==istice) then if (j>=1 .and. h2osoi_liq(c,j) > 0._r8) then ! ice layer with meltwater ! melting corresponds to a negative ice flux - qflx_glcice_melt_diags(c) = qflx_glcice_melt_diags(c) + h2osoi_liq(c,j)/dtime - qflx_glcice_diags(c) = qflx_glcice_daigs(c) - h2osoi_liq(c,j)/dtime + qflx_glcice_melt_diag(c) = qflx_glcice_melt_diag(c) + h2osoi_liq(c,j)/dtime + qflx_glcice_diag(c) = qflx_glcice_diag(c) - h2osoi_liq(c,j)/dtime endif ! liquid water is present endif ! istice_mec diff --git a/components/elm/src/data_types/ColumnDataType.F90 b/components/elm/src/data_types/ColumnDataType.F90 index c6a54b3314c5..4abfc4360108 100644 --- a/components/elm/src/data_types/ColumnDataType.F90 +++ b/components/elm/src/data_types/ColumnDataType.F90 @@ -506,7 +506,7 @@ module ColumnDataType real(r8), pointer :: qflx_glcice_melt (:) => null() ! ice melt (positive definite) (mm H2O/s) real(r8), pointer :: qflx_glcice_diag (:) => null() ! net flux of new glacial ice (growth - melt) (mm H2O/s), passed to GLC real(r8), pointer :: qflx_glcice_frz_diag (:) => null() ! ice growth (positive definite) (mm H2O/s) - real(r8), pointer :: qflx_glcice_melt_daig(:) => null() ! ice melt (positive definite) (mm H2O/s) + real(r8), pointer :: qflx_glcice_melt_diag(:) => null() ! ice melt (positive definite) (mm H2O/s) real(r8), pointer :: qflx_drain_vr (:,:) => null() ! liquid water lost as drainage (m /time step) real(r8), pointer :: qflx_h2osfc2topsoi (:) => null() ! liquid water coming from surface standing water top soil (mm H2O/s) real(r8), pointer :: qflx_snow2topsoi (:) => null() ! liquid water coming from residual snow to topsoil (mm H2O/s) @@ -5732,7 +5732,7 @@ subroutine col_wf_init(this, begc, endc) allocate(this%qflx_glcice (begc:endc)) ; this%qflx_glcice (:) = spval allocate(this%qflx_glcice_frz (begc:endc)) ; this%qflx_glcice_frz (:) = spval allocate(this%qflx_glcice_melt (begc:endc)) ; this%qflx_glcice_melt (:) = spval - allocate(this%qflx_glcice_diag (begc:endc)) ; this%qflx_glcice_daig (:) = spval + allocate(this%qflx_glcice_diag (begc:endc)) ; this%qflx_glcice_diag (:) = spval allocate(this%qflx_glcice_frz_diag (begc:endc)) ; this%qflx_glcice_frz_diag (:) = spval allocate(this%qflx_glcice_melt_diag (begc:endc)) ; this%qflx_glcice_melt_diag(:) = spval allocate(this%qflx_drain_vr (begc:endc,1:nlevgrnd)) ; this%qflx_drain_vr (:,:) = spval From 7aff30e2b5729d01508608ae083d5635289e2b2d Mon Sep 17 00:00:00 2001 From: Chloe Date: Mon, 7 Oct 2024 18:00:53 -0700 Subject: [PATCH 32/39] debugging QICE --- components/elm/src/biogeophys/HydrologyDrainageMod.F90 | 3 +++ components/elm/src/biogeophys/SoilTemperatureMod.F90 | 4 ++++ components/elm/src/main/elm_driver.F90 | 3 +++ 3 files changed, 10 insertions(+) diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index 7b58508f399f..c7a69b7d3194 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -217,6 +217,7 @@ subroutine HydrologyDrainage(bounds, & do c = bounds%begc,bounds%endc qflx_glcice_frz(c) = 0._r8 + qflx_glcice_frz_diag(c) = 0._r8 end do do fc = 1,num_do_smb_c c = filter_do_smb_c(fc) @@ -233,6 +234,8 @@ subroutine HydrologyDrainage(bounds, & if (lun_pp%itype(l)==istice) then qflx_glcice_frz_diag(c) = qflx_snwcp_ice(c) qflx_glcice_diag(c) = qflx_glcice_diag(c) + qflx_glcice_frz_diag(c) + !write(iulog,*) 'CAW lun_pp%itype(l)==istice',lun_pp%itype(l)==istice + !write(iulog,*) 'qflx_snwcp_ice(c)',qflx_snwcp_ice(c) endif end do diff --git a/components/elm/src/biogeophys/SoilTemperatureMod.F90 b/components/elm/src/biogeophys/SoilTemperatureMod.F90 index bab7cfc61aee..9508e54849c6 100644 --- a/components/elm/src/biogeophys/SoilTemperatureMod.F90 +++ b/components/elm/src/biogeophys/SoilTemperatureMod.F90 @@ -1397,6 +1397,7 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & qflx_snofrz_lyr(c,-nlevsno+1:0) = 0._r8 qflx_snofrz_col(c) = 0._r8 qflx_glcice_melt(c) = 0._r8 + qflx_glcice_melt_diag(c) = 0._r8 qflx_snow_melt(c) = 0._r8 end do @@ -1667,6 +1668,9 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & ! melting corresponds to a negative ice flux qflx_glcice_melt_diag(c) = qflx_glcice_melt_diag(c) + h2osoi_liq(c,j)/dtime qflx_glcice_diag(c) = qflx_glcice_diag(c) - h2osoi_liq(c,j)/dtime + !write(iulog,*) 'CAW lun_pp%itype(l)==istice',lun_pp%itype(l)==istice + !write(iulog,*) 'CAW j',j + !write(iulog,*) 'CAW h2osoi_liq(c,j) ',h2osoi_liq(c,j) endif ! liquid water is present endif ! istice_mec diff --git a/components/elm/src/main/elm_driver.F90 b/components/elm/src/main/elm_driver.F90 index 51f08bf235c3..18ecae88e37f 100644 --- a/components/elm/src/main/elm_driver.F90 +++ b/components/elm/src/main/elm_driver.F90 @@ -1604,6 +1604,8 @@ subroutine elm_drv_init(bounds, & qflx_glcice => col_wf%qflx_glcice , & ! Output: [real(r8) (:) ] flux of new glacier ice (mm H2O/s) [+ = ice grows] + qflx_glcice_diag => col_wf%qflx_glcice_diag , & ! Output: [real(r8) (:) ] flux of new glacier ice (mm H2O/s) [+ = ice grows] + eflx_bot => col_ef%eflx_bot , & ! Output: [real(r8) (:) ] heat flux from beneath soil/ice column (W/m**2) cisun_z => photosyns_vars%cisun_z_patch , & ! Output: [real(r8) (:) ] intracellular sunlit leaf CO2 (Pa) @@ -1637,6 +1639,7 @@ subroutine elm_drv_init(bounds, & ! Initialize qflx_glcice everywhere, to zero. qflx_glcice(c) = 0._r8 + qflx_glcice_diag(c) = 0._r8 end do From 78897e7fdf7f964a9304d6cba0484ab9f8650dda Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 8 Oct 2024 09:31:37 -0700 Subject: [PATCH 33/39] removed debugging write statements --- components/elm/src/biogeophys/HydrologyDrainageMod.F90 | 2 -- components/elm/src/biogeophys/SoilTemperatureMod.F90 | 3 --- 2 files changed, 5 deletions(-) diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index c7a69b7d3194..41a646705ef4 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -234,8 +234,6 @@ subroutine HydrologyDrainage(bounds, & if (lun_pp%itype(l)==istice) then qflx_glcice_frz_diag(c) = qflx_snwcp_ice(c) qflx_glcice_diag(c) = qflx_glcice_diag(c) + qflx_glcice_frz_diag(c) - !write(iulog,*) 'CAW lun_pp%itype(l)==istice',lun_pp%itype(l)==istice - !write(iulog,*) 'qflx_snwcp_ice(c)',qflx_snwcp_ice(c) endif end do diff --git a/components/elm/src/biogeophys/SoilTemperatureMod.F90 b/components/elm/src/biogeophys/SoilTemperatureMod.F90 index 9508e54849c6..72797e719e6e 100644 --- a/components/elm/src/biogeophys/SoilTemperatureMod.F90 +++ b/components/elm/src/biogeophys/SoilTemperatureMod.F90 @@ -1668,9 +1668,6 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & ! melting corresponds to a negative ice flux qflx_glcice_melt_diag(c) = qflx_glcice_melt_diag(c) + h2osoi_liq(c,j)/dtime qflx_glcice_diag(c) = qflx_glcice_diag(c) - h2osoi_liq(c,j)/dtime - !write(iulog,*) 'CAW lun_pp%itype(l)==istice',lun_pp%itype(l)==istice - !write(iulog,*) 'CAW j',j - !write(iulog,*) 'CAW h2osoi_liq(c,j) ',h2osoi_liq(c,j) endif ! liquid water is present endif ! istice_mec From f779abd0aa6f4c29e991aa3338871ec175b1d8d0 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 15 Oct 2024 10:18:58 -0700 Subject: [PATCH 34/39] move QICE_FRZ diag calc outside of do_smb loop --- components/elm/src/biogeophys/HydrologyDrainageMod.F90 | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index 41a646705ef4..170dbe0ab59b 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -218,6 +218,11 @@ subroutine HydrologyDrainage(bounds, & do c = bounds%begc,bounds%endc qflx_glcice_frz(c) = 0._r8 qflx_glcice_frz_diag(c) = 0._r8 + + if (lun_pp%itype(l)==istice .and. qflx_snwcp_ice(c) > 0.0_r8) then + qflx_glcice_frz_diag(c) = qflx_snwcp_ice(c) + qflx_glcice_diag(c) = qflx_glcice_diag(c) + qflx_glcice_frz_diag(c) + endif end do do fc = 1,num_do_smb_c c = filter_do_smb_c(fc) From 06c6a629a5765caee6bd24b27cbddc788e2efa1b Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 15 Oct 2024 13:15:17 -0700 Subject: [PATCH 35/39] remove old QICE diag field calc from do_smb loop --- components/elm/src/biogeophys/HydrologyDrainageMod.F90 | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index 170dbe0ab59b..bc19e924af2d 100755 --- a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 +++ b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 @@ -236,10 +236,10 @@ subroutine HydrologyDrainage(bounds, & if (glc_dyn_runoff_routing(g)) qflx_snwcp_ice(c) = 0._r8 end if - if (lun_pp%itype(l)==istice) then - qflx_glcice_frz_diag(c) = qflx_snwcp_ice(c) - qflx_glcice_diag(c) = qflx_glcice_diag(c) + qflx_glcice_frz_diag(c) - endif + !if (lun_pp%itype(l)==istice) then + ! qflx_glcice_frz_diag(c) = qflx_snwcp_ice(c) + ! qflx_glcice_diag(c) = qflx_glcice_diag(c) + qflx_glcice_frz_diag(c) + !endif end do From 93624bbbac3d01e05a5d284b82c8aa7467f4e183 Mon Sep 17 00:00:00 2001 From: Chloe Date: Sat, 19 Oct 2024 10:44:53 -0700 Subject: [PATCH 36/39] added PR #6682 to maint branch for deep firn simulations From f4688baf7f2cb4f991431795031efc967104fca4 Mon Sep 17 00:00:00 2001 From: Chloe Date: Tue, 21 Jan 2025 16:53:39 -0800 Subject: [PATCH 37/39] minor changes for PR --- components/elm/src/biogeophys/SoilTemperatureMod.F90 | 3 +-- components/elm/src/data_types/ColumnDataType.F90 | 2 +- 2 files changed, 2 insertions(+), 3 deletions(-) diff --git a/components/elm/src/biogeophys/SoilTemperatureMod.F90 b/components/elm/src/biogeophys/SoilTemperatureMod.F90 index 72797e719e6e..6826b3e22e63 100644 --- a/components/elm/src/biogeophys/SoilTemperatureMod.F90 +++ b/components/elm/src/biogeophys/SoilTemperatureMod.F90 @@ -1628,8 +1628,7 @@ subroutine Phasechange_beta (bounds, num_nolakec, filter_nolakec, dhsdT, & if (imelt(c,j) == 1 .AND. j < 1) then qflx_snomelt(c) = qflx_snomelt(c) + max(0._r8,(wice0(c,j)-h2osoi_ice(c,j)))/dtime - qflx_snomelt_lyr(c,j) = qflx_snomelt(c) - + qflx_snomelt_lyr(c,j) = max(0._r8,(wice0(c,j)-h2osoi_ice(c,j)))/dtime endif diff --git a/components/elm/src/data_types/ColumnDataType.F90 b/components/elm/src/data_types/ColumnDataType.F90 index 4abfc4360108..33a7d1d9c27e 100644 --- a/components/elm/src/data_types/ColumnDataType.F90 +++ b/components/elm/src/data_types/ColumnDataType.F90 @@ -5692,7 +5692,7 @@ subroutine col_wf_init(this, begc, endc) allocate(this%qflx_evap_grnd (begc:endc)) ; this%qflx_evap_grnd (:) = spval allocate(this%qflx_snwcp_liq (begc:endc)) ; this%qflx_snwcp_liq (:) = spval allocate(this%qflx_snwcp_ice (begc:endc)) ; this%qflx_snwcp_ice (:) = spval - allocate(this%qflx_ice_runoff_xs (begc:endc)) ; this%qflx_ice_runoff_xs (:) = 0._r8 + allocate(this%qflx_ice_runoff_xs (begc:endc)) ; this%qflx_ice_runoff_xs (:) = spval allocate(this%qflx_tran_veg (begc:endc)) ; this%qflx_tran_veg (:) = spval allocate(this%qflx_dew_snow (begc:endc)) ; this%qflx_dew_snow (:) = spval allocate(this%qflx_dew_grnd (begc:endc)) ; this%qflx_dew_grnd (:) = spval From 9b7b29c8080e7f0c00abb8ee75c6df3bb8f90b5c Mon Sep 17 00:00:00 2001 From: Chloe Date: Wed, 22 Jan 2025 09:52:08 -0800 Subject: [PATCH 38/39] minor PR changes & testing --- components/elm/src/data_types/ColumnDataType.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/components/elm/src/data_types/ColumnDataType.F90 b/components/elm/src/data_types/ColumnDataType.F90 index 33a7d1d9c27e..4abfc4360108 100644 --- a/components/elm/src/data_types/ColumnDataType.F90 +++ b/components/elm/src/data_types/ColumnDataType.F90 @@ -5692,7 +5692,7 @@ subroutine col_wf_init(this, begc, endc) allocate(this%qflx_evap_grnd (begc:endc)) ; this%qflx_evap_grnd (:) = spval allocate(this%qflx_snwcp_liq (begc:endc)) ; this%qflx_snwcp_liq (:) = spval allocate(this%qflx_snwcp_ice (begc:endc)) ; this%qflx_snwcp_ice (:) = spval - allocate(this%qflx_ice_runoff_xs (begc:endc)) ; this%qflx_ice_runoff_xs (:) = spval + allocate(this%qflx_ice_runoff_xs (begc:endc)) ; this%qflx_ice_runoff_xs (:) = 0._r8 allocate(this%qflx_tran_veg (begc:endc)) ; this%qflx_tran_veg (:) = spval allocate(this%qflx_dew_snow (begc:endc)) ; this%qflx_dew_snow (:) = spval allocate(this%qflx_dew_grnd (begc:endc)) ; this%qflx_dew_grnd (:) = spval From 16383cb081a12e5eebae864de2205a308dd6accd Mon Sep 17 00:00:00 2001 From: Chloe Date: Wed, 22 Jan 2025 14:16:52 -0800 Subject: [PATCH 39/39] remove MALI fully coupled B case 20thTR config --- cime_config/allactive/config_compsets.xml | 4 ---- 1 file changed, 4 deletions(-) diff --git a/cime_config/allactive/config_compsets.xml b/cime_config/allactive/config_compsets.xml index f474f4932ebc..d5a52a8bac38 100755 --- a/cime_config/allactive/config_compsets.xml +++ b/cime_config/allactive/config_compsets.xml @@ -459,10 +459,6 @@ 1850_EAM%CMIP6_ELM%SPBC_MPASSI_MPASO_MOSART_MALI_SWAV - - BGWCYCL20TR - 20TRSOI_EAM%CMIP6_ELM%CNPRDCTCBCTOP_MPASSI_MPASO_MOSART_MALI_SWAV -