Skip to content

Commit

Permalink
Merge pull request #15 from climbfuji/move_number_concentration_bugfi…
Browse files Browse the repository at this point in the history 
…x_reorg_interstitials

Move number concentration bug fix, reorganization of  interstitial code
  • Loading branch information
@DomHeinzeller
DomHeinzeller authored Feb 18, 2020
2 parents ba242de + cebdfa4 commit e9621ef
Show file tree
Hide file tree
Showing 5 changed files with 116 additions and 159 deletions.
117 changes: 53 additions & 64 deletions physics/GFS_suite_interstitial.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -463,13 +463,13 @@ end subroutine GFS_suite_interstitial_3_finalize
subroutine GFS_suite_interstitial_3_run (im, levs, nn, cscnv, &
satmedmf, trans_trac, do_shoc, ltaerosol, ntrac, ntcw, &
ntiw, ntclamt, ntrw, ntsw, ntrnc, ntsnc, ntgl, ntgnc, &
xlat, gq0, imp_physics, imp_physics_mg, &
xlat, gt0, gq0, imp_physics, imp_physics_mg, &
imp_physics_zhao_carr, imp_physics_zhao_carr_pdf, &
imp_physics_gfdl, imp_physics_thompson, &
imp_physics_wsm6, imp_physics_fer_hires, prsi, &
prsl, prslk, rhcbot,rhcpbl, rhctop, rhcmax, islmsk, &
work1, work2, kpbl, kinver,clw, rhc, save_qc, save_qi, &
errmsg, errflg)
save_tcp, errmsg, errflg)

use machine, only: kind_phys

Expand All @@ -487,11 +487,13 @@ subroutine GFS_suite_interstitial_3_run (im, levs, nn, cscnv, &
real(kind=kind_phys), dimension(im, levs), intent(in) :: prsl, prslk
real(kind=kind_phys), dimension(im, levs+1), intent(in) :: prsi
real(kind=kind_phys), dimension(im), intent(in) :: xlat
real(kind=kind_phys), dimension(im, levs), intent(in) :: gt0
real(kind=kind_phys), dimension(im, levs, ntrac), intent(in) :: gq0

real(kind=kind_phys), dimension(im, levs), intent(inout) :: rhc, save_qc
! save_qi is not allocated for Zhao-Carr MP
real(kind=kind_phys), dimension(:, :), intent(inout) :: save_qi
real(kind=kind_phys), dimension(:, :), intent(inout) :: save_tcp ! ONLY ALLOCATE FOR THOMPSON! TODO
real(kind=kind_phys), dimension(im, levs, nn), intent(inout) :: clw

character(len=*), intent(out) :: errmsg
Expand All @@ -512,33 +514,6 @@ subroutine GFS_suite_interstitial_3_run (im, levs, nn, cscnv, &
errmsg = ''
errflg = 0

!GF* The following section (initializing convective variables) is already executed in GFS_typedefs%interstitial_phys_reset
! do k=1,levs
! do i=1,im
! clw(i,k,1) = 0.0
! clw(i,k,2) = -999.9
! enddo
! enddo
! if (Model%imfdeepcnv >= 0 .or. Model%imfshalcnv > 0 .or. &
! (Model%npdf3d == 3 .and. Model%num_p3d == 4) .or. &
! (Model%npdf3d == 0 .and. Model%ncnvcld3d == 1) ) then
! do k=1,levs
! do i=1,im
! cnvc(i,k) = 0.0
! cnvw(i,k) = 0.0
! enddo
! enddo
! endif
! if(imp_physics == 8) then
! if(Model%ltaerosol) then
! ice00 (:,:) = 0.0
! liq0 (:,:) = 0.0
! else
! ice00 (:,:) = 0.0
! endif
! endif
!*GF

if (cscnv .or. satmedmf .or. trans_trac ) then
tracers = 2
do n=2,ntrac
Expand Down Expand Up @@ -596,6 +571,8 @@ subroutine GFS_suite_interstitial_3_run (im, levs, nn, cscnv, &
enddo
enddo
endif
else
rhc(:,:) = 1.0
endif

if (imp_physics == imp_physics_zhao_carr .or. imp_physics == imp_physics_zhao_carr_pdf) then ! zhao-carr microphysics
Expand All @@ -615,8 +592,9 @@ subroutine GFS_suite_interstitial_3_run (im, levs, nn, cscnv, &
elseif (imp_physics == imp_physics_thompson) then
do k=1,levs
do i=1,im
clw(i,k,1) = gq0(i,k,ntiw) ! ice
clw(i,k,2) = gq0(i,k,ntcw) ! water
clw(i,k,1) = gq0(i,k,ntiw) ! ice
clw(i,k,2) = gq0(i,k,ntcw) ! water
save_tcp(i,k) = gt0(i,k)
enddo
enddo
if(ltaerosol) then
Expand All @@ -632,15 +610,7 @@ subroutine GFS_suite_interstitial_3_run (im, levs, nn, cscnv, &
clw(i,k,2) = gq0(i,k,ntcw) ! water
enddo
enddo
else ! if_ntcw
!GF* never executed unless imp_physics = imp_physics_zhao_carr or imp_physics_zhao_carr_pdf
! do i=1,im
! psautco_l(i) = Model%psautco(1)*work1(i) + Model%psautco(2)*work2(i)
! prautco_l(i) = Model%prautco(1)*work1(i) + Model%prautco(2)*work2(i)
! enddo
!*GF
rhc(:,:) = 1.0
endif ! end if_ntcw
endif

end subroutine GFS_suite_interstitial_3_run

Expand All @@ -662,9 +632,10 @@ end subroutine GFS_suite_interstitial_4_finalize
subroutine GFS_suite_interstitial_4_run (im, levs, ltaerosol, cplchm, tracers_total, ntrac, ntcw, ntiw, ntclamt, &
ntrw, ntsw, ntrnc, ntsnc, ntgl, ntgnc, ntlnc, ntinc, nn, imp_physics, imp_physics_gfdl, imp_physics_thompson, &
imp_physics_zhao_carr, imp_physics_zhao_carr_pdf, dtf, save_qc, save_qi, con_pi, &
gq0, clw, dqdti, imfdeepcnv, imfdeepcnv_gf, errmsg, errflg)
gq0, clw, prsl, save_tcp, con_rd, nwfa, spechum, dqdti, imfdeepcnv, imfdeepcnv_gf, errmsg, errflg)

use machine, only: kind_phys
use module_mp_thompson_make_number_concentrations, only: make_IceNumber, make_DropletNumber

implicit none

Expand All @@ -683,6 +654,11 @@ subroutine GFS_suite_interstitial_4_run (im, levs, ltaerosol, cplchm, tracers_to

real(kind=kind_phys), dimension(im,levs,ntrac), intent(inout) :: gq0
real(kind=kind_phys), dimension(im,levs,nn), intent(inout) :: clw
real(kind=kind_phys), dimension(im,levs), intent(in) :: prsl
real(kind=kind_phys), intent(in) :: con_rd
real(kind=kind_phys), dimension(:,:), intent(in) :: nwfa, save_tcp
real(kind=kind_phys), dimension(im,levs), intent(in) :: spechum

! dqdti may not be allocated
real(kind=kind_phys), dimension(:,:), intent(inout) :: dqdti

Expand All @@ -693,10 +669,12 @@ subroutine GFS_suite_interstitial_4_run (im, levs, ltaerosol, cplchm, tracers_to
! local variables
integer :: i,k,n,tracers

real(kind=kind_phys) :: liqm, icem

liqm = 4./3.*con_pi*1.e-12
icem = 4./3.*con_pi*3.2768*1.e-14*890.
real(kind=kind_phys), dimension(im,levs) :: rho_dryair
real(kind=kind_phys), dimension(im,levs) :: qv_mp !< kg kg-1 (dry mixing ratio)
real(kind=kind_phys), dimension(im,levs) :: qc_mp !< kg kg-1 (dry mixing ratio)
real(kind=kind_phys), dimension(im,levs) :: qi_mp !< kg kg-1 (dry mixing ratio)
real(kind=kind_phys), dimension(im,levs) :: nc_mp !< kg-1 (dry mixing ratio)
real(kind=kind_phys), dimension(im,levs) :: ni_mp !< kg-1 (dry mixing ratio)

! Initialize CCPP error handling variables
errmsg = ''
Expand Down Expand Up @@ -729,32 +707,42 @@ subroutine GFS_suite_interstitial_4_run (im, levs, ltaerosol, cplchm, tracers_to
imp_physics == imp_physics_zhao_carr_pdf .or. &
imp_physics == imp_physics_gfdl) then
gq0(1:im,:,ntcw) = clw(1:im,:,1) + clw(1:im,:,2)

elseif (ntiw > 0) then
do k=1,levs
do i=1,im
gq0(i,k,ntiw) = clw(i,k,1) ! ice
gq0(i,k,ntcw) = clw(i,k,2) ! water
enddo
enddo
! if (imp_physics == imp_physics_thompson) then
if (imp_physics == imp_physics_thompson .and. imfdeepcnv /= imfdeepcnv_gf) then
if (ltaerosol) then
do k=1,levs
do i=1,im
gq0(i,k,ntlnc) = gq0(i,k,ntlnc) &
+ max(0.0, (clw(i,k,2)-save_qc(i,k))) / liqm
gq0(i,k,ntinc) = gq0(i,k,ntinc) &
+ max(0.0, (clw(i,k,1)-save_qi(i,k))) / icem
enddo
enddo
else
do k=1,levs
do i=1,im
gq0(i,k,ntinc) = gq0(i,k,ntinc) &
+ max(0.0, (clw(i,k,1)-save_qi(i,k))) / icem
enddo
enddo
endif

if (imp_physics == imp_physics_thompson .and. (ntlnc>0 .or. ntinc>0)) then
do k=1,levs
do i=1,im
!> - Density of air in kg m-3
rho_dryair(i,k) = prsl(i,k)/(con_rd*save_tcp(i,k))
!> - Convert specific humidity to dry mixing ratio
qv_mp(i,k) = spechum(i,k)/(1.0_kind_phys-spechum(i,k))
if (ntlnc>0) then
!> - Convert moist mixing ratio to dry mixing ratio
qc_mp(i,k) = save_qc(i,k)/(1.0_kind_phys-spechum(i,k))
!> - Convert number concentration from moist to dry
nc_mp(i,k) = gq0(i,k,ntlnc)/(1.0_kind_phys-spechum(i,k))
nc_mp(i,k) = nc_mp(i,k) + max(0.0, make_DropletNumber(qc_mp(i,k) * rho_dryair(i,k), nwfa(i,k)) * (1.0/rho_dryair(i,k)))
!> - Convert number concentrations from dry to moist
gq0(i,k,ntlnc) = nc_mp(i,k)/(1.0_kind_phys+qv_mp(i,k))
endif
if (ntinc>0) then
!> - Convert moist mixing ratio to dry mixing ratio
qi_mp(i,k) = save_qi(i,k)/(1.0_kind_phys-spechum(i,k))
!> - Convert number concentration from moist to dry
ni_mp(i,k) = gq0(i,k,ntinc)/(1.0_kind_phys-spechum(i,k))
ni_mp(i,k) = ni_mp(i,k) + max(0.0, make_IceNumber(qi_mp(i,k) * rho_dryair(i,k), save_tcp(i,k)) * (1.0/rho_dryair(i,k)))
!> - Convert number concentrations from dry to moist
gq0(i,k,ntinc) = ni_mp(i,k)/(1.0_kind_phys+qv_mp(i,k))
endif
enddo
enddo
endif

else
Expand All @@ -764,6 +752,7 @@ subroutine GFS_suite_interstitial_4_run (im, levs, ltaerosol, cplchm, tracers_to
enddo
enddo
endif ! end if_ntiw

else
do k=1,levs
do i=1,im
Expand Down
63 changes: 63 additions & 0 deletions physics/GFS_suite_interstitial.meta
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1218,6 +1218,15 @@
kind = kind_phys
intent = in
optional = F
[gt0]
standard_name = air_temperature_updated_by_physics
long_name = temperature updated by physics
units = K
dimensions = (horizontal_dimension,vertical_dimension)
type = real
kind = kind_phys
intent = in
optional = F
[gq0]
standard_name = tracer_concentration_updated_by_physics
long_name = tracer concentration updated by physics
Expand Down Expand Up @@ -1432,6 +1441,15 @@
kind = kind_phys
intent = inout
optional = F
[save_tcp]
standard_name = air_temperature_save_from_cumulus_paramterization
long_name = air temperature after cumulus parameterization
units = K
dimensions = (horizontal_dimension,vertical_dimension)
type = real
kind = kind_phys
intent = inout
optional = F
[errmsg]
standard_name = ccpp_error_message
long_name = error message for error handling in CCPP
Expand Down Expand Up @@ -1692,6 +1710,51 @@
kind = kind_phys
intent = inout
optional = F
[prsl]
standard_name = air_pressure
long_name = mean layer pressure
units = Pa
dimensions = (horizontal_dimension,vertical_dimension)
type = real
kind = kind_phys
intent = in
optional = F
[save_tcp]
standard_name = air_temperature_save_from_cumulus_paramterization
long_name = air temperature after cumulus parameterization
units = K
dimensions = (horizontal_dimension,vertical_dimension)
type = real
kind = kind_phys
intent = in
optional = F
[con_rd]
standard_name = gas_constant_dry_air
long_name = ideal gas constant for dry air
units = J kg-1 K-1
dimensions = ()
type = real
kind = kind_phys
intent = in
optional = F
[nwfa]
standard_name = water_friendly_aerosol_number_concentration
long_name = number concentration of water-friendly aerosols
units = kg-1
dimensions = (horizontal_dimension,vertical_dimension)
type = real
kind = kind_phys
intent = in
optional = F
[spechum]
standard_name = water_vapor_specific_humidity
long_name = water vapor specific humidity
units = kg kg-1
dimensions = (horizontal_dimension,vertical_dimension)
type = real
kind = kind_phys
intent = inout
optional = F
[dqdti]
standard_name = instantaneous_water_vapor_specific_humidity_tendency_due_to_convection
long_name = instantaneous moisture tendency due to convection
Expand Down
26 changes: 0 additions & 26 deletions physics/cu_gf_driver.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,6 @@ module cu_gf_driver
use machine , only: kind_phys
use cu_gf_deep, only: cu_gf_deep_run,neg_check,autoconv,aeroevap,fct1d3
use cu_gf_sh , only: cu_gf_sh_run
use module_mp_thompson_make_number_concentrations, only: make_IceNumber, make_DropletNumber

implicit none

Expand Down Expand Up @@ -74,7 +73,6 @@ subroutine cu_gf_driver_run(ntracer,garea,im,ix,km,dt,cactiv, &
us,vs,t2di,w,qv2di_spechum,p2di,psuri, &
hbot,htop,kcnv,xland,hfx2,qfx2,cliw,clcw, &
pbl,ud_mf,dd_mf,dt_mf,cnvw_moist,cnvc,imfshalcnv, &
nwfa,con_rd,gq0,ntinc,ntlnc,imp_physics,imp_physics_thompson, &
errmsg,errflg)
!-------------------------------------------------------------
implicit none
Expand Down Expand Up @@ -126,12 +124,6 @@ subroutine cu_gf_driver_run(ntracer,garea,im,ix,km,dt,cactiv, &
real(kind=kind_phys), dimension( im ),intent(in) :: garea
real(kind=kind_phys), intent(in ) :: dt

! additional variables for number concentrations
real(kind=kind_phys), intent(in) :: nwfa(1:im,1:km)
real(kind=kind_phys), intent(in) :: con_rd
real(kind=kind_phys), dimension(im,km,ntracer), intent(inout) :: gq0
integer, intent(in) :: imp_physics,imp_physics_thompson,ntlnc,ntinc

integer, intent(in ) :: imfshalcnv
character(len=*), intent(out) :: errmsg
integer, intent(out) :: errflg
Expand Down Expand Up @@ -826,26 +818,8 @@ subroutine cu_gf_driver_run(ntracer,garea,im,ix,km,dt,cactiv, &
cliw(i,k) = max(0.,cliw(i,k) + tem)
endif

!
!> calculate cloud water and cloud ice number concentrations
!
rho_dryar(i,k) = p2di(i,k)/(con_rd*t(i,k)) ! Density of dry air in kg m-3
if (imp_physics == imp_physics_thompson) then
if ((tem*tem1)>1.e-5) then
gq0(i,k,ntinc) = max(0., gq0(i,k,ntinc) + &
make_IceNumber(tem*tem1*rho_dryar(i,k), t(i,k)) * &
(1/rho_dryar(i,k)))
end if
if ((tem*(1-tem1))>1.e-5) then
gq0(i,k,ntlnc) = max(0., gq0(i,k,ntlnc) + &
make_DropletNumber(tem*(1-tem1)*rho_dryar(i,k), nwfa(i,k)) &
* (1/rho_dryar(i,k)))
end if
end if

enddo


gdc(i,1,10)=forcing(i,1)
gdc(i,2,10)=forcing(i,2)
gdc(i,3,10)=forcing(i,3)
Expand Down
Loading

0 comments on commit e9621ef

Please sign in to comment.