Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

fixed stratosphere warm bias and code optimization for MERRA2 #88

Merged
merged 2 commits into from
Jul 24, 2023
Merged

fixed stratosphere warm bias and code optimization for MERRA2 #88

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
bb2af8c
optimized the code and decrease the stratosphere warm bias for mraero…
AnningCheng-NOAA Jun 29, 2023
a69c8a3
Merge branch 'ufs/dev' into strat_warm_bias_fix_cheng
grantfirl Jul 19, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
7 changes: 1 addition & 6 deletions physics/GFS_rrtmgp_cloud_mp.F90
View file
Open in desktop
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
20 changes: 12 additions & 8 deletions physics/module_mp_thompson.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3400,8 +3400,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 +3654,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 +3744,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 +4234,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
View file
Open in desktop
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
Loading