Skip to content

Commit

Permalink
Merge branch 'ufs/dev' into srw300_scidoc_noimages
Browse files Browse the repository at this point in the history
  • Loading branch information
grantfirl committed Aug 18, 2023
2 parents 816f607 + af890d4 commit e2818ec
Show file tree
Hide file tree
Showing 3 changed files with 46 additions and 47 deletions.
7 changes: 1 addition & 6 deletions physics/GFS_rrtmgp_cloud_mp.F90
Original file line number Diff line number Diff line change
Expand Up @@ -875,17 +875,12 @@ subroutine cmp_reff_Thompson(nLev, nCol, i_cldliq, i_cldice, i_cldsnow, i_cldice
qi_mp(iCol,iLay) = tracer(iCol,iLay,i_cldice) / (1.-q_lay(iCol,iLay))
qs_mp(iCol,iLay) = tracer(iCol,iLay,i_cldsnow) / (1.-q_lay(iCol,iLay))
ni_mp(iCol,iLay) = tracer(iCol,iLay,i_cldice_nc) / (1.-q_lay(iCol,iLay))
if (ltaerosol) then
if (ltaerosol .or. mraerosol) then
nc_mp(iCol,iLay) = tracer(iCol,iLay,i_cldliq_nc) / (1.-q_lay(iCol,iLay))
nwfa(iCol,iLay) = tracer(iCol,iLay,i_twa)
if (qc_mp(iCol,iLay) > 1.e-12 .and. nc_mp(iCol,iLay) < 100.) then
nc_mp(iCol,iLay) = make_DropletNumber(qc_mp(iCol,iLay)*rho, nwfa(iCol,iLay)*rho) * orho
endif
elseif (mraerosol) then
nc_mp(iCol,iLay) = tracer(iCol,iLay,i_cldliq_nc) / (1.-q_lay(iCol,iLay))
if (qc_mp(iCol,iLay) > 1.e-12 .and. nc_mp(iCol,iLay) < 100.) then
nc_mp(iCol,iLay) = make_DropletNumber(qc_mp(iCol,iLay)*rho, nwfa(iCol,iLay)*rho) * orho
endif
else
if (nint(lsmask(iCol)) == 1) then !land
nc_mp(iCol,iLay) = nt_c_l*orho
Expand Down
66 changes: 34 additions & 32 deletions physics/module_mp_thompson.F90
Original file line number Diff line number Diff line change
Expand Up @@ -708,9 +708,9 @@ SUBROUTINE thompson_init(is_aerosol_aware_in, &
dtc(n) = (Dc(n) - Dc(n-1))
enddo

!> - Create bins of cloud ice (from min diameter up to 5x min snow size)
!> - Create bins of cloud ice (from min diameter up to 2x min snow size)
xDx(1) = D0i*1.0d0
xDx(nbi+1) = 5.0d0*D0s
xDx(nbi+1) = 2.0d0*D0s
do n = 2, nbi
xDx(n) = DEXP(DFLOAT(n-1)/DFLOAT(nbi) &
*DLOG(xDx(nbi+1)/xDx(1)) +DLOG(xDx(1)))
Expand Down Expand Up @@ -2822,7 +2822,7 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, &
prr_rcg(k) = MIN(DBLE(rg(k)*odts), prr_rcg(k))
prg_rcg(k) = -prr_rcg(k)
!> - Put in explicit drop break-up due to collisions.
pnr_rcg(k) = -5.*tnr_gacr(idx_g1,idx_g,idx_r1,idx_r) ! RAIN2M
pnr_rcg(k) = -1.5*tnr_gacr(idx_g1,idx_g,idx_r1,idx_r) ! RAIN2M
endif
endif
endif
Expand Down Expand Up @@ -3053,34 +3053,32 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, &
if (prr_sml(k) .gt. 0.) then
prr_sml(k) = prr_sml(k) + 4218.*olfus*tempc &
* (prr_rcs(k)+prs_scw(k))
endif
prr_sml(k) = MIN(DBLE(rs(k)*odts), MAX(0.D0, prr_sml(k)))
pnr_sml(k) = smo0(k)/rs(k)*prr_sml(k) * 10.0**(-0.25*tempc) ! RAIN2M
pnr_sml(k) = MIN(DBLE(smo0(k)*odts), pnr_sml(k))

if (ssati(k).lt. 0.) then
prs_sde(k) = C_cube*t1_subl*diffu(k)*ssati(k)*rvs &
* (t1_qs_sd*smo1(k) &
+ t2_qs_sd*rhof2(k)*vsc2(k)*smof(k))
prs_sde(k) = MAX(DBLE(-rs(k)*odts), prs_sde(k))
prr_sml(k) = MIN(DBLE(rs(k)*odts), prr_sml(k))
pnr_sml(k) = smo0(k)/rs(k)*prr_sml(k) * 10.0**(-0.25*tempc) ! RAIN2M
pnr_sml(k) = MIN(DBLE(smo0(k)*odts), pnr_sml(k))
elseif (ssati(k).lt. 0.) then
prr_sml(k) = 0.0
prs_sde(k) = C_cube*t1_subl*diffu(k)*ssati(k)*rvs &
* (t1_qs_sd*smo1(k) &
+ t2_qs_sd*rhof2(k)*vsc2(k)*smof(k))
prs_sde(k) = MAX(DBLE(-rs(k)*odts), prs_sde(k))
endif
endif

if (L_qg(k)) then
prr_gml(k) = (tempc*tcond(k)-lvap0*diffu(k)*delQvs(k)) &
* N0_g(k)*(t1_qg_me*ilamg(k)**cge(10) &
+ t2_qg_me*rhof2(k)*vsc2(k)*ilamg(k)**cge(11))
!-GT prr_gml(k) = prr_gml(k) + 4218.*olfus*tempc &
!-GT * (prr_rcg(k)+prg_gcw(k))
prr_gml(k) = MIN(DBLE(rg(k)*odts), MAX(0.D0, prr_gml(k)))
pnr_gml(k) = N0_g(k)*cgg(2)*ilamg(k)**cge(2) / rg(k) & ! RAIN2M
* prr_gml(k) * 10.0**(-0.5*tempc)

if (ssati(k).lt. 0.) then
prg_gde(k) = C_cube*t1_subl*diffu(k)*ssati(k)*rvs &
* N0_g(k) * (t1_qg_sd*ilamg(k)**cge(10) &
+ t2_qg_sd*vsc2(k)*rhof2(k)*ilamg(k)**cge(11))
prg_gde(k) = MAX(DBLE(-rg(k)*odts), prg_gde(k))
if (prr_gml(k) .gt. 0.) then
prr_gml(k) = MIN(DBLE(rg(k)*odts), prr_gml(k))
pnr_gml(k) = N0_g(k)*cgg(2)*ilamg(k)**cge(2) / rg(k) & ! RAIN2M
* prr_gml(k) * 10.0**(-0.5*tempc)
elseif (ssati(k).lt. 0.) then
prr_gml(k) = 0.0
prg_gde(k) = C_cube*t1_subl*diffu(k)*ssati(k)*rvs &
* N0_g(k) * (t1_qg_sd*ilamg(k)**cge(10) &
+ t2_qg_sd*vsc2(k)*rhof2(k)*ilamg(k)**cge(11))
prg_gde(k) = MAX(DBLE(-rg(k)*odts), prg_gde(k))
endif
endif

Expand Down Expand Up @@ -3400,8 +3398,8 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, &
tcond(k) = (5.69 + 0.0168*tempc)*1.0E-5 * 418.936
ocp(k) = 1./(Cp*(1.+0.887*qv(k)))
lvt2(k)=lvap(k)*lvap(k)*ocp(k)*oRv*otemp*otemp

nwfa(k) = MAX(11.1E6*rho(k), (nwfa1d(k) + nwfaten(k)*DT)*rho(k))
if (is_aerosol_aware) &
nwfa(k) = MAX(11.1E6*rho(k), (nwfa1d(k) + nwfaten(k)*DT)*rho(k))
enddo

do k = kts, kte
Expand Down Expand Up @@ -3654,7 +3652,8 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, &
qvten(k) = qvten(k) - prw_vcd(k)
qcten(k) = qcten(k) + prw_vcd(k)
ncten(k) = ncten(k) + pnc_wcd(k)
nwfaten(k) = nwfaten(k) - pnc_wcd(k)
if (is_aerosol_aware) &
nwfaten(k) = nwfaten(k) - pnc_wcd(k)
tten(k) = tten(k) + lvap(k)*ocp(k)*prw_vcd(k)*(1-IFDRY)
rc(k) = MAX(R1, (qc1d(k) + DT*qcten(k))*rho(k))
if (rc(k).eq.R1) L_qc(k) = .false.
Expand Down Expand Up @@ -3743,7 +3742,8 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, &
qrten(k) = qrten(k) - prv_rev(k)
qvten(k) = qvten(k) + prv_rev(k)
nrten(k) = nrten(k) - pnr_rev(k)
nwfaten(k) = nwfaten(k) + pnr_rev(k)
if (is_aerosol_aware) &
nwfaten(k) = nwfaten(k) + pnr_rev(k)
tten(k) = tten(k) - lvap(k)*ocp(k)*prv_rev(k)*(1-IFDRY)

rr(k) = MAX(R1, (qr1d(k) + DT*qrten(k))*rho(k))
Expand Down Expand Up @@ -4232,10 +4232,12 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, &
qv1d(k) = MAX(1.E-10, qv1d(k) + qvten(k)*DT)
qc1d(k) = qc1d(k) + qcten(k)*DT
nc1d(k) = MAX(2./rho(k), MIN(nc1d(k) + ncten(k)*DT, Nt_c_max))
nwfa1d(k) = MAX(11.1E6, MIN(9999.E6, &
(nwfa1d(k)+nwfaten(k)*DT)))
nifa1d(k) = MAX(naIN1*0.01, MIN(9999.E6, &
(nifa1d(k)+nifaten(k)*DT)))
if (is_aerosol_aware) then
nwfa1d(k) = MAX(11.1E6, MIN(9999.E6, &
(nwfa1d(k)+nwfaten(k)*DT)))
nifa1d(k) = MAX(naIN1*0.01, MIN(9999.E6, &
(nifa1d(k)+nifaten(k)*DT)))
end if
if (qc1d(k) .le. R1) then
qc1d(k) = 0.0
nc1d(k) = 0.0
Expand Down
20 changes: 11 additions & 9 deletions physics/mp_thompson.F90
Original file line number Diff line number Diff line change
Expand Up @@ -151,6 +151,10 @@ subroutine mp_thompson_init(ncol, nlev, con_g, con_rd, con_eps, &
!> - Convert specific humidity to water vapor mixing ratio.
!> - Also, hydrometeor variables are mass or number mixing ratio
!> - either kg of species per kg of dry air, or per kg of (dry + vapor).
if (merra2_aerosol_aware) then
call get_niwfa(aerfld, nifa, nwfa, ncol, nlev)
end if


qv = spechum/(1.0_kind_phys-spechum)

Expand All @@ -163,7 +167,7 @@ subroutine mp_thompson_init(ncol, nlev, con_g, con_rd, con_eps, &

ni = ni/(1.0_kind_phys-spechum)
nr = nr/(1.0_kind_phys-spechum)
if (is_aerosol_aware) then
if (is_aerosol_aware .or. merra2_aerosol_aware) then
nc = nc/(1.0_kind_phys-spechum)
nwfa = nwfa/(1.0_kind_phys-spechum)
nifa = nifa/(1.0_kind_phys-spechum)
Expand Down Expand Up @@ -208,8 +212,6 @@ subroutine mp_thompson_init(ncol, nlev, con_g, con_rd, con_eps, &
nwfa(i,k) = naCCN1+naCCN0*exp(-((hgt(i,k)-hgt(i,1))/1000.)*niCCN3)
enddo
enddo
else if (merra2_aerosol_aware) then
call get_niwfa(aerfld, nifa, nwfa, ncol, nlev)
else
if (mpirank==mpiroot) write(*,*) ' Apparently initial CCN aerosols are present.'
if (MAXVAL(nwfa2d) .lt. eps) then
Expand Down Expand Up @@ -555,6 +557,9 @@ subroutine mp_thompson_run(ncol, nlev, con_g, con_rd, &
else
dtstep = dtp
end if
if (merra2_aerosol_aware) then
call get_niwfa(aerfld, nifa, nwfa, ncol, nlev)
end if

!> - Convert specific humidity to water vapor mixing ratio.
!> - Also, hydrometeor variables are mass or number mixing ratio
Expand All @@ -574,7 +579,7 @@ subroutine mp_thompson_run(ncol, nlev, con_g, con_rd, &

ni = ni/(1.0_kind_phys-spechum)
nr = nr/(1.0_kind_phys-spechum)
if (is_aerosol_aware) then
if (is_aerosol_aware .or. merra2_aerosol_aware) then
nc = nc/(1.0_kind_phys-spechum)
nwfa = nwfa/(1.0_kind_phys-spechum)
nifa = nifa/(1.0_kind_phys-spechum)
Expand Down Expand Up @@ -681,9 +686,6 @@ subroutine mp_thompson_run(ncol, nlev, con_g, con_rd, &
ncten3 => diag3d(:,:,36:36)
qcten3 => diag3d(:,:,37:37)
end if set_extended_diagnostic_pointers
if (merra2_aerosol_aware) then
call get_niwfa(aerfld, nifa, nwfa, ncol, nlev)
end if
!> - Call mp_gt_driver() with or without aerosols, with or without effective radii, ...
if (is_aerosol_aware .or. merra2_aerosol_aware) then
call mp_gt_driver(qv=qv, qc=qc, qr=qr, qi=qi, qs=qs, qg=qg, ni=ni, nr=nr, &
Expand Down Expand Up @@ -921,8 +923,8 @@ subroutine get_niwfa(aerfld, nifa, nwfa, ncol, nlev)
aerfld(:,:,4)/1011.5142+ aerfld(:,:,5)/5683.3501)*1.e15

nwfa=((aerfld(:,:,6)/0.0045435214+aerfld(:,:,7)/0.2907854+aerfld(:,:,8)/12.91224+ &
aerfld(:,:,9)/206.2216+ aerfld(:,:,10)/4326.23)*1.+aerfld(:,:,11)/0.3053104*5+ &
aerfld(:,:,15)/0.3232698*1)*1.e15
aerfld(:,:,9)/206.2216+ aerfld(:,:,10)/4326.23)*9.+aerfld(:,:,11)/0.3053104*5+ &
aerfld(:,:,15)/0.3232698*8)*1.e15
end subroutine get_niwfa

end module mp_thompson

0 comments on commit e2818ec

Please sign in to comment.