Merge pull request #4 from ericaligo-NOAA/winterwx

Add rhonewsn1 to nuopc and to the NOAH MP driver and source code

drivers/ccpp/noahmpdrv.F90

@@ -138,7 +138,7 @@ subroutine noahmpdrv_run                                       &
       idveg, iopt_crs, iopt_btr, iopt_run, iopt_sfc, iopt_frz,   &
       iopt_inf, iopt_rad, iopt_alb, iopt_snf, iopt_tbot,         &
       iopt_stc, iopt_trs,xlatin, xcoszin, iyrlen, julian, garea, &
-      rainn_mp, rainc_mp, snow_mp, graupel_mp, ice_mp,           &
+      rainn_mp, rainc_mp, snow_mp, graupel_mp, ice_mp,rhonewsn1, &
       con_hvap, con_cp, con_jcal, rhoh2o, con_eps, con_epsm1,    &
       con_fvirt, con_rd, con_hfus, thsfc_loc,                    &
 
@@ -262,6 +262,7 @@ subroutine noahmpdrv_run                                       &
   real(kind=kind_phys), dimension(:)     , intent(in)    :: snow_mp    ! microphysics snow [mm]
   real(kind=kind_phys), dimension(:)     , intent(in)    :: graupel_mp ! microphysics graupel [mm]
   real(kind=kind_phys), dimension(:)     , intent(in)    :: ice_mp     ! microphysics ice/hail [mm]
+  real(kind=kind_phys), dimension(:)     , intent(in)    :: rhonewsn1  ! precipitation ice density (kg/m^3)
   real(kind=kind_phys)                   , intent(in)    :: con_hvap   ! latent heat condensation [J/kg]
   real(kind=kind_phys)                   , intent(in)    :: con_cp     ! specific heat air [J/kg/K] 
   real(kind=kind_phys)                   , intent(in)    :: con_jcal   ! joules per calorie (not used)
@@ -757,7 +758,7 @@ subroutine noahmpdrv_run                                       &
 
       call transfer_mp_parameters(vegetation_category, soil_category, &
                         slope_category, soil_color_category, crop_type,parameters)
-
+      parameters%prcpiceden = rhonewsn1(i)
       call noahmp_options(idveg ,iopt_crs, iopt_btr , iopt_run, iopt_sfc,  &
                                  iopt_frz, iopt_inf , iopt_rad, iopt_alb,  &
                                  iopt_snf, iopt_tbot, iopt_stc, iopt_rsf,  &

drivers/ccpp/noahmpdrv.meta

@@ -533,6 +533,14 @@
   type = real
   kind = kind_phys
   intent = in
+[rhonewsn1]
+  standard_name = lwe_density_of_precip_ice
+  long_name = density of precipitation ice
+  units = kg m-3
+  dimensions = (horizontal_loop_extent)
+  type = real
+  kind = kind_phys
+  intent = in
 [con_hvap]
   standard_name = latent_heat_of_vaporization_of_water_at_0C
   long_name = latent heat of evaporation/sublimation

drivers/nuopc/lnd_comp_driver.F90

@@ -574,6 +574,7 @@ subroutine drv_run(gcomp, noahmp, rc)
          noahmp%model%julian   , noahmp%domain%garea    , &
          noahmp%model%rainn_mp , noahmp%model%rainc_mp  , &
          noahmp%model%snow_mp  , noahmp%model%graupel_mp, noahmp%model%ice_mp    , &
+         noahmp%model%rhonewsn1, &
          con_hvap              , con_cp                 , con_jcal               , &
          rhoh2o                , con_eps                , con_epsm1              , &
          con_fvirt             , con_rd                 , con_hfus               , &

drivers/nuopc/lnd_comp_types.F90

@@ -190,6 +190,7 @@ module lnd_comp_types
      real(kind=kp), allocatable :: edir       (:)   ! soil surface evaporation (mm/s)
      real(kind=kp), allocatable :: snowc      (:)   ! fractional snow cover 
      real(kind=kp), allocatable :: stm        (:)   ! total soil column moisture content (m)
+     real(kind=kp), allocatable :: rhonewsn1  (:)   ! precipitation ice density (kg/m^3) 
      real(kind=kp), allocatable :: snohf      (:)   ! snow/freezing-rain latent heat flux (W/m^2)
      real(kind=kp), allocatable :: smcwlt2    (:)   ! dry soil moisture threshold
      real(kind=kp), allocatable :: smcref2    (:)   ! soil moisture threshold
@@ -469,6 +470,7 @@ subroutine InitializeAllocate(this, begl, endl, km, lsnowl)
     allocate(this%model%edir       (begl:endl))
     allocate(this%model%snowc      (begl:endl))
     allocate(this%model%stm        (begl:endl))
+    allocate(this%model%rhonewsn1  (begl:endl))
     allocate(this%model%snohf      (begl:endl))
     allocate(this%model%smcwlt2    (begl:endl))
     allocate(this%model%smcref2    (begl:endl))
@@ -635,6 +637,7 @@ subroutine InitializeDefault(this)
     this%model%edir        = 0.0_kp
     this%model%snowc       = 0.0_kp
     this%model%stm         = 0.0_kp
+    this%model%rhonewsn1   = 0.0_kp
     this%model%snohf       = 0.0_kp
     this%model%smcwlt2     = 0.0_kp
     this%model%smcref2     = 0.0_kp

src/module_sf_noahmplsm.f90

@@ -218,6 +218,7 @@ module module_sf_noahmplsm
     real (kind=kind_phys) :: saim(12)           !< monthly stem area index, one-sided
     real (kind=kind_phys) :: laim(12)           !< monthly leaf area index, one-sided
     real (kind=kind_phys) :: sla                !< single-side leaf area per kg [m2/kg]
+    real (kind=kind_phys) :: prcpiceden         !< precipitation ice density [kg/m^3]
     real (kind=kind_phys) :: dilefc             !< coeficient for leaf stress death [1/s]
     real (kind=kind_phys) :: dilefw             !< coeficient for leaf stress death [1/s]
     real (kind=kind_phys) :: fragr              !< fraction of growth respiration  !original was 0.3 
@@ -1067,13 +1068,14 @@ subroutine atm (parameters,sfcprs  ,sfctmp   ,q2      ,
 ! fresh snow density
 
      bdfall = min(120.,67.92+51.25*exp((sfctmp-tfrz)/2.59))       !mb/an: change to min  
-     if(opt_snf == 4) then
+     if(opt_snf == 4 .or. opt_snf == 5) then
         prcp_frozen = prcpsnow + prcpgrpl + prcphail
         if(prcpnonc > 0. .and. prcp_frozen > 0.) then
 	  fpice = min(1.0,prcp_frozen/prcpnonc)
 	  fpice = max(0.0,fpice)
-	  bdfall = bdfall*(prcpsnow/prcp_frozen) + rho_grpl*(prcpgrpl/prcp_frozen) + &
+	  if(opt_snf==4) bdfall = bdfall*(prcpsnow/prcp_frozen) + rho_grpl*(prcpgrpl/prcp_frozen) + &
 	             rho_hail*(prcphail/prcp_frozen)
+          if(opt_snf==5) bdfall = parameters%prcpiceden
 	else
 	  fpice = 0.0
         endif