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.

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

dtc/develop: update from EMC 2019/11/27 #360

Merged
2 changes: 1 addition & 1 deletion physics/GFS_surface_composites.F90
Original file line number Diff line number Diff line change
Expand Up @@ -123,7 +123,7 @@ subroutine GFS_surface_composites_pre_run (im, frac_grid, flag_cice, cplflx, lan
if (cice(i) < one) then
wet(i) = .true.
! tsfco(i) = tgice
tsfco(i) = max(tisfc(i), tgice)
if (.not. cplflx) tsfco(i) = max(tisfc(i), tgice)
! tsfco(i) = max((tsfc(i) - cice(i)*tisfc(i)) &
! / (one - cice(i)), tgice)
endif
Expand Down
44 changes: 16 additions & 28 deletions physics/drag_suite.F90
Original file line number Diff line number Diff line change
Expand Up @@ -332,11 +332,11 @@ subroutine drag_suite_run( &
& hpbl(im), &
& slmsk(im)
real(kind=kind_phys), dimension(im) :: govrth,xland
real(kind=kind_phys), dimension(im,km) :: dz2
!real(kind=kind_phys), dimension(im,km) :: dz2
real(kind=kind_phys) :: tauwavex0,tauwavey0, &
& XNBV,density,tvcon,hpbl2
integer :: kpbl2,kvar
real(kind=kind_phys), dimension(im,km+1) :: zq ! = PHII/g
!real(kind=kind_phys), dimension(im,km+1) :: zq ! = PHII/g
real(kind=kind_phys), dimension(im,km) :: zl ! = PHIL/g

!SPP
Expand Down Expand Up @@ -413,10 +413,10 @@ subroutine drag_suite_run( &
! local variables
!
integer :: i,j,k,lcap,lcapp1,nwd,idir, &
klcap,kp1,ikount,kk
klcap,kp1
!
real(kind=kind_phys) :: rcs,rclcs,csg,fdir,cleff,cleff_ss,cs, &
rcsks,wdir,ti,rdz,temp,tem2,dw2,shr2, &
real(kind=kind_phys) :: rcs,csg,fdir,cleff,cleff_ss,cs, &
rcsks,wdir,ti,rdz,tem2,dw2,shr2, &
bvf2,rdelks,wtkbj,tem,gfobnv,hd,fro, &
rim,temc,tem1,efact,temv,dtaux,dtauy, &
dtauxb,dtauyb,eng0,eng1
Expand All @@ -442,7 +442,6 @@ subroutine drag_suite_run( &
coefm(im),coefm_ss(im)
!
integer :: kbl(im),klowtop(im)
logical :: iope
integer,parameter :: mdir=8
!integer :: nwdir(mdir)
!data nwdir/6,7,5,8,2,3,1,4/
Expand Down Expand Up @@ -661,6 +660,17 @@ subroutine drag_suite_run( &
enddo
enddo
!
! calculate mid-layer height (zl), interface height (zq), and layer depth (dz2).
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Not a big deal, but this comment no longer reflects what the code is doing.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I actually left it on purpose (and just commented out the parts that are no longer calculated, zq and dz2) in case Mike wants to put it back later.

!
!zq=0.
do k = kts,km
do i = its,im
!zq(i,k+1) = PHII(i,k+1)*g_inv
!dz2(i,k) = (PHII(i,k+1)-PHII(i,k))*g_inv
zl(i,k) = PHIL(i,k)*g_inv
enddo
enddo
!
! determine reference level: maximum of 2*var and pbl heights
!
do i = its,im
Expand Down Expand Up @@ -895,7 +905,6 @@ subroutine drag_suite_run( &
density=1.2
utendwave=0.
vtendwave=0.
zq=0.
!
IF ( (gwd_opt_ss .EQ. 1).and.(ss_taper.GT.1.E-02) ) THEN
if (me==master) print *,"in Drag Suite: Running small-scale gravity wave drag"
Expand All @@ -914,14 +923,6 @@ subroutine drag_suite_run( &
thvx(i,k) = thx(i,k)*tvcon
enddo
enddo
! Calculate mid-layer height (zl), interface height (zq), and layer depth (dz2).
do k = kts,km
do i = its,im
zq(i,k+1) = PHII(i,k+1)*g_inv
dz2(i,k) = (PHII(i,k+1)-PHII(i,k))*g_inv
zl(i,k) = PHIL(i,k)*g_inv
enddo
enddo

do i=its,im
hpbl2 = hpbl(i)+10.
Expand Down Expand Up @@ -1027,19 +1028,6 @@ subroutine drag_suite_run( &

utendform=0.
vtendform=0.
zq=0.

IF ( (gwd_opt_ss .NE. 1).and.(ss_taper.GT.1.E-02) ) THEN
! Defining mid-layer height (zl), interface height (zq), and layer depth (dz2).
! This is already done above if the small-scale GWD is activated.
do k = kts,km
do i = its,im
zq(i,k+1) = PHII(i,k+1)*g_inv
dz2(i,k) = (PHII(i,k+1)-PHII(i,k))*g_inv
zl(i,k) = PHIL(i,k)*g_inv
enddo
enddo
ENDIF

DO i=its,im
IF ((xland(i)-1.5) .le. 0.) then
Expand Down
8 changes: 4 additions & 4 deletions physics/module_gfdl_cloud_microphys.F90
Original file line number Diff line number Diff line change
Expand Up @@ -4729,7 +4729,7 @@ subroutine cloud_diagnosis (is, ie, ks, ke, den, delp, lsm, qmw, qmi, qmr, qms,
real :: n0r = 8.0e6, n0s = 3.0e6, n0g = 4.0e6
real :: alphar = 0.8, alphas = 0.25, alphag = 0.5
real :: gammar = 17.837789, gammas = 8.2850630, gammag = 11.631769
real :: qmin = 1.0e-12, beta = 1.22
real :: qmin = 1.0e-12, beta = 1.22, qmin1 = 9.e-6

do k = ks, ke
do i = is, ie
Expand Down Expand Up @@ -4759,7 +4759,7 @@ subroutine cloud_diagnosis (is, ie, ks, ke, den, delp, lsm, qmw, qmi, qmr, qms,
! cloud ice (Heymsfield and Mcfarquhar, 1996)
! -----------------------------------------------------------------------

if (qmi (i, k) .gt. qmin) then
if (qmi (i, k) .gt. qmin1) then
qci (i, k) = dpg * qmi (i, k) * 1.0e3
rei_fac = log (1.0e3 * qmi (i, k) * den (i, k))
if (t (i, k) - tice .lt. - 50) then
Expand All @@ -4785,7 +4785,7 @@ subroutine cloud_diagnosis (is, ie, ks, ke, den, delp, lsm, qmw, qmi, qmr, qms,
! cloud ice (Wyser, 1998)
! -----------------------------------------------------------------------

if (qmi (i, k) .gt. qmin) then
if (qmi (i, k) .gt. qmin1) then
qci (i, k) = dpg * qmi (i, k) * 1.0e3
bw = - 2. + 1.e-3 * log10 (den (i, k) * qmi (i, k) / rho_0) * max (0.0, tice - t (i, k)) ** 1.5
rei (i, k) = 377.4 + bw * (203.3 + bw * (37.91 + 2.3696 * bw))
Expand Down Expand Up @@ -4815,7 +4815,7 @@ subroutine cloud_diagnosis (is, ie, ks, ke, den, delp, lsm, qmw, qmi, qmr, qms,
! snow (Lin et al., 1983)
! -----------------------------------------------------------------------

if (qms (i, k) .gt. qmin) then
if (qms (i, k) .gt. qmin1) then
qcs (i, k) = dpg * qms (i, k) * 1.0e3
lambdas = exp (0.25 * log (pi * rhos * n0s / qms (i, k) / den (i, k)))
res (i, k) = 0.5 * exp (log (gammas / 6) / alphas) / lambdas * 1.0e6
Expand Down
16 changes: 8 additions & 8 deletions physics/module_mp_thompson.F90
Original file line number Diff line number Diff line change
Expand Up @@ -1302,7 +1302,7 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, &
kmax_qc = k
qc_max = qc1d(k)
elseif (qc1d(k) .lt. 0.0) then
write(*,*) 'WARNING, negative qc ', qc1d(k), &
write(*,'(a,e16.7,a,3i8)') 'WARNING, negative qc ', qc1d(k), &
' at i,j,k=', i,j,k
endif
if (qr1d(k) .gt. qr_max) then
Expand All @@ -1311,7 +1311,7 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, &
kmax_qr = k
qr_max = qr1d(k)
elseif (qr1d(k) .lt. 0.0) then
write(*,*) 'WARNING, negative qr ', qr1d(k), &
write(*,'(a,e16.7,a,3i8)') 'WARNING, negative qr ', qr1d(k), &
' at i,j,k=', i,j,k
endif
if (nr1d(k) .gt. nr_max) then
Expand All @@ -1320,7 +1320,7 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, &
kmax_nr = k
nr_max = nr1d(k)
elseif (nr1d(k) .lt. 0.0) then
write(*,*) 'WARNING, negative nr ', nr1d(k), &
write(*,'(a,e16.7,a,3i8)') 'WARNING, negative nr ', nr1d(k), &
' at i,j,k=', i,j,k
endif
if (qs1d(k) .gt. qs_max) then
Expand All @@ -1329,7 +1329,7 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, &
kmax_qs = k
qs_max = qs1d(k)
elseif (qs1d(k) .lt. 0.0) then
write(*,*) 'WARNING, negative qs ', qs1d(k), &
write(*,'(a,e16.7,a,3i8)') 'WARNING, negative qs ', qs1d(k), &
' at i,j,k=', i,j,k
endif
if (qi1d(k) .gt. qi_max) then
Expand All @@ -1338,7 +1338,7 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, &
kmax_qi = k
qi_max = qi1d(k)
elseif (qi1d(k) .lt. 0.0) then
write(*,*) 'WARNING, negative qi ', qi1d(k), &
write(*,'(a,e16.7,a,3i8)') 'WARNING, negative qi ', qi1d(k), &
' at i,j,k=', i,j,k
endif
if (qg1d(k) .gt. qg_max) then
Expand All @@ -1347,7 +1347,7 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, &
kmax_qg = k
qg_max = qg1d(k)
elseif (qg1d(k) .lt. 0.0) then
write(*,*) 'WARNING, negative qg ', qg1d(k), &
write(*,'(a,e16.7,a,3i8)') 'WARNING, negative qg ', qg1d(k), &
' at i,j,k=', i,j,k
endif
if (ni1d(k) .gt. ni_max) then
Expand All @@ -1356,11 +1356,11 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, &
kmax_ni = k
ni_max = ni1d(k)
elseif (ni1d(k) .lt. 0.0) then
write(*,*) 'WARNING, negative ni ', ni1d(k), &
write(*,'(a,e16.7,a,3i8)') 'WARNING, negative ni ', ni1d(k), &
' at i,j,k=', i,j,k
endif
if (qv1d(k) .lt. 0.0) then
write(*,*) 'WARNING, negative qv ', qv1d(k), &
write(*,'(a,e16.7,a,3i8)') 'WARNING, negative qv ', qv1d(k), &
' at i,j,k=', i,j,k
if (k.lt.kte-2 .and. k.gt.kts+1) then
write(*,*) ' below and above are: ', qv(i,k-1,j), qv(i,k+1,j)
Expand Down
17 changes: 16 additions & 1 deletion physics/module_mp_thompson_make_number_concentrations.F90
Original file line number Diff line number Diff line change
Expand Up @@ -79,6 +79,11 @@ elemental real function make_IceNumber (Q_ice, temp)
161.503, 168.262, 175.248, 182.473, 189.952, 197.699, &
205.728, 214.055, 222.694, 231.661, 240.971, 250.639 /)

if (Q_ice == 0) then
make_IceNumber = 0
return
end if

!+---+-----------------------------------------------------------------+
!..From the model 3D temperature field, subtract 179K for which
!.. index value of retab as a start. Value of corr is for
Expand Down Expand Up @@ -133,6 +138,11 @@ elemental real function make_DropletNumber (Q_cloud, qnwfa)
real:: q_nwfa, x1, xDc
integer:: nu_c

if (Q_cloud == 0) then
make_DropletNumber = 0
return
end if

!+---+

q_nwfa = MAX(99.E6, MIN(qnwfa,5.E10))
Expand Down Expand Up @@ -160,6 +170,11 @@ elemental real function make_RainNumber (Q_rain, temp)
!real, parameter:: PI = 3.1415926536
real, parameter:: am_r = PI*1000./6.

if (Q_rain == 0) then
make_RainNumber = 0
return
end if

!+---+-----------------------------------------------------------------+
!.. Not thrilled with it, but set Y-intercept parameter to Marshal-Palmer value
!.. that basically assumes melting snow becomes typical rain. However, for
Expand All @@ -172,7 +187,7 @@ elemental real function make_RainNumber (Q_rain, temp)
N0 = 8.E6

if (temp .le. 271.15) then
N0 = 8.E8
N0 = 8.E8
elseif (temp .gt. 271.15 .and. temp.lt.273.15) then
N0 = 8. * 10**(279.15-temp)
endif
Expand Down
10 changes: 6 additions & 4 deletions physics/mp_thompson_post.F90
Original file line number Diff line number Diff line change
Expand Up @@ -67,7 +67,7 @@ end subroutine mp_thompson_post_init
!!
#endif
subroutine mp_thompson_post_run(ncol, nlev, tgrs_save, tgrs, prslk, dtp, &
mpicomm, mpirank, mpiroot, errmsg, errflg)
kdt, mpicomm, mpirank, mpiroot, errmsg, errflg)

implicit none

Expand All @@ -78,6 +78,7 @@ subroutine mp_thompson_post_run(ncol, nlev, tgrs_save, tgrs, prslk, dtp, &
real(kind_phys), dimension(1:ncol,1:nlev), intent(inout) :: tgrs
real(kind_phys), dimension(1:ncol,1:nlev), intent(in) :: prslk
real(kind_phys), intent(in) :: dtp
integer, intent(in) :: kdt
! MPI information
integer, intent(in ) :: mpicomm
integer, intent(in ) :: mpirank
Expand Down Expand Up @@ -115,8 +116,8 @@ subroutine mp_thompson_post_run(ncol, nlev, tgrs_save, tgrs, prslk, dtp, &

if (tgrs_save(i,k) + mp_tend(i,k)*prslk(i,k) .ne. tgrs(i,k)) then
#ifdef DEBUG
write(0,*) "mp_thompson_post_run mp_tend limiter: i, k, t_old, t_new, t_lim:", &
& i, k, tgrs_save(i,k), tgrs(i,k), tgrs_save(i,k) + mp_tend(i,k)*prslk(i,k)
write(0,'(a,3i6,3e16.7)') "mp_thompson_post_run mp_tend limiter: kdt, i, k, t_old, t_new, t_lim:", &
& kdt, i, k, tgrs_save(i,k), tgrs(i,k), tgrs_save(i,k) + mp_tend(i,k)*prslk(i,k)
#endif
events = events + 1
end if
Expand All @@ -125,7 +126,8 @@ subroutine mp_thompson_post_run(ncol, nlev, tgrs_save, tgrs, prslk, dtp, &
end do

if (events > 0) then
write(0,'(a,i0,a,i0,a)') "mp_thompson_post_run: mp_tend_lim applied ", events, "/", nlev*ncol, " times"
write(0,'(a,i0,a,i0,a,i0)') "mp_thompson_post_run: mp_tend_lim applied ", events, "/", nlev*ncol, &
& " times at timestep ", kdt
end if

end subroutine mp_thompson_post_run
Expand Down
8 changes: 8 additions & 0 deletions physics/mp_thompson_post.meta
Original file line number Diff line number Diff line change
Expand Up @@ -92,6 +92,14 @@
kind = kind_phys
intent = in
optional = F
[kdt]
standard_name = index_of_time_step
long_name = current forecast iteration
units = index
dimensions = ()
type = integer
intent = in
optional = F
[mpicomm]
standard_name = mpi_comm
long_name = MPI communicator
Expand Down
Loading