INTERFACE:
subroutine SnowAge_grain(lbc, ubc, num_snowc, filter_snowc, num_nosnowc, filter_nosnowc)DESCRIPTION:
! Updates the snow effective grain size (radius). ! Contributions to grain size evolution are from: ! 1. vapor redistribution (dry snow) ! 2. liquid water redistribution (wet snow) ! 3. re-freezing of liquid water ! ! Vapor redistribution: Method is to retrieve 3 best-bit parameters that ! depend on snow temperature, temperature gradient, and density, ! that are derived from the microphysical model described in: ! Flanner and Zender (2006), Linking snowpack microphysics and albedo ! evolution, J. Geophys. Res., 111, D12208, doi:10.1029/2005JD006834. ! The parametric equation has the form: ! dr/dt = drdt_0*(tau/(dr_fresh+tau))^(1/kappa), where: ! r is the effective radius, ! tau and kappa are best-fit parameters, ! drdt_0 is the initial rate of change of effective radius, and ! dr_fresh is the difference between the current and fresh snow states ! (r_current - r_fresh). ! ! Liquid water redistribution: Apply the grain growth function from: ! Brun, E. (1989), Investigation of wet-snow metamorphism in respect of ! liquid-water content, Annals of Glaciology, 13, 22-26. ! There are two parameters that describe the grain growth rate as ! a function of snow liquid water content (LWC). The "LWC=0" parameter ! is zeroed here because we are accounting for dry snowing with a ! different representation ! ! Re-freezing of liquid water: Assume that re-frozen liquid water clumps ! into an arbitrarily large effective grain size (snw_rds_refrz). ! The phenomenon is observed (Grenfell), but so far unquantified, as far as ! I am aware. ! !
USES:
use clmtype use clm_time_manager , only : get_step_size, get_nstep use clm_varpar , only : nlevsno use clm_varcon , only : spval use abortutils , only : endrun use shr_const_mod , only : SHR_CONST_RHOICE, SHR_CONST_PIARGUMENTS:
implicit none integer, intent(in) :: lbc, ubc ! column bounds integer, intent(in) :: num_snowc ! number of column snow points in column filter integer, intent(in) :: filter_snowc(ubc-lbc+1) ! column filter for snow points integer, intent(in) :: num_nosnowc ! number of column non-snow points in column filter integer, intent(in) :: filter_nosnowc(ubc-lbc+1) ! column filter for non-snow pointsCALLED FROM:
LOCAL VARIABLES:
! local pointers to implicit arguments real(r8), pointer :: t_soisno(:,:) ! soil and snow temperature (col,lyr) [K] integer, pointer :: snl(:) ! negative number of snow layers (col) [nbr] real(r8), pointer :: t_grnd(:) ! ground temperature (col) [K] real(r8), pointer :: dz(:,:) ! layer thickness (col,lyr) [m] real(r8), pointer :: h2osno(:) ! snow water (col) [mm H2O] real(r8), pointer :: snw_rds(:,:) ! effective grain radius (col,lyr) [microns, m-6] real(r8), pointer :: snw_rds_top(:) ! effective grain radius, top layer (col) [microns, m-6] real(r8), pointer :: sno_liq_top(:) ! liquid water fraction (mass) in top snow layer (col) [frc] real(r8), pointer :: h2osoi_liq(:,:) ! liquid water content (col,lyr) [kg m-2] real(r8), pointer :: h2osoi_ice(:,:) ! ice content (col,lyr) [kg m-2] real(r8), pointer :: snot_top(:) ! snow temperature in top layer (col) [K] real(r8), pointer :: dTdz_top(:) ! temperature gradient in top layer (col) [K m-1] real(r8), pointer :: qflx_snow_grnd_col(:) ! snow on ground after interception (col) [kg m-2 s-1] real(r8), pointer :: qflx_snwcp_ice(:) ! excess precipitation due to snow capping [kg m-2 s-1] real(r8), pointer :: qflx_snofrz_lyr(:,:) ! snow freezing rate (col,lyr) [kg m-2 s-1] logical , pointer :: do_capsnow(:) ! true => do snow capping ! !OTHER LOCAL VARIABLES: integer :: snl_top ! top snow layer index [idx] integer :: snl_btm ! bottom snow layer index [idx] integer :: i ! layer index [idx] integer :: c_idx ! column index [idx] integer :: fc ! snow column filter index [idx] integer :: T_idx ! snow aging lookup table temperature index [idx] integer :: Tgrd_idx ! snow aging lookup table temperature gradient index [idx] integer :: rhos_idx ! snow aging lookup table snow density index [idx] real(r8) :: t_snotop ! temperature at upper layer boundary [K] real(r8) :: t_snobtm ! temperature at lower layer boundary [K] real(r8) :: dTdz(lbc:ubc,-nlevsno:0) ! snow temperature gradient (col,lyr) [K m-1] real(r8) :: bst_tau ! snow aging parameter retrieved from lookup table [hour] real(r8) :: bst_kappa ! snow aging parameter retrieved from lookup table [unitless] real(r8) :: bst_drdt0 ! snow aging parameter retrieved from lookup table [um hr-1] real(r8) :: dr ! incremental change in snow effective radius [um] real(r8) :: dr_wet ! incremental change in snow effective radius from wet growth [um] real(r8) :: dr_fresh ! difference between fresh snow r_e and current r_e [um] real(r8) :: newsnow ! fresh snowfall [kg m-2] real(r8) :: refrzsnow ! re-frozen snow [kg m-2] real(r8) :: frc_newsnow ! fraction of layer mass that is new snow [frc] real(r8) :: frc_oldsnow ! fraction of layer mass that is old snow [frc] real(r8) :: frc_refrz ! fraction of layer mass that is re-frozen snow [frc] real(r8) :: frc_liq ! fraction of layer mass that is liquid water[frc] real(r8) :: dtime ! land model time step [sec] real(r8) :: rhos ! snow density [kg m-3] real(r8) :: h2osno_lyr ! liquid + solid H2O in snow layer [kg m-2] ! Assign local pointers to derived subtypes components (column-level) t_soisno => clm3%g%l%c%ces%t_soisno snl => clm3%g%l%c%cps%snl t_grnd => clm3%g%l%c%ces%t_grnd dz => clm3%g%l%c%cps%dz h2osno => clm3%g%l%c%cws%h2osno snw_rds => clm3%g%l%c%cps%snw_rds h2osoi_liq => clm3%g%l%c%cws%h2osoi_liq h2osoi_ice => clm3%g%l%c%cws%h2osoi_ice snot_top => clm3%g%l%c%cps%snot_top dTdz_top => clm3%g%l%c%cps%dTdz_top snw_rds_top => clm3%g%l%c%cps%snw_rds_top sno_liq_top => clm3%g%l%c%cps%sno_liq_top qflx_snow_grnd_col => clm3%g%l%c%cwf%pwf_a%qflx_snow_grnd qflx_snwcp_ice => clm3%g%l%c%cwf%pwf_a%qflx_snwcp_ice qflx_snofrz_lyr => clm3%g%l%c%cwf%qflx_snofrz_lyr do_capsnow => clm3%g%l%c%cps%do_capsnow ! set timestep and step interval dtime = get_step_size() ! loop over columns that have at least one snow layer do fc = 1, num_snowc c_idx = filter_snowc(fc) snl_btm = 0 snl_top = snl(c_idx) + 1 ! loop over snow layers do i=snl_top,snl_btm,1 !********** 1. DRY SNOW AGING *********** h2osno_lyr = h2osoi_liq(c_idx,i) + h2osoi_ice(c_idx,i) ! temperature gradient if (i == snl_top) then ! top layer t_snotop = t_grnd(c_idx) t_snobtm = (t_soisno(c_idx,i+1)*dz(c_idx,i) + t_soisno(c_idx,i)*dz(c_idx,i+1)) / (dz(c_idx,i)+dz(c_idx,i+1)) else t_snotop = (t_soisno(c_idx,i-1)*dz(c_idx,i) + t_soisno(c_idx,i)*dz(c_idx,i-1)) / (dz(c_idx,i)+dz(c_idx,i-1)) t_snobtm = (t_soisno(c_idx,i+1)*dz(c_idx,i) + t_soisno(c_idx,i)*dz(c_idx,i+1)) / (dz(c_idx,i)+dz(c_idx,i+1)) endif dTdz(c_idx,i) = abs((t_snotop - t_snobtm) / dz(c_idx,i)) ! snow density rhos = (h2osoi_liq(c_idx,i)+h2osoi_ice(c_idx,i)) / dz(c_idx,i) ! best-fit table indecies T_idx = nint((t_soisno(c_idx,i)-223) / 5) + 1 Tgrd_idx = nint(dTdz(c_idx,i) / 10) + 1 rhos_idx = nint((rhos-50) / 50) + 1 ! boundary check: if (T_idx < idx_T_min) then T_idx = idx_T_min endif if (T_idx > idx_T_max) then T_idx = idx_T_max endif if (Tgrd_idx < idx_Tgrd_min) then Tgrd_idx = idx_Tgrd_min endif if (Tgrd_idx > idx_Tgrd_max) then Tgrd_idx = idx_Tgrd_max endif if (rhos_idx < idx_rhos_min) then rhos_idx = idx_rhos_min endif if (rhos_idx > idx_rhos_max) then rhos_idx = idx_rhos_max endif ! best-fit parameters bst_tau = snowage_tau(rhos_idx,Tgrd_idx,T_idx) bst_kappa = snowage_kappa(rhos_idx,Tgrd_idx,T_idx) bst_drdt0 = snowage_drdt0(rhos_idx,Tgrd_idx,T_idx) ! change in snow effective radius, using best-fit parameters dr_fresh = snw_rds(c_idx,i)-snw_rds_min dr = (bst_drdt0*(bst_tau/(dr_fresh+bst_tau))**(1/bst_kappa)) * (dtime/3600) !********** 2. WET SNOW AGING *********** ! We are assuming wet and dry evolution occur simultaneously, and ! the contributions from both can be summed. ! This is justified by setting the linear offset constant C1_liq_Brun89 to zero [Brun, 1989] ! liquid water faction frc_liq = min(0.1_r8, (h2osoi_liq(c_idx,i) / (h2osoi_liq(c_idx,i)+h2osoi_ice(c_idx,i)))) !dr_wet = 1E6_r8*(dtime*(C1_liq_Brun89 + C2_liq_Brun89*(frc_liq**(3))) / (4*SHR_CONST_PI*(snw_rds(c_idx,i)/1E6)**(2))) !simplified, units of microns: dr_wet = 1E18_r8*(dtime*(C2_liq_Brun89*(frc_liq**(3))) / (4*SHR_CONST_PI*snw_rds(c_idx,i)**(2))) dr = dr + dr_wet !********** 3. SNOWAGE SCALING (TURNED OFF BY DEFAULT) ************* ! Multiply rate of change of effective radius by some constant, xdrdt if (flg_snoage_scl) then dr = dr*xdrdt endif !********** 4. INCREMENT EFFECTIVE RADIUS, ACCOUNTING FOR: *********** ! DRY AGING ! WET AGING ! FRESH SNOW ! RE-FREEZING ! new snowfall [kg/m2] if (do_capsnow(c_idx)) then newsnow = max(0._r8, (qflx_snwcp_ice(c_idx)*dtime)) else newsnow = max(0._r8, (qflx_snow_grnd_col(c_idx)*dtime)) endif ! snow that has re-frozen [kg/m2] refrzsnow = max(0._r8, (qflx_snofrz_lyr(c_idx,i)*dtime)) ! fraction of layer mass that is re-frozen frc_refrz = refrzsnow / h2osno_lyr ! fraction of layer mass that is new snow if (i == snl_top) then frc_newsnow = newsnow / h2osno_lyr else frc_newsnow = 0._r8 endif if ((frc_refrz + frc_newsnow) > 1._r8) then frc_refrz = frc_refrz / (frc_refrz + frc_newsnow) frc_newsnow = 1._r8 - frc_refrz frc_oldsnow = 0._r8 else frc_oldsnow = 1._r8 - frc_refrz - frc_newsnow endif ! mass-weighted mean of fresh snow, old snow, and re-frozen snow effective radius snw_rds(c_idx,i) = (snw_rds(c_idx,i)+dr)*frc_oldsnow + snw_rds_min*frc_newsnow + snw_rds_refrz*frc_refrz !********** 5. CHECK BOUNDARIES *********** ! boundary check if (snw_rds(c_idx,i) < snw_rds_min) then snw_rds(c_idx,i) = snw_rds_min endif if (snw_rds(c_idx,i) > snw_rds_max) then snw_rds(c_idx,i) = snw_rds_max end if ! set top layer variables for history files if (i == snl_top) then snot_top(c_idx) = t_soisno(c_idx,i) dTdz_top(c_idx) = dTdz(c_idx,i) snw_rds_top(c_idx) = snw_rds(c_idx,i) sno_liq_top(c_idx) = h2osoi_liq(c_idx,i) / (h2osoi_liq(c_idx,i)+h2osoi_ice(c_idx,i)) endif enddo enddo ! Special case: snow on ground, but not enough to have defined a snow layer: ! set snw_rds to fresh snow grain size: do fc = 1, num_nosnowc c_idx = filter_nosnowc(fc) if (h2osno(c_idx) > 0._r8) then snw_rds(c_idx,0) = snw_rds_min endif enddo end subroutine SnowAge_grain subroutine SnowOptics_init( ) use fileutils , only : getfil use CLM_varctl , only : fsnowoptics use spmdMod , only : masterproc use ncdio_pio , only : file_desc_t, ncd_io, ncd_pio_openfile, ncd_pio_closefile type(file_desc_t) :: ncid ! netCDF file id character(len=256) :: locfn ! local filename character(len= 32) :: subname = 'SnowOptics_init' ! subroutine name integer :: ier ! error status ! Open optics file: if(masterproc) write(iulog,*) 'Attempting to read snow optical properties .....' call getfil (fsnowoptics, locfn, 0) call ncd_pio_openfile(ncid, locfn, 0) if(masterproc) write(iulog,*) subname,trim(fsnowoptics) ! direct-beam snow Mie parameters: call ncd_io('ss_alb_ice_drc', ss_alb_snw_drc, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'asm_prm_ice_drc',asm_prm_snw_drc, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'ext_cff_mss_ice_drc', ext_cff_mss_snw_drc, 'read', ncid, posNOTonfile=.true.) ! diffuse snow Mie parameters call ncd_io( 'ss_alb_ice_dfs', ss_alb_snw_dfs, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'asm_prm_ice_dfs', asm_prm_snw_dfs, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'ext_cff_mss_ice_dfs', ext_cff_mss_snw_dfs, 'read', ncid, posNOTonfile=.true.) ! BC species 1 Mie parameters call ncd_io( 'ss_alb_bcphil', ss_alb_bc1, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'asm_prm_bcphil', asm_prm_bc1, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'ext_cff_mss_bcphil', ext_cff_mss_bc1, 'read', ncid, posNOTonfile=.true.) ! BC species 2 Mie parameters call ncd_io( 'ss_alb_bcphob', ss_alb_bc2, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'asm_prm_bcphob', asm_prm_bc2, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'ext_cff_mss_bcphob', ext_cff_mss_bc2, 'read', ncid, posNOTonfile=.true.) ! OC species 1 Mie parameters call ncd_io( 'ss_alb_ocphil', ss_alb_oc1, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'asm_prm_ocphil', asm_prm_oc1, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'ext_cff_mss_ocphil', ext_cff_mss_oc1, 'read', ncid, posNOTonfile=.true.) ! OC species 2 Mie parameters call ncd_io( 'ss_alb_ocphob', ss_alb_oc2, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'asm_prm_ocphob', asm_prm_oc2, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'ext_cff_mss_ocphob', ext_cff_mss_oc2, 'read', ncid, posNOTonfile=.true.) ! dust species 1 Mie parameters call ncd_io( 'ss_alb_dust01', ss_alb_dst1, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'asm_prm_dust01', asm_prm_dst1, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'ext_cff_mss_dust01', ext_cff_mss_dst1, 'read', ncid, posNOTonfile=.true.) ! dust species 2 Mie parameters call ncd_io( 'ss_alb_dust02', ss_alb_dst2, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'asm_prm_dust02', asm_prm_dst2, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'ext_cff_mss_dust02', ext_cff_mss_dst2, 'read', ncid, posNOTonfile=.true.) ! dust species 3 Mie parameters call ncd_io( 'ss_alb_dust03', ss_alb_dst3, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'asm_prm_dust03', asm_prm_dst3, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'ext_cff_mss_dust03', ext_cff_mss_dst3, 'read', ncid, posNOTonfile=.true.) ! dust species 4 Mie parameters call ncd_io( 'ss_alb_dust04', ss_alb_dst4, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'asm_prm_dust04', asm_prm_dst4, 'read', ncid, posNOTonfile=.true.) call ncd_io( 'ext_cff_mss_dust04', ext_cff_mss_dst4, 'read', ncid, posNOTonfile=.true.) call ncd_pio_closefile(ncid) if (masterproc) then write(iulog,*) 'Successfully read snow optical properties' ! print some diagnostics: write (iulog,*) 'SNICAR: Mie single scatter albedos for direct-beam ice, rds=100um: ', & ss_alb_snw_drc(71,1), ss_alb_snw_drc(71,2), ss_alb_snw_drc(71,3), & ss_alb_snw_drc(71,4), ss_alb_snw_drc(71,5) write (iulog,*) 'SNICAR: Mie single scatter albedos for diffuse ice, rds=100um: ', & ss_alb_snw_dfs(71,1), ss_alb_snw_dfs(71,2), ss_alb_snw_dfs(71,3), & ss_alb_snw_dfs(71,4), ss_alb_snw_dfs(71,5) if (DO_SNO_OC) then write (iulog,*) 'SNICAR: Including OC aerosols from snow radiative transfer calculations' else write (iulog,*) 'SNICAR: Excluding OC aerosols from snow radiative transfer calculations' endif write (iulog,*) 'SNICAR: Mie single scatter albedos for hydrophillic BC: ', & ss_alb_bc1(1), ss_alb_bc1(2), ss_alb_bc1(3), ss_alb_bc1(4), ss_alb_bc1(5) write (iulog,*) 'SNICAR: Mie single scatter albedos for hydrophobic BC: ', & ss_alb_bc2(1), ss_alb_bc2(2), ss_alb_bc2(3), ss_alb_bc2(4), ss_alb_bc2(5) if (DO_SNO_OC) then write (iulog,*) 'SNICAR: Mie single scatter albedos for hydrophillic OC: ', & ss_alb_oc1(1), ss_alb_oc1(2), ss_alb_oc1(3), ss_alb_oc1(4), ss_alb_oc1(5) write (iulog,*) 'SNICAR: Mie single scatter albedos for hydrophobic OC: ', & ss_alb_oc2(1), ss_alb_oc2(2), ss_alb_oc2(3), ss_alb_oc2(4), ss_alb_oc2(5) endif write (iulog,*) 'SNICAR: Mie single scatter albedos for dust species 1: ', & ss_alb_dst1(1), ss_alb_dst1(2), ss_alb_dst1(3), ss_alb_dst1(4), ss_alb_dst1(5) write (iulog,*) 'SNICAR: Mie single scatter albedos for dust species 2: ', & ss_alb_dst2(1), ss_alb_dst2(2), ss_alb_dst2(3), ss_alb_dst2(4), ss_alb_dst2(5) write (iulog,*) 'SNICAR: Mie single scatter albedos for dust species 3: ', & ss_alb_dst3(1), ss_alb_dst3(2), ss_alb_dst3(3), ss_alb_dst3(4), ss_alb_dst3(5) write (iulog,*) 'SNICAR: Mie single scatter albedos for dust species 4: ', & ss_alb_dst4(1), ss_alb_dst4(2), ss_alb_dst4(3), ss_alb_dst4(4), ss_alb_dst4(5) write(iulog,*) end if end subroutine SnowOptics_init subroutine SnowAge_init( ) use CLM_varctl , only : fsnowaging use fileutils , only : getfil use spmdMod , only : masterproc use ncdio_pio , only : file_desc_t, ncd_io, ncd_pio_openfile, ncd_pio_closefile type(file_desc_t) :: ncid ! netCDF file id character(len=256) :: locfn ! local filename character(len= 32) :: subname = 'SnowOptics_init' ! subroutine name integer :: varid ! netCDF id's integer :: ier ! error status ! Open snow aging (effective radius evolution) file: allocate(snowage_tau(idx_rhos_max,idx_Tgrd_max,idx_T_max)) allocate(snowage_kappa(idx_rhos_max,idx_Tgrd_max,idx_T_max)) allocate(snowage_drdt0(idx_rhos_max,idx_Tgrd_max,idx_T_max)) if(masterproc) write(iulog,*) 'Attempting to read snow aging parameters .....' call getfil (fsnowaging, locfn, 0) call ncd_pio_openfile(ncid, locfn, 0) if(masterproc) write(iulog,*) subname,trim(fsnowaging) ! snow aging parameters call ncd_io('tau', snowage_tau, 'read', ncid, posNOTonfile=.true.) call ncd_io('kappa', snowage_kappa, 'read', ncid, posNOTonfile=.true.) call ncd_io('drdsdt0', snowage_drdt0, 'read', ncid, posNOTonfile=.true.) call ncd_pio_closefile(ncid) if (masterproc) then write(iulog,*) 'Successfully read snow aging properties' ! print some diagnostics: write (iulog,*) 'SNICAR: snowage tau for T=263K, dTdz = 100 K/m, rhos = 150 kg/m3: ', snowage_tau(3,11,9) write (iulog,*) 'SNICAR: snowage kappa for T=263K, dTdz = 100 K/m, rhos = 150 kg/m3: ', snowage_kappa(3,11,9) write (iulog,*) 'SNICAR: snowage dr/dt_0 for T=263K, dTdz = 100 K/m, rhos = 150 kg/m3: ', snowage_drdt0(3,11,9) endif end subroutine SnowAge_init end module SNICARMod \markboth{Left}{Source File: SnowHydrologyMod.F90, Date: Wed Jun 15 14:32:23 MDT 2011 } module SnowHydrologyMod ----------------------------------------------------------------------- %///////////////////////////////////////////////////////////// \mbox{}\hrulefill\ \subsection{Fortran: Module Interface SnowHydrologyMod (Source File: SnowHydrologyMod.F90)} Calculate snow hydrology. \bigskip{\em USES:} \begin{verbatim} use shr_kind_mod, only: r8 => shr_kind_r8 use clm_varpar , only : nlevsnoPUBLIC TYPES:
implicit none savePUBLIC MEMBER FUNCTIONS:
public :: SnowWater ! Change of snow mass and the snow water onto soil public :: SnowCompaction ! Change in snow layer thickness due to compaction public :: CombineSnowLayers ! Combine snow layers less than a min thickness public :: DivideSnowLayers ! Subdivide snow layers if they exceed maximum thickness public :: BuildSnowFilter ! Construct snow/no-snow filtersPRIVATE MEMBER FUNCTIONS:
private :: Combo ! Returns the combined variables: dz, t, wliq, wice.REVISION HISTORY:
Created by Mariana Vertenstein