Skip to content

Commit

Permalink
Merge branch 'main' of https://github.com/NCAR/ccpp-physics into ccpp…
Browse files Browse the repository at this point in the history 
…_error_code_in_prebuild
  • Loading branch information
@climbfuji
climbfuji committed Jan 21, 2022
2 parents d1527d5 + 899df7f commit f6e185d
Show file tree
Hide file tree
Showing 4 changed files with 95 additions and 117 deletions.
52 changes: 18 additions & 34 deletions physics/GFS_rrtmg_pre.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -183,7 +183,8 @@ subroutine GFS_rrtmg_pre_run (im, levs, lm, lmk, lmp, n_var_lndp, &

integer :: i, j, k, k1, k2, lsk, lv, n, itop, ibtc, LP1, lla, llb, lya,lyb

real(kind=kind_phys) :: es, qs, delt, tem0d, gridkm, pfac
real(kind=kind_phys) :: es, qs, delt, tem0d, pfac
real(kind=kind_phys), dimension(im) :: gridkm

real(kind=kind_phys), dimension(im) :: cvt1, cvb1, tem1d, tskn, xland

Expand All @@ -194,9 +195,9 @@ subroutine GFS_rrtmg_pre_run (im, levs, lm, lmk, lmp, n_var_lndp, &
effrl, effri, effrr, effrs, rho, orho, plyrpa

! for Thompson MP
real(kind=kind_phys), dimension(im,lm+LTP) :: &
re_cloud, re_ice, re_snow, qv_mp, qc_mp, &
qi_mp, qs_mp, nc_mp, ni_mp, nwfa
real(kind=kind_phys), dimension(im,lm+LTP) :: &
qv_mp, qc_mp, qi_mp, qs_mp, &
nc_mp, ni_mp, nwfa
real (kind=kind_phys), dimension(lm) :: cldfra1d, qv1d, &
& qc1d, qi1d, qs1d, dz1d, p1d, t1d

Expand Down Expand Up @@ -235,16 +236,14 @@ subroutine GFS_rrtmg_pre_run (im, levs, lm, lmk, lmp, n_var_lndp, &

LP1 = LM + 1 ! num of in/out levels


gridkm = sqrt(2.0)*sqrt(dx(1)*0.001*dx(1)*0.001)

if (imp_physics == imp_physics_thompson) then
max_relh = 1.5
else
max_relh = 1.1
endif

do i = 1, IM
gridkm(i) = dx(i)*0.001
lwp_ex(i) = 0.0
iwp_ex(i) = 0.0
lwp_fc(i) = 0.0
Expand Down Expand Up @@ -817,30 +816,18 @@ subroutine GFS_rrtmg_pre_run (im, levs, lm, lmk, lmp, n_var_lndp, &
! it will raise the low limit from 5 to 10, but the high limit will remain 125.
call calc_effectRad (tlyr(i,:), plyr(i,:)*100., qv_mp(i,:), qc_mp(i,:), &
nc_mp(i,:), qi_mp(i,:), ni_mp(i,:), qs_mp(i,:), &
re_cloud(i,:), re_ice(i,:), re_snow(i,:), 1, lm )
effrl(i,:), effri(i,:), effrs(i,:), 1, lm )
! Scale Thompson's effective radii from meter to micron
do k=1,lm
re_cloud(i,k) = MAX(re_qc_min, MIN(re_cloud(i,k), re_qc_max))
re_ice(i,k) = MAX(re_qi_min, MIN(re_ice(i,k), re_qi_max))
re_snow(i,k) = MAX(re_qs_min, MIN(re_snow(i,k), re_qs_max))
end do
end do
! Scale Thompson's effective radii from meter to micron
do k=1,lm
do i=1,im
re_cloud(i,k) = re_cloud(i,k)*1.e6
re_ice(i,k) = re_ice(i,k)*1.e6
re_snow(i,k) = re_snow(i,k)*1.e6
effrl(i,k) = MAX(re_qc_min, MIN(effrl(i,k), re_qc_max))*1.e6
effri(i,k) = MAX(re_qi_min, MIN(effri(i,k), re_qi_max))*1.e6
effrs(i,k) = MAX(re_qs_min, MIN(effrs(i,k), re_qs_max))*1.e6
end do
effrl(i,lmk) = re_qc_min*1.e6
effri(i,lmk) = re_qi_min*1.e6
effrs(i,lmk) = re_qs_min*1.e6
end do
do k=1,lm
k1 = k + kd
do i=1,im
effrl(i,k1) = re_cloud (i,k)
effri(i,k1) = re_ice (i,k)
effrr(i,k1) = 1000. ! rrain_def=1000.
effrs(i,k1) = re_snow(i,k)
enddo
enddo
effrr(:,:) = 1000. ! rrain_def=1000.
! Update global arrays
do k=1,lm
k1 = k + kd
Expand Down Expand Up @@ -972,8 +959,7 @@ subroutine GFS_rrtmg_pre_run (im, levs, lm, lmk, lmp, n_var_lndp, &
ntrac-1, ntcw-1,ntiw-1,ntrw-1, &
ntsw-1,ntgl-1, &
im, lm, lmp, uni_cld, lmfshal, lmfdeep2, &
cldcov(:,1:LM), effrl_inout, &
effri_inout, effrs_inout, &
cldcov(:,1:LM), effrl, effri, effrs, &
lwp_ex, iwp_ex, lwp_fc, iwp_fc, &
dzb, xlat_d, julian, yearlen, gridkm, &
clouds, cldsa, mtopa ,mbota, de_lgth, alpha) ! --- outputs
Expand Down Expand Up @@ -1006,8 +992,7 @@ subroutine GFS_rrtmg_pre_run (im, levs, lm, lmk, lmp, n_var_lndp, &
ntrac-1, ntcw-1,ntiw-1,ntrw-1, &
ntsw-1,ntgl-1, &
im, lm, lmp, uni_cld, lmfshal, lmfdeep2, &
cldcov(:,1:LM), effrl_inout, &
effri_inout, effrs_inout, &
cldcov(:,1:LM), effrl, effri, effrs, &
lwp_ex, iwp_ex, lwp_fc, iwp_fc, &
dzb, xlat_d, julian, yearlen, gridkm, &
clouds, cldsa, mtopa ,mbota, de_lgth, alpha) ! --- outputs
Expand All @@ -1018,8 +1003,7 @@ subroutine GFS_rrtmg_pre_run (im, levs, lm, lmk, lmp, n_var_lndp, &
ntrac-1, ntcw-1,ntiw-1,ntrw-1, &
ntsw-1,ntgl-1, &
im, lmk, lmp, uni_cld, lmfshal, lmfdeep2, &
cldcov(:,1:LMK), cnvw, effrl_inout, &
effri_inout, effrs_inout, &
cldcov(:,1:LMK), cnvw, effrl, effri, effrs,&
lwp_ex, iwp_ex, lwp_fc, iwp_fc, &
dzb, xlat_d, julian, yearlen, &
clouds, cldsa, mtopa ,mbota, de_lgth, alpha) ! --- outputs
Expand Down
24 changes: 7 additions & 17 deletions physics/module_mp_thompson.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1252,16 +1252,6 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, &
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

Expand Down Expand Up @@ -2188,15 +2178,15 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, &
ni(k) = MAX(R2, ni1d(k)*rho(k))
if (ni(k).le. R2) then
lami = cie(2)/5.E-6
ni(k) = MIN(999.D3, cig(1)*oig2*ri(k)/am_i*lami**bm_i)
ni(k) = MIN(499.D3, cig(1)*oig2*ri(k)/am_i*lami**bm_i)
endif
L_qi(k) = .true.
lami = (am_i*cig(2)*oig1*ni(k)/ri(k))**obmi
ilami = 1./lami
xDi = (bm_i + mu_i + 1.) * ilami
if (xDi.lt. 5.E-6) then
lami = cie(2)/5.E-6
ni(k) = MIN(999.D3, cig(1)*oig2*ri(k)/am_i*lami**bm_i)
ni(k) = MIN(499.D3, cig(1)*oig2*ri(k)/am_i*lami**bm_i)
elseif (xDi.gt. 300.E-6) then
lami = cie(2)/300.E-6
ni(k) = cig(1)*oig2*ri(k)/am_i*lami**bm_i
Expand Down Expand Up @@ -2901,7 +2891,7 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, &

!> - Freezing of aqueous aerosols based on Koop et al (2001, Nature)
xni = smo0(k)+ni(k) + (pni_rfz(k)+pni_wfz(k)+pni_inu(k))*dtsave
if (is_aerosol_aware .AND. homogIce .AND. (xni.le.999.E3) &
if (is_aerosol_aware .AND. homogIce .AND. (xni.le.499.E3) &
& .AND.(temp(k).lt.238).AND.(ssati(k).ge.0.4) ) then
xnc = iceKoop(temp(k),qv(k),qvs(k),nwfa(k), dtsave)
pni_iha(k) = xnc*odts
Expand Down Expand Up @@ -3237,7 +3227,7 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, &
xDi = (bm_i + mu_i + 1.) * ilami
if (xDi.lt. 5.E-6) then
lami = cie(2)/5.E-6
xni = MIN(999.D3, cig(1)*oig2*xri/am_i*lami**bm_i)
xni = MIN(499.D3, cig(1)*oig2*xri/am_i*lami**bm_i)
niten(k) = (xni-ni1d(k)*rho(k))*odts*orho
elseif (xDi.gt. 300.E-6) then
lami = cie(2)/300.E-6
Expand All @@ -3248,8 +3238,8 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, &
niten(k) = -ni1d(k)*odts
endif
xni=MAX(0.,(ni1d(k) + niten(k)*dtsave)*rho(k))
if (xni.gt.999.E3) &
niten(k) = (999.E3-ni1d(k)*rho(k))*odts*orho
if (xni.gt.499.E3) &
niten(k) = (499.E3-ni1d(k)*rho(k))*odts*orho

!> - Rain tendency
qrten(k) = qrten(k) + (prr_wau(k) + prr_rcw(k) &
Expand Down Expand Up @@ -4187,7 +4177,7 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, &
lami = cie(2)/300.E-6
endif
ni1d(k) = MIN(cig(1)*oig2*qi1d(k)/am_i*lami**bm_i, &
999.D3/rho(k))
499.D3/rho(k))
endif
qr1d(k) = qr1d(k) + qrten(k)*DT
nr1d(k) = MAX(R2/rho(k), nr1d(k) + nrten(k)*DT)
Expand Down
68 changes: 39 additions & 29 deletions physics/mp_thompson.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -91,7 +91,7 @@ subroutine mp_thompson_init(ncol, nlev, con_g, con_rd, con_eps, &
real(kind_phys) :: orho(1:ncol,1:nlev) ! m3 kg-1
real(kind_phys) :: nc_local(1:ncol,1:nlev) ! needed because nc is only allocated if is_aerosol_aware is true
!
real (kind=kind_phys) :: h_01, airmass, niIN3, niCCN3
real (kind=kind_phys) :: h_01, z1, niIN3, niCCN3
integer :: i, k

! Initialize the CCPP error handling variables
Expand Down Expand Up @@ -192,8 +192,8 @@ subroutine mp_thompson_init(ncol, nlev, con_g, con_rd, con_eps, &
endif
niCCN3 = -1.0*ALOG(naCCN1/naCCN0)/h_01
nwfa(i,1) = naCCN1+naCCN0*exp(-((hgt(i,2)-hgt(i,1))/1000.)*niCCN3)
airmass = 1./orho(i,1) * (hgt(i,2)-hgt(i,1))*area(i) ! kg
nwfa2d(i) = nwfa(i,1) * 0.000196 * (airmass*5.E-11)
z1 = hgt(i,2)-hgt(i,1)
nwfa2d(i) = nwfa(i,1) * 0.000196 * (50./z1)
do k = 2, nlev
nwfa(i,k) = naCCN1+naCCN0*exp(-((hgt(i,k)-hgt(i,1))/1000.)*niCCN3)
enddo
Expand All @@ -212,8 +212,8 @@ subroutine mp_thompson_init(ncol, nlev, con_g, con_rd, con_eps, &
!+---+-----------------------------------------------------------------+
if (mpirank==mpiroot) write(*,*) ' Apparently there are no initial CCN aerosol surface emission rates.'
do i = 1, ncol
airmass = 1./orho(i,1) * (hgt(i,2)-hgt(i,1))*area(i) ! kg
nwfa2d(i) = nwfa(i,1) * 0.000196 * (airmass*5.E-11)
z1 = hgt(i,2)-hgt(i,1)
nwfa2d(i) = nwfa(i,1) * 0.000196 * (50./z1)
enddo
else
if (mpirank==mpiroot) write(*,*) ' Apparently initial CCN aerosol surface emission rates are present.'
Expand Down Expand Up @@ -377,6 +377,7 @@ subroutine mp_thompson_run(ncol, nlev, con_g, con_rd, &

! Reduced time step if subcycling is used
real(kind_phys) :: dtstep
integer :: ndt
! Air density
real(kind_phys) :: rho(1:ncol,1:nlev) !< kg m-3
! Water vapor mixing ratio (instead of specific humidity)
Expand Down Expand Up @@ -456,11 +457,39 @@ subroutine mp_thompson_run(ncol, nlev, con_g, con_rd, &
errmsg = ''
errflg = 0

! Check initialization state
if (.not.is_initialized) then
write(errmsg, fmt='((a))') 'mp_thompson_run called before mp_thompson_init'
errflg = 1
return
if (first_time_step .and. istep==1 .and. blkno==1) then
! Check initialization state
if (.not.is_initialized) then
write(errmsg, fmt='((a))') 'mp_thompson_run called before mp_thompson_init'
errflg = 1
return
end if
! Check forr optional arguments of aerosol-aware microphysics
if (is_aerosol_aware .and. .not. (present(nc) .and. &
present(nwfa) .and. &
present(nifa) .and. &
present(nwfa2d) .and. &
present(nifa2d) )) then
write(errmsg,fmt='(*(a))') 'Logic error in mp_thompson_run:', &
' aerosol-aware microphysics require all of the', &
' following optional arguments:', &
' nc, nwfa, nifa, nwfa2d, nifa2d'
errflg = 1
return
end if
! Consistency cheecks - subcycling and inner loop at the same time are not supported
if (nsteps>1 .and. dt_inner < dtp) then
write(errmsg,'(*(a))') "Logic error: Subcycling and inner loop cannot be used at the same time"
errflg = 1
return
else if (mpirank==mpiroot .and. nsteps>1) then
write(*,'(a,i0,a,a,f6.2,a)') 'Thompson MP is using ', nsteps, ' substep(s) per time step with an ', &
'effective time step of ', dtp/real(nsteps, kind=kind_phys), ' seconds'
else if (mpirank==mpiroot .and. dt_inner < dtp) then
ndt = max(nint(dtp/dt_inner),1)
write(*,'(a,i0,a,a,f6.2,a)') 'Thompson MP is using ', ndt, ' inner loops per time step with an ', &
'effective time step of ', dtp/real(ndt, kind=kind_phys), ' seconds'
end if
end if

! Set reduced time step if subcycling is used
Expand All @@ -469,25 +498,6 @@ subroutine mp_thompson_run(ncol, nlev, con_g, con_rd, &
else
dtstep = dtp
end if
if (first_time_step .and. istep==1 .and. mpirank==mpiroot .and. blkno==1) then
write(*,'(a,i0,a,a,f8.2,a)') 'Thompson MP is using ', nsteps, ' substep(s) per time step', &
' with an effective time step of ', dtstep, ' seconds'
end if

if (first_time_step .and. istep==1) then
if (is_aerosol_aware .and. .not. (present(nc) .and. &
present(nwfa) .and. &
present(nifa) .and. &
present(nwfa2d) .and. &
present(nifa2d) )) then
write(errmsg,fmt='(*(a))') 'Logic error in mp_thompson_run:', &
' aerosol-aware microphysics require all of the', &
' following optional arguments:', &
' nc, nwfa, nifa, nwfa2d, nifa2d'
errflg = 1
return
end if
end if

!> - Convert specific humidity to water vapor mixing ratio.
!> - Also, hydrometeor variables are mass or number mixing ratio
Expand Down
Loading

0 comments on commit f6e185d

Please sign in to comment.