diff --git a/components/elm/src/biogeophys/HydrologyDrainageMod.F90 b/components/elm/src/biogeophys/HydrologyDrainageMod.F90 index 850b7eedac52..fbd4d8e2b6f3 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 @@ -120,7 +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 @@ -215,6 +217,12 @@ 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) @@ -222,10 +230,10 @@ 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 + .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 end do 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)) diff --git a/components/elm/src/biogeophys/SoilTemperatureMod.F90 b/components/elm/src/biogeophys/SoilTemperatureMod.F90 index d4a1074bf4a7..451b3cbccf22 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 @@ -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) eflx_snomelt => col_ef%eflx_snomelt , & ! Output: [real(r8) (:) ] snow melt heat flux (W/m**2) @@ -1393,6 +1395,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 @@ -1643,8 +1646,8 @@ 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_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 @@ -1656,6 +1659,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_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 + end do ! end of column-loop enddo ! end of level-loop diff --git a/components/elm/src/data_types/ColumnDataType.F90 b/components/elm/src/data_types/ColumnDataType.F90 index 724cefb5cfdf..fd1a17d74d63 100644 --- a/components/elm/src/data_types/ColumnDataType.F90 +++ b/components/elm/src/data_types/ColumnDataType.F90 @@ -502,6 +502,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_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) @@ -5725,6 +5728,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_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 allocate(this%qflx_h2osfc2topsoi (begc:endc)) ; this%qflx_h2osfc2topsoi (:) = spval allocate(this%qflx_snow2topsoi (begc:endc)) ; this%qflx_snow2topsoi (:) = spval @@ -5842,23 +5848,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', default='inactive') - + 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') + 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 + ! 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 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