diff --git a/physics/GFS_radiation_surface.F90 b/physics/GFS_radiation_surface.F90 index 2481af163..7da6e5df7 100644 --- a/physics/GFS_radiation_surface.F90 +++ b/physics/GFS_radiation_surface.F90 @@ -88,7 +88,7 @@ subroutine GFS_radiation_surface_run ( & alnwf, facsf, facwf, & semis_lnd, semis_ice, snoalb character(len=3) , dimension(:), intent(in) :: lndp_var_list - real(kind=kind_phys), dimension(:), intent(in) :: albdvis_lnd, albdnir_lnd, & + real(kind=kind_phys), dimension(:), intent(inout) :: albdvis_lnd, albdnir_lnd, & albivis_lnd, albinir_lnd real(kind=kind_phys), dimension(:), intent(in) :: albdvis_ice, albdnir_ice, & albivis_ice, albinir_ice diff --git a/physics/GFS_radiation_surface.meta b/physics/GFS_radiation_surface.meta index c38ffe2a3..82faaa865 100644 --- a/physics/GFS_radiation_surface.meta +++ b/physics/GFS_radiation_surface.meta @@ -411,7 +411,7 @@ dimensions = (horizontal_loop_extent) type = real kind = kind_phys - intent = in + intent = inout optional = F [albdnir_lnd] standard_name = surface_albedo_direct_NIR_over_land @@ -420,7 +420,7 @@ dimensions = (horizontal_loop_extent) type = real kind = kind_phys - intent = in + intent = inout optional = F [albivis_lnd] standard_name = surface_albedo_diffuse_visible_over_land @@ -429,7 +429,7 @@ dimensions = (horizontal_loop_extent) type = real kind = kind_phys - intent = in + intent = inout optional = F [albinir_lnd] standard_name = surface_albedo_diffuse_NIR_over_land @@ -438,7 +438,7 @@ dimensions = (horizontal_loop_extent) type = real kind = kind_phys - intent = in + intent = inout optional = F [albdvis_ice] standard_name = surface_albedo_direct_visible_over_ice diff --git a/physics/module_mp_thompson.F90 b/physics/module_mp_thompson.F90 index b1301d744..f05aa8ba2 100644 --- a/physics/module_mp_thompson.F90 +++ b/physics/module_mp_thompson.F90 @@ -969,7 +969,7 @@ END SUBROUTINE thompson_init SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, & nwfa, nifa, nwfa2d, nifa2d, & tt, th, pii, & - p, w, dz, dt_in, & + p, w, dz, dt_in, dt_inner, & RAINNC, RAINNCV, & SNOWNC, SNOWNCV, & ICENC, ICENCV, & @@ -1046,7 +1046,7 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, & REAL, DIMENSION(ims:ime, kms:kme, jms:jme), OPTIONAL, INTENT(INOUT):: & vt_dbz_wt LOGICAL, INTENT(IN) :: first_time_step - REAL, INTENT(IN):: dt_in + REAL, INTENT(IN):: dt_in, dt_inner ! To support subcycling: current step and maximum number of steps INTEGER, INTENT (IN) :: istep, nsteps LOGICAL, INTENT (IN) :: reset_dBZ @@ -1107,6 +1107,7 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, & LOGICAL, OPTIONAL, INTENT(IN) :: diagflag INTEGER, OPTIONAL, INTENT(IN) :: do_radar_ref logical :: melti = .false. + INTEGER :: ndt, it ! CCPP error handling character(len=*), optional, intent( out) :: errmsg @@ -1236,7 +1237,30 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, & ! j_end = jte ! endif - dt = dt_in +! dt = dt_in + RAINNC(:,:) = 0.0 + SNOWNC(:,:) = 0.0 + ICENC(:,:) = 0.0 + GRAUPELNC(:,:) = 0.0 + pcp_ra(:,:) = 0.0 + pcp_sn(:,:) = 0.0 + pcp_gr(:,:) = 0.0 + pcp_ic(:,:) = 0.0 + ndt = max(nint(dt_in/dt_inner),1) + dt = dt_in/ndt + if(dt_in .le. dt_inner) dt= dt_in + if(nsteps>1 .and. ndt>1) then + if (present(errmsg) .and. present(errflg)) then + write(errmsg, '(a)') 'Logic error in mp_gt_driver: inner loop cannot be used with subcycling' + errflg = 1 + return + else + write(*,'(a)') 'Warning: inner loop cannot be used with subcycling, resetting ndt=1' + ndt = 1 + endif + endif + + do it = 1, ndt qc_max = 0. qr_max = 0. @@ -1412,10 +1436,10 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, & tten1, qvten1, qrten1, qsten1, & qgten1, qiten1, niten1, nrten1, ncten1, qcten1) - pcp_ra(i,j) = pptrain - pcp_sn(i,j) = pptsnow - pcp_gr(i,j) = pptgraul - pcp_ic(i,j) = pptice + pcp_ra(i,j) = pcp_ra(i,j) + pptrain + pcp_sn(i,j) = pcp_sn(i,j) + pptsnow + pcp_gr(i,j) = pcp_gr(i,j) + pptgraul + pcp_ic(i,j) = pcp_ic(i,j) + pptice RAINNCV(i,j) = pptrain + pptsnow + pptgraul + pptice RAINNC(i,j) = RAINNC(i,j) + pptrain + pptsnow + pptgraul + pptice IF ( PRESENT(snowncv) .AND. PRESENT(snownc) ) THEN @@ -1594,9 +1618,26 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, & enddo endif assign_extended_diagnostics + if (ndt>1 .and. it==ndt) then + + SR(i,j) = (pcp_sn(i,j) + pcp_gr(i,j) + pcp_ic(i,j))/(RAINNC(i,j)+1.e-12) + RAINNCV(i,j) = RAINNC(i,j) + IF ( PRESENT (snowncv) ) THEN + SNOWNCV(i,j) = SNOWNC(i,j) + ENDIF + IF ( PRESENT (icencv) ) THEN + ICENCV(i,j) = ICENC(i,j) + ENDIF + IF ( PRESENT (graupelncv) ) THEN + GRAUPELNCV(i,j) = GRAUPELNC(i,j) + ENDIF + endif + ! Diagnostic calculations only for last step ! if Thompson MP is called multiple times - last_step_only: IF (istep == nsteps) THEN + last_step_only: IF ((ndt>1 .and. it==ndt) .or. & + (nsteps>1 .and. istep==nsteps) .or. & + (nsteps==1 .and. ndt==1)) THEN !> - Call calc_refl10cm() @@ -1656,6 +1697,7 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, & ! 'ni: ', ni_max, '(', imax_ni, ',', jmax_ni, ',', kmax_ni, ')', & ! 'nr: ', nr_max, '(', imax_nr, ',', jmax_nr, ',', kmax_nr, ')' ! END DEBUG - GT + enddo ! end of nt loop ! These are always allocated !deallocate (vtsk1) diff --git a/physics/module_sf_noahmplsm.f90 b/physics/module_sf_noahmplsm.f90 index e9e17f399..9fcb7edf8 100644 --- a/physics/module_sf_noahmplsm.f90 +++ b/physics/module_sf_noahmplsm.f90 @@ -1904,6 +1904,10 @@ subroutine energy (parameters,ice ,vegtyp ,ist ,nsnow ,nsoil , & !in ! ground snow cover fraction [niu and yang, 2007, jgr] fsno = 0. + if(snowh <= 1.e-6 .or. sneqv <= 1.e-3) then + snowh = 0.0 + sneqv = 0.0 + end if if(snowh.gt.0.) then bdsno = sneqv / snowh fmelt = (bdsno/100.)**parameters%mfsno @@ -7070,7 +7074,7 @@ subroutine snowh2o (parameters,nsnow ,nsoil ,dt ,qsnfro ,qsnsub , & !in end if end if - if(snowh <= 1.e-8 .or. sneqv <= 1.e-6) then + if(snowh <= 1.e-6 .or. sneqv <= 1.e-3) then snowh = 0.0 sneqv = 0.0 end if diff --git a/physics/mp_thompson.F90 b/physics/mp_thompson.F90 index 6fb039b9d..c31d90b09 100644 --- a/physics/mp_thompson.F90 +++ b/physics/mp_thompson.F90 @@ -334,7 +334,7 @@ subroutine mp_thompson_run(ncol, nlev, con_g, con_rd, & spechum, qc, qr, qi, qs, qg, ni, nr, & is_aerosol_aware, nc, nwfa, nifa, & nwfa2d, nifa2d, & - tgrs, prsl, phii, omega, & + tgrs, prsl, phii, omega, dt_inner, & dtp, first_time_step, istep, nsteps, & prcp, rain, graupel, ice, snow, sr, & refl_10cm, reset_dBZ, do_radar_ref, & @@ -389,6 +389,7 @@ subroutine mp_thompson_run(ncol, nlev, con_g, con_rd, & ! Radar reflectivity real(kind_phys), intent( out) :: refl_10cm(:,:) logical, optional, intent(in ) :: do_radar_ref + real, intent(in ) :: dt_inner ! Cloud effective radii real(kind_phys), optional, intent( out) :: re_cloud(:,:) real(kind_phys), optional, intent( out) :: re_ice(:,:) @@ -673,7 +674,7 @@ subroutine mp_thompson_run(ncol, nlev, con_g, con_rd, & if (do_effective_radii) then call mp_gt_driver(qv=qv, qc=qc, qr=qr, qi=qi, qs=qs, qg=qg, ni=ni, nr=nr, & nc=nc, nwfa=nwfa, nifa=nifa, nwfa2d=nwfa2d, nifa2d=nifa2d, & - tt=tgrs, p=prsl, w=w, dz=dz, dt_in=dtp, & + tt=tgrs, p=prsl, w=w, dz=dz, dt_in=dtstep, dt_inner=dt_inner, & rainnc=rain_mp, rainncv=delta_rain_mp, & snownc=snow_mp, snowncv=delta_snow_mp, & icenc=ice_mp, icencv=delta_ice_mp, & @@ -713,7 +714,7 @@ subroutine mp_thompson_run(ncol, nlev, con_g, con_rd, & else call mp_gt_driver(qv=qv, qc=qc, qr=qr, qi=qi, qs=qs, qg=qg, ni=ni, nr=nr, & nc=nc, nwfa=nwfa, nifa=nifa, nwfa2d=nwfa2d, nifa2d=nifa2d, & - tt=tgrs, p=prsl, w=w, dz=dz, dt_in=dtp, & + tt=tgrs, p=prsl, w=w, dz=dz, dt_in=dtstep, dt_inner=dt_inner, & rainnc=rain_mp, rainncv=delta_rain_mp, & snownc=snow_mp, snowncv=delta_snow_mp, & icenc=ice_mp, icencv=delta_ice_mp, & @@ -753,7 +754,7 @@ subroutine mp_thompson_run(ncol, nlev, con_g, con_rd, & else if (do_effective_radii) then call mp_gt_driver(qv=qv, qc=qc, qr=qr, qi=qi, qs=qs, qg=qg, ni=ni, nr=nr, & - tt=tgrs, p=prsl, w=w, dz=dz, dt_in=dtp, & + tt=tgrs, p=prsl, w=w, dz=dz, dt_in=dtstep, dt_inner=dt_inner, & rainnc=rain_mp, rainncv=delta_rain_mp, & snownc=snow_mp, snowncv=delta_snow_mp, & icenc=ice_mp, icencv=delta_ice_mp, & @@ -792,7 +793,7 @@ subroutine mp_thompson_run(ncol, nlev, con_g, con_rd, & qcten3=qcten3) else call mp_gt_driver(qv=qv, qc=qc, qr=qr, qi=qi, qs=qs, qg=qg, ni=ni, nr=nr, & - tt=tgrs, p=prsl, w=w, dz=dz, dt_in=dtp, & + tt=tgrs, p=prsl, w=w, dz=dz, dt_in=dtstep, dt_inner=dt_inner, & rainnc=rain_mp, rainncv=delta_rain_mp, & snownc=snow_mp, snowncv=delta_snow_mp, & icenc=ice_mp, icencv=delta_ice_mp, & diff --git a/physics/mp_thompson.meta b/physics/mp_thompson.meta index 1ab496c25..15605fcbe 100644 --- a/physics/mp_thompson.meta +++ b/physics/mp_thompson.meta @@ -562,6 +562,15 @@ kind = kind_phys intent = in optional = F +[dt_inner] + standard_name = time_step_for_inner_loop + long_name = time step for inner loop + units = s + dimensions = () + type = real + kind = kind_phys + intent = in + optional = F [first_time_step] standard_name = flag_for_first_time_step long_name = flag for first time step for time integration loop (cold/warmstart) diff --git a/physics/radiation_surface.f b/physics/radiation_surface.f index 3ec34513c..78525c1de 100644 --- a/physics/radiation_surface.f +++ b/physics/radiation_surface.f @@ -411,13 +411,15 @@ subroutine setalb & real (kind=kind_phys), dimension(:), intent(in) :: & & slmsk, snowf, zorlf, coszf, tsknf, tairf, hprif, & & alvsf, alnsf, alvwf, alnwf, facsf, facwf, fice, tisfc, & - & lsmalbdvis, lsmalbdnir, lsmalbivis, lsmalbinir, & & icealbdvis, icealbdnir, icealbivis, icealbinir, & & sncovr, sncovr_ice, snoalb, albPpert ! sfc-perts, mgehne real (kind=kind_phys), intent(in) :: pertalb ! sfc-perts, mgehne real (kind=kind_phys), intent(in) :: min_seaice real (kind=kind_phys), dimension(:), intent(in) :: & & fracl, fraco, fraci + real (kind=kind_phys), dimension(:),intent(inout) :: & + & lsmalbdvis, lsmalbdnir, lsmalbivis, lsmalbinir + logical, dimension(:), intent(in) :: & & icy @@ -537,6 +539,10 @@ subroutine setalb & alndnd = alnwf(i)*flnd + snoalb(i) * fsno alndvb = ab2bm *flnd + snoalb(i) * fsno alndvd = alvwf(i)*flnd + snoalb(i) * fsno + lsmalbdnir(i) = min(0.99,max(0.01,alndnb)) + lsmalbinir(i) = min(0.99,max(0.01,alndnd)) + lsmalbdvis(i) = min(0.99,max(0.01,alndvb)) + lsmalbivis(i) = min(0.99,max(0.01,alndvd)) else !-- fill in values for land albedo alndnb = 0. diff --git a/physics/sfc_drv.f b/physics/sfc_drv.f index 2ecb26469..817897fe7 100644 --- a/physics/sfc_drv.f +++ b/physics/sfc_drv.f @@ -94,10 +94,12 @@ end subroutine lsm_noah_finalize ! call sfc_drv ! ! --- inputs: ! ! ( im, km, ps, t1, q1, soiltyp, vegtype, sigmaf, ! -! sfcemis, dlwflx, dswsfc, snet, delt, tg3, cm, ch, ! +! sfcemis, dlwflx, dswsfc, delt, tg3, cm, ch, ! ! prsl1, prslki, zf, land, wind, slopetyp, ! ! shdmin, shdmax, snoalb, sfalb, flag_iter, flag_guess, ! ! lheatstrg, isot, ivegsrc, ! +! albdvis_lnd, albdnir_lnd, albivis_lnd, albinir_lnd, ! +! adjvisbmd, adjnirbmd, adjvisdfd, adjnirdfd, ! ! --- in/outs: ! ! weasd, snwdph, tskin, tprcp, srflag, smc, stc, slc, ! ! canopy, trans, tsurf, zorl, ! @@ -133,7 +135,6 @@ end subroutine lsm_noah_finalize ! sfcemis - real, sfc lw emissivity ( fraction ) im ! ! dlwflx - real, total sky sfc downward lw flux ( w/m**2 ) im ! ! dswflx - real, total sky sfc downward sw flux ( w/m**2 ) im ! -! snet - real, total sky sfc netsw flx into ground(w/m**2) im ! ! delt - real, time interval (second) 1 ! ! tg3 - real, deep soil temperature (k) im ! ! cm - real, surface exchange coeff for momentum (m/s) im ! @@ -148,6 +149,14 @@ end subroutine lsm_noah_finalize ! shdmax - real, max fractnl cover of green veg (not used) im ! ! snoalb - real, upper bound on max albedo over deep snow im ! ! sfalb - real, mean sfc diffused sw albedo (fractional) im ! +! albdvis_lnd - real, sfc vis direct sw albedo over land im ! +! albdnir_lnd - real, sfc nir direct sw albedo over land im ! +! albivis_lnd - real, sfc vis diffused sw albedo over land im ! +! albinir_lnd - real, sfc nir diffused sw albedo over land im ! +! adjvisbmd - real, t adj sfc uv+vis-beam sw dnward flux im ! +! adjnirbmd - real, t adj sfc nir-beam sw dnward flux im ! +! adjvisdfd - real, t adj sfc uv+vis-diff sw dnward flux im ! +! adjnirdfd - real, t adj sfc nir-diff sw dnward flux im ! ! flag_iter- logical, im ! ! flag_guess-logical, im ! ! lheatstrg- logical, flag for canopy heat storage 1 ! @@ -203,11 +212,13 @@ end subroutine lsm_noah_finalize subroutine lsm_noah_run & & ( im, km, grav, cp, hvap, rd, eps, epsm1, rvrdm1, ps, & ! --- inputs: & t1, q1, soiltyp, vegtype, sigmaf, & - & sfcemis, dlwflx, dswsfc, snet, delt, tg3, cm, ch, & + & sfcemis, dlwflx, dswsfc, delt, tg3, cm, ch, & & prsl1, prslki, zf, land, wind, slopetyp, & & shdmin, shdmax, snoalb, sfalb, flag_iter, flag_guess, & & lheatstrg, isot, ivegsrc, & & bexppert, xlaipert, vegfpert,pertvegf, & ! sfc perts, mgehne + & albdvis_lnd, albdnir_lnd, albivis_lnd, albinir_lnd, & + & adjvisbmd, adjnirbmd, adjvisdfd, adjnirdfd, & ! --- in/outs: & weasd, snwdph, tskin, tprcp, srflag, smc, stc, slc, & & canopy, trans, tsurf, zorl, & @@ -247,10 +258,12 @@ subroutine lsm_noah_run & integer, dimension(:), intent(in) :: soiltyp, vegtype, slopetyp real (kind=kind_phys), dimension(:), intent(in) :: ps, & - & t1, q1, sigmaf, sfcemis, dlwflx, dswsfc, snet, tg3, cm, & + & t1, q1, sigmaf, sfcemis, dlwflx, dswsfc, tg3, cm, & & ch, prsl1, prslki, wind, shdmin, shdmax, & & snoalb, sfalb, zf, & - & bexppert, xlaipert, vegfpert + & bexppert, xlaipert, vegfpert, & + & albdvis_lnd, albdnir_lnd, albivis_lnd, albinir_lnd, & + & adjvisbmd, adjnirbmd, adjvisdfd, adjnirdfd real (kind=kind_phys), intent(in) :: delt @@ -400,7 +413,12 @@ subroutine lsm_noah_run & lwdn = dlwflx(i) !..downward lw flux at sfc in w/m2 swdn = dswsfc(i) !..downward sw flux at sfc in w/m2 - solnet = snet(i) !..net sw rad flx (dn-up) at sfc in w/m2 +!net sw rad flx (dn-up) at sfc in w/m2 + solnet = adjvisbmd(i)*(1-albdvis_lnd(i)) & + & +adjnirbmd(i)*(1-albdnir_lnd(i)) & + & +adjvisdfd(i)*(1-albivis_lnd(i)) & + & +adjnirdfd(i)*(1-albinir_lnd(i)) + sfcems = sfcemis(i) sfcprs = prsl1(i) diff --git a/physics/sfc_drv.meta b/physics/sfc_drv.meta index c68102e7e..7812505c5 100644 --- a/physics/sfc_drv.meta +++ b/physics/sfc_drv.meta @@ -273,15 +273,6 @@ kind = kind_phys intent = in optional = F -[snet] - standard_name = surface_net_downwelling_shortwave_flux - long_name = total sky surface net shortwave flux - units = W m-2 - dimensions = (horizontal_loop_extent) - type = real - kind = kind_phys - intent = in - optional = F [delt] standard_name = time_step_for_dynamics long_name = dynamics time step @@ -482,6 +473,78 @@ kind = kind_phys intent = in optional = F +[albdvis_lnd] + standard_name = surface_albedo_direct_visible_over_land + long_name = direct surface albedo visible band over land + units = frac + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in + optional = F +[albdnir_lnd] + standard_name = surface_albedo_direct_NIR_over_land + long_name = direct surface albedo NIR band over land + units = frac + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in + optional = F +[albivis_lnd] + standard_name = surface_albedo_diffuse_visible_over_land + long_name = diffuse surface albedo visible band over land + units = frac + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in + optional = F +[albinir_lnd] + standard_name = surface_albedo_diffuse_NIR_over_land + long_name = diffuse surface albedo NIR band over land + units = frac + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in + optional = F +[adjvisbmd] + standard_name = surface_downwelling_direct_ultraviolet_and_visible_shortwave_flux + long_name = surface downwelling beam ultraviolet plus visible shortwave flux at current time + units = W m-2 + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in + optional = F +[adjnirbmd] + standard_name = surface_downwelling_direct_near_infrared_shortwave_flux + long_name = surface downwelling beam near-infrared shortwave flux at current time + units = W m-2 + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in + optional = F +[adjvisdfd] + standard_name = surface_downwelling_diffuse_ultraviolet_and_visible_shortwave_flux + long_name = surface downwelling diffuse ultraviolet plus visible shortwave flux at current time + units = W m-2 + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in + optional = F +[adjnirdfd] + standard_name = surface_downwelling_diffuse_near_infrared_shortwave_flux + long_name = surface downwelling diffuse near-infrared shortwave flux at current time + units = W m-2 + dimensions = (horizontal_loop_extent) + type = real + kind = kind_phys + intent = in + optional = F [weasd] standard_name = water_equivalent_accumulated_snow_depth_over_land long_name = water equiv of acc snow depth over land