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

Update to the latest GFDL dev/emc #70

Merged
merged 16 commits into from
Mar 30, 2022
Merged

Update to the latest GFDL dev/emc #70

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
16 commits
Select commit Hold shift + click to select a range
7b763a2
Add get_nth_domain_info subroutine
DusanJovic-NOAA Sep 30, 2021
80ebd27
Deallocate local arrays in fv_dyn_bundle_setup
DusanJovic-NOAA Nov 3, 2021
fa86482
adding back a line that was mistakenly deleted in a previous commit (…
laurenchilutti Nov 30, 2021
fba2d78
Merge remote-tracking branch 'origin/dev/emc' into multiple_output_grids
DusanJovic-NOAA Jan 4, 2022
9af3191
Do not print debug messages to stderr
DusanJovic-NOAA Jan 5, 2022
da619f3
Merge remote-tracking branch 'gfdl/dev/emc' into fms_mixedmode
MinsukJi-NOAA Jan 31, 2022
5193c6b
fix for 4diau with iau_filter_increments=T (#167)
jswhit Feb 4, 2022
7ec6a47
Merge remote-tracking branch 'origin/dev/emc' into multiple_output_grids
DusanJovic-NOAA Feb 10, 2022
9d1fabc
Use mpp_error instead of write statements in model/fv_regional_bc.F90
DusanJovic-NOAA Feb 11, 2022
348cee7
Merge pull request #176 from DusanJovic-NOAA/multiple_output_grids
laurenchilutti Feb 15, 2022
7ce7aa9
Attempt at integrating fixes on top of dev/emc branch. (#173)
MatthewPyle-NOAA Feb 25, 2022
45cbd3d
Merge remote-tracking branch 'emc/fms_mixedmode' into fms_mixedmode
MinsukJi-NOAA Mar 2, 2022
9f08a54
Merge remote-tracking branch 'gfdl/dev/emc' into fms_mixedmode
MinsukJi-NOAA Mar 2, 2022
43f7ed3
Support for cloud microphysics hail species (#171)
MicroTed Mar 4, 2022
b5ab47d
Merge remote-tracking branch 'emc/fms_mixedmode' into fms_mixedmode
MinsukJi-NOAA Mar 30, 2022
c184369
Merge remote-tracking branch 'gfdl/dev/emc' into fms_mixedmode
MinsukJi-NOAA Mar 30, 2022
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
2 changes: 1 addition & 1 deletion model/fv_arrays.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -56,7 +56,7 @@ module fv_arrays_mod
integer :: id_ws, id_te, id_amdt, id_mdt, id_divg, id_aam
logical :: initialized = .false.
real sphum, liq_wat, ice_wat ! GFDL physics
real rainwat, snowwat, graupel
real rainwat, snowwat, graupel, hailwat

real :: efx(max_step), efx_sum, efx_nest(max_step), efx_sum_nest, mtq(max_step), mtq_sum
integer :: steps
Expand Down
66 changes: 58 additions & 8 deletions model/fv_dynamics.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -95,6 +95,10 @@ module fv_dynamics_mod
! <td>neg_adj2</td>
! </tr>
! <tr>
! <td>fv_sg_mod</td>
! <td>neg_adj4</td>
! </tr>
! <tr>
! <td>fv_timing_mod</td>
! <td>timing_on, timing_off</td>
! </tr>
Expand All @@ -115,6 +119,10 @@ module fv_dynamics_mod
! <td>neg_adj3</td>
! </tr>
! <tr>
! <td>fv_sg_mod</td>
! <td>neg_adj4</td>
! </tr>
! <tr>
! <td>tracer_manager_mod</td>
! <td>get_tracer_index</td>
! </tr>
Expand All @@ -140,7 +148,7 @@ module fv_dynamics_mod
use mpp_mod, only: mpp_pe
use field_manager_mod, only: MODEL_ATMOS
use tracer_manager_mod, only: get_tracer_index
use fv_sg_mod, only: neg_adj3, neg_adj2
use fv_sg_mod, only: neg_adj3, neg_adj2, neg_adj4
use fv_nesting_mod, only: setup_nested_grid_BCs
use fv_regional_mod, only: regional_boundary_update, set_regional_BCs
use fv_regional_mod, only: dump_field, H_STAGGER, U_STAGGER, V_STAGGER
Expand Down Expand Up @@ -283,7 +291,7 @@ subroutine fv_dynamics(npx, npy, npz, nq_tot, ng, bdt, consv_te, fill,
integer:: kord_tracer(ncnst)
integer :: i,j,k, n, iq, n_map, nq, nr, nwat, k_split
integer :: sphum, liq_wat = -999, ice_wat = -999 ! GFDL physics
integer :: rainwat = -999, snowwat = -999, graupel = -999, cld_amt = -999
integer :: rainwat = -999, snowwat = -999, graupel = -999, hailwat = -999, cld_amt = -999
integer :: theta_d = -999
logical used, do_omega
integer, parameter :: max_packs=13
Expand Down Expand Up @@ -399,6 +407,7 @@ subroutine fv_dynamics(npx, npy, npz, nq_tot, ng, bdt, consv_te, fill,
rainwat = get_tracer_index (MODEL_ATMOS, 'rainwat')
snowwat = get_tracer_index (MODEL_ATMOS, 'snowwat')
graupel = get_tracer_index (MODEL_ATMOS, 'graupel')
hailwat = get_tracer_index (MODEL_ATMOS, 'hailwat')
cld_amt = get_tracer_index (MODEL_ATMOS, 'cld_amt')
endif

Expand All @@ -424,12 +433,12 @@ subroutine fv_dynamics(npx, npy, npz, nq_tot, ng, bdt, consv_te, fill,
#else
!$OMP parallel do default(none) shared(is,ie,js,je,isd,ied,jsd,jed,npz,dp1,zvir,nwat,q,q_con,sphum,liq_wat, &
#endif
!$OMP rainwat,ice_wat,snowwat,graupel) private(cvm,i,j,k)
!$OMP rainwat,ice_wat,snowwat,graupel,hailwat) private(cvm,i,j,k)
do k=1,npz
do j=js,je
#ifdef USE_COND
call moist_cp(is,ie,isd,ied,jsd,jed, npz, j, k, nwat, sphum, liq_wat, rainwat, &
ice_wat, snowwat, graupel, q, q_con(is:ie,j,k), cvm)
ice_wat, snowwat, graupel, hailwat, q, q_con(is:ie,j,k), cvm)
#endif
do i=is,ie
dp1(i,j,k) = zvir*q(i,j,k,sphum)
Expand All @@ -442,7 +451,7 @@ subroutine fv_dynamics(npx, npy, npz, nq_tot, ng, bdt, consv_te, fill,
#else
!$OMP parallel do default(none) shared(is,ie,js,je,isd,ied,jsd,jed,npz,dp1,zvir,q,q_con,sphum,liq_wat, &
#endif
!$OMP rainwat,ice_wat,snowwat,graupel,pkz,flagstruct, &
!$OMP rainwat,ice_wat,snowwat,graupel,hailwat,pkz,flagstruct, &
#ifdef MULTI_GASES
!$OMP kapad, &
#endif
Expand All @@ -457,7 +466,7 @@ subroutine fv_dynamics(npx, npy, npz, nq_tot, ng, bdt, consv_te, fill,
do j=js,je
#ifdef MOIST_CAPPA
call moist_cv(is,ie,isd,ied,jsd,jed, npz, j, k, nwat, sphum, liq_wat, rainwat, &
ice_wat, snowwat, graupel, q, q_con(is:ie,j,k), cvm)
ice_wat, snowwat, graupel, hailwat, q, q_con(is:ie,j,k), cvm)
#endif
do i=is,ie
#ifdef MULTI_GASES
Expand Down Expand Up @@ -522,7 +531,7 @@ subroutine fv_dynamics(npx, npy, npz, nq_tot, ng, bdt, consv_te, fill,
gridstruct%rsin2, gridstruct%cosa_s, &
zvir, cp_air, rdgas, hlv, te_2d, ua, va, teq, &
flagstruct%moist_phys, nwat, sphum, liq_wat, rainwat, &
ice_wat, snowwat, graupel, hydrostatic, idiag%id_te)
ice_wat, snowwat, graupel, hailwat, hydrostatic, idiag%id_te)
if( idiag%id_te>0 ) then
used = send_data(idiag%id_te, teq, fv_time)
! te_den=1.E-9*g_sum(teq, is, ie, js, je, ng, area, 0)/(grav*4.*pi*radius**2)
Expand Down Expand Up @@ -728,6 +737,9 @@ subroutine fv_dynamics(npx, npy, npz, nq_tot, ng, bdt, consv_te, fill,
call fill2D(is, ie, js, je, ng, npz, q(isd,jsd,1,snowwat), delp, gridstruct%area, domain, gridstruct%bounded_domain, npx, npy)
if ( graupel > 0 ) &
call fill2D(is, ie, js, je, ng, npz, q(isd,jsd,1,graupel), delp, gridstruct%area, domain, gridstruct%bounded_domain, npx, npy)
if ( hailwat > 0 ) &
call fill2D(is, ie, js, je, ng, npz, q(isd,jsd,1,hailwat), delp, gridstruct%area, domain, gridstruct%bounded_domain, npx, npy)

call timing_off('Fill2D')
endif
#endif
Expand Down Expand Up @@ -786,6 +798,7 @@ subroutine fv_dynamics(npx, npy, npz, nq_tot, ng, bdt, consv_te, fill,
if (ice_wat > 0) call prt_mxm('ice_wat_dyn', q(isd,jsd,1,ice_wat), is, ie, js, je, ng, npz, 1.,gridstruct%area_64, domain)
if (snowwat > 0) call prt_mxm('snowwat_dyn', q(isd,jsd,1,snowwat), is, ie, js, je, ng, npz, 1.,gridstruct%area_64, domain)
if (graupel > 0) call prt_mxm('graupel_dyn', q(isd,jsd,1,graupel), is, ie, js, je, ng, npz, 1.,gridstruct%area_64, domain)
if (hailwat > 0) call prt_mxm('hailwat_dyn', q(isd,jsd,1,hailwat), is, ie, js, je, ng, npz, 1.,gridstruct%area_64, domain)
endif
#ifdef AVEC_TIMERS
call avec_timer_stop(6)
Expand Down Expand Up @@ -853,7 +866,44 @@ subroutine fv_dynamics(npx, npy, npz, nq_tot, ng, bdt, consv_te, fill,
used = send_data(idiag%id_mdt, dtdt_m, fv_time)
endif

if( nwat==6 ) then
if( nwat==7 ) then
if (cld_amt > 0) then
call neg_adj4(is, ie, js, je, ng, npz, &
flagstruct%hydrostatic, &
peln, delz, &
pt, delp, q(isd,jsd,1,sphum), &
q(isd,jsd,1,liq_wat), &
q(isd,jsd,1,rainwat), &
q(isd,jsd,1,ice_wat), &
q(isd,jsd,1,snowwat), &
q(isd,jsd,1,graupel), &
q(isd,jsd,1,hailwat), &
q(isd,jsd,1,cld_amt), flagstruct%check_negative)
else
call neg_adj4(is, ie, js, je, ng, npz, &
flagstruct%hydrostatic, &
peln, delz, &
pt, delp, q(isd,jsd,1,sphum), &
q(isd,jsd,1,liq_wat), &
q(isd,jsd,1,rainwat), &
q(isd,jsd,1,ice_wat), &
q(isd,jsd,1,snowwat), &
q(isd,jsd,1,graupel), &
q(isd,jsd,1,hailwat), check_negative=flagstruct%check_negative)
endif
if ( flagstruct%fv_debug ) then
call prt_mxm('T_dyn_a3', pt, is, ie, js, je, ng, npz, 1., gridstruct%area_64, domain)
call prt_mxm('SPHUM_dyn', q(isd,jsd,1,sphum ), is, ie, js, je, ng, npz, 1.,gridstruct%area_64, domain)
call prt_mxm('liq_wat_dyn', q(isd,jsd,1,liq_wat), is, ie, js, je, ng, npz, 1.,gridstruct%area_64, domain)
call prt_mxm('rainwat_dyn', q(isd,jsd,1,rainwat), is, ie, js, je, ng, npz, 1.,gridstruct%area_64, domain)
call prt_mxm('ice_wat_dyn', q(isd,jsd,1,ice_wat), is, ie, js, je, ng, npz, 1.,gridstruct%area_64, domain)
call prt_mxm('snowwat_dyn', q(isd,jsd,1,snowwat), is, ie, js, je, ng, npz, 1.,gridstruct%area_64, domain)
call prt_mxm('graupel_dyn', q(isd,jsd,1,graupel), is, ie, js, je, ng, npz, 1.,gridstruct%area_64, domain)
IF ( hailwat > 0 ) call prt_mxm('hailwat_dyn', q(isd,jsd,1,hailwat), is, ie, js, je, ng, npz, 1.,gridstruct%area_64, domain)
endif
endif

if( nwat == 6 ) then
if (cld_amt > 0) then
call neg_adj3(is, ie, js, je, ng, npz, &
flagstruct%hydrostatic, &
Expand Down
68 changes: 52 additions & 16 deletions model/fv_mapz.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -226,7 +226,7 @@ subroutine Lagrangian_to_Eulerian(last_step, consv, ps, pe, delp, pkz, pk, &
real rcp, rg, rrg, bkh, dtmp, k1k
integer:: i,j,k
integer:: kdelz
integer:: nt, liq_wat, ice_wat, rainwat, snowwat, cld_amt, graupel, ccn_cm3, iq, n, kp, k_next
integer:: nt, liq_wat, ice_wat, rainwat, snowwat, cld_amt, graupel, hailwat, ccn_cm3, iq, n, kmp, kp, k_next
integer :: ierr

ccpp_associate: associate( fast_mp_consv => CCPP_interstitial%fast_mp_consv, &
Expand All @@ -242,6 +242,7 @@ subroutine Lagrangian_to_Eulerian(last_step, consv, ps, pe, delp, pkz, pk, &
rainwat = get_tracer_index (MODEL_ATMOS, 'rainwat')
snowwat = get_tracer_index (MODEL_ATMOS, 'snowwat')
graupel = get_tracer_index (MODEL_ATMOS, 'graupel')
hailwat = get_tracer_index (MODEL_ATMOS, 'hailwat')
cld_amt = get_tracer_index (MODEL_ATMOS, 'cld_amt')
ccn_cm3 = get_tracer_index (MODEL_ATMOS, 'ccn_cm3')

Expand All @@ -251,7 +252,7 @@ subroutine Lagrangian_to_Eulerian(last_step, consv, ps, pe, delp, pkz, pk, &

!$OMP parallel do default(none) shared(is,ie,js,je,km,pe,ptop,kord_tm,hydrostatic, &
!$OMP pt,pk,rg,peln,q,nwat,liq_wat,rainwat,ice_wat,snowwat, &
!$OMP graupel,q_con,sphum,cappa,r_vir,rcp,k1k,delp, &
!$OMP graupel,hailwat,q_con,sphum,cappa,r_vir,rcp,k1k,delp, &
!$OMP delz,akap,pkz,te,u,v,ps, gridstruct, last_step, &
!$OMP ak,bk,nq,isd,ied,jsd,jed,kord_tr,fill, adiabatic, &
#ifdef MULTI_GASES
Expand Down Expand Up @@ -294,7 +295,7 @@ subroutine Lagrangian_to_Eulerian(last_step, consv, ps, pe, delp, pkz, pk, &
do k=1,km
#ifdef MOIST_CAPPA
call moist_cv(is,ie,isd,ied,jsd,jed, km, j, k, nwat, sphum, liq_wat, rainwat, &
ice_wat, snowwat, graupel, q, gz, cvm)
ice_wat, snowwat, graupel, hailwat, q, gz, cvm)
do i=is,ie
q_con(i,j,k) = gz(i)
#ifdef MULTI_GASES
Expand Down Expand Up @@ -545,7 +546,7 @@ subroutine Lagrangian_to_Eulerian(last_step, consv, ps, pe, delp, pkz, pk, &
do k=1,km
#ifdef MOIST_CAPPA
call moist_cv(is,ie,isd,ied,jsd,jed, km, j, k, nwat, sphum, liq_wat, rainwat, &
ice_wat, snowwat, graupel, q, gz, cvm)
ice_wat, snowwat, graupel, hailwat, q, gz, cvm)
do i=is,ie
q_con(i,j,k) = gz(i)
#ifdef MULTI_GASES
Expand Down Expand Up @@ -667,7 +668,7 @@ subroutine Lagrangian_to_Eulerian(last_step, consv, ps, pe, delp, pkz, pk, &
!$OMP parallel default(none) shared(is,ie,js,je,km,ptop,u,v,pe,ua,va,isd,ied,jsd,jed,kord_mt, &
!$OMP te_2d,te,delp,hydrostatic,hs,rg,pt,peln, adiabatic, &
!$OMP cp,delz,nwat,rainwat,liq_wat,ice_wat,snowwat, &
!$OMP graupel,q_con,r_vir,sphum,w,pk,pkz,last_step,consv, &
!$OMP graupel,hailwat,q_con,r_vir,sphum,w,pk,pkz,last_step,consv, &
!$OMP do_adiabatic_init,zsum1,zsum0,te0_2d,domain, &
!$OMP ng,gridstruct,E_Flux,pdt,dtmp,reproduce_sum,q, &
!$OMP mdt,cld_amt,cappa,dtdt,out_dt,rrg,akap,do_sat_adj, &
Expand All @@ -682,7 +683,7 @@ subroutine Lagrangian_to_Eulerian(last_step, consv, ps, pe, delp, pkz, pk, &
!$OMP parallel default(none) shared(is,ie,js,je,km,kmp,ptop,u,v,pe,ua,va,isd,ied,jsd,jed,kord_mt, &
!$OMP te_2d,te,delp,hydrostatic,hs,rg,pt,peln, adiabatic, &
!$OMP cp,delz,nwat,rainwat,liq_wat,ice_wat,snowwat, &
!$OMP graupel,q_con,r_vir,sphum,w,pk,pkz,last_step,consv, &
!$OMP graupel,hailwat,q_con,r_vir,sphum,w,pk,pkz,last_step,consv, &
!$OMP do_adiabatic_init,zsum1,zsum0,te0_2d,domain, &
!$OMP ng,gridstruct,E_Flux,pdt,dtmp,reproduce_sum,q, &
!$OMP mdt,cld_amt,cappa,dtdt,out_dt,rrg,akap,do_sat_adj, &
Expand Down Expand Up @@ -744,7 +745,7 @@ subroutine Lagrangian_to_Eulerian(last_step, consv, ps, pe, delp, pkz, pk, &
do k=1,km
#ifdef MOIST_CAPPA
call moist_cv(is,ie,isd,ied,jsd,jed, km, j, k, nwat, sphum, liq_wat, rainwat, &
ice_wat, snowwat, graupel, q, gz, cvm)
ice_wat, snowwat, graupel, hailwat, q, gz, cvm)
do i=is,ie
! KE using 3D winds:
q_con(i,j,k) = gz(i)
Expand Down Expand Up @@ -874,11 +875,20 @@ subroutine Lagrangian_to_Eulerian(last_step, consv, ps, pe, delp, pkz, pk, &
pt(i,j,k) = (pt(i,j,k)+dtmp*pkz(i,j,k))/ virq(q(i,j,k,1:num_gas))
#else
pt(i,j,k) = (pt(i,j,k)+dtmp*pkz(i,j,k))/((1.+r_vir*q(i,j,k,sphum))*(1.-gz(i)))
#endif
enddo
elseif ( nwat==7 ) then
do i=is,ie
gz(i) = q(i,j,k,liq_wat)+q(i,j,k,rainwat)+q(i,j,k,ice_wat)+q(i,j,k,snowwat)+q(i,j,k,graupel)+q(i,j,k,hailwat)
#ifdef MULTI_GASES
pt(i,j,k) = (pt(i,j,k)+dtmp*pkz(i,j,k))/ virq(q(i,j,k,1:num_gas))
#else
pt(i,j,k) = (pt(i,j,k)+dtmp*pkz(i,j,k))/((1.+r_vir*q(i,j,k,sphum))*(1.-gz(i)))
#endif
enddo
else
call moist_cv(is,ie,isd,ied,jsd,jed, km, j, k, nwat, sphum, liq_wat, rainwat, &
ice_wat, snowwat, graupel, q, gz, cvm)
ice_wat, snowwat, graupel, hailwat, q, gz, cvm)
do i=is,ie
#ifdef MULTI_GASES
pt(i,j,k) = (pt(i,j,k)+dtmp*pkz(i,j,k)) / virq(q(i,j,k,1:num_gas))
Expand Down Expand Up @@ -925,13 +935,13 @@ subroutine compute_total_energy(is, ie, js, je, isd, ied, jsd, jed, km, &
u, v, w, delz, pt, delp, q, qc, pe, peln, hs, &
rsin2_l, cosa_s_l, &
r_vir, cp, rg, hlv, te_2d, ua, va, teq, &
moist_phys, nwat, sphum, liq_wat, rainwat, ice_wat, snowwat, graupel, hydrostatic, id_te)
moist_phys, nwat, sphum, liq_wat, rainwat, ice_wat, snowwat, graupel, hailwat, hydrostatic, id_te)
!------------------------------------------------------
! Compute vertically integrated total energy per column
!------------------------------------------------------
! !INPUT PARAMETERS:
integer, intent(in):: km, is, ie, js, je, isd, ied, jsd, jed, id_te
integer, intent(in):: sphum, liq_wat, ice_wat, rainwat, snowwat, graupel, nwat
integer, intent(in):: sphum, liq_wat, ice_wat, rainwat, snowwat, graupel, hailwat, nwat
real, intent(inout), dimension(isd:ied,jsd:jed,km):: ua, va
real, intent(in), dimension(isd:ied,jsd:jed,km):: pt, delp
real, intent(in), dimension(isd:ied,jsd:jed,km,*):: q
Expand Down Expand Up @@ -966,7 +976,7 @@ subroutine compute_total_energy(is, ie, js, je, isd, ied, jsd, jed, km, &
#ifdef MULTI_GASES
!$OMP num_gas, &
#endif
!$OMP q,nwat,liq_wat,rainwat,ice_wat,snowwat,graupel,sphum) &
!$OMP q,nwat,liq_wat,rainwat,ice_wat,snowwat,graupel,hailwat,sphum) &
!$OMP private(phiz, tv, cvm, qd)
do j=js,je

Expand Down Expand Up @@ -1012,7 +1022,7 @@ subroutine compute_total_energy(is, ie, js, je, isd, ied, jsd, jed, km, &
do k=1,km
#ifdef MOIST_CAPPA
call moist_cv(is,ie,isd,ied,jsd,jed, km, j, k, nwat, sphum, liq_wat, rainwat, &
ice_wat, snowwat, graupel, q, qd, cvm)
ice_wat, snowwat, graupel, hailwat, q, qd, cvm)
#endif
do i=is,ie
#ifdef MOIST_CAPPA
Expand Down Expand Up @@ -3408,9 +3418,9 @@ end subroutine mappm
!! including the heating capacity of water vapor and condensates.
!>@details See \cite emanuel1994atmospheric for information on variable heat capacities.
subroutine moist_cv(is,ie, isd,ied, jsd,jed, km, j, k, nwat, sphum, liq_wat, rainwat, &
ice_wat, snowwat, graupel, q, qd, cvm, t1)
ice_wat, snowwat, graupel, hailwat, q, qd, cvm, t1)
integer, intent(in):: is, ie, isd,ied, jsd,jed, km, nwat, j, k
integer, intent(in):: sphum, liq_wat, rainwat, ice_wat, snowwat, graupel
integer, intent(in):: sphum, liq_wat, rainwat, ice_wat, snowwat, graupel, hailwat
#ifdef MULTI_GASES
real, intent(in), dimension(isd:ied,jsd:jed,km,num_gas):: q
#else
Expand Down Expand Up @@ -3499,6 +3509,18 @@ subroutine moist_cv(is,ie, isd,ied, jsd,jed, km, j, k, nwat, sphum, liq_wat, rai
cvm(i) = (1.-(qv(i)+qd(i)))*cv_air*vicvqd(q(i,j,k,1:num_gas)) + qv(i)*cv_vap + ql(i)*c_liq + qs(i)*c_ice
#else
cvm(i) = (1.-(qv(i)+qd(i)))*cv_air + qv(i)*cv_vap + ql(i)*c_liq + qs(i)*c_ice
#endif
enddo
case(7)
do i=is,ie
qv(i) = q(i,j,k,sphum)
ql(i) = q(i,j,k,liq_wat) + q(i,j,k,rainwat)
qs(i) = q(i,j,k,ice_wat) + q(i,j,k,snowwat) + q(i,j,k,graupel) + q(i,j,k,hailwat)
qd(i) = ql(i) + qs(i)
#ifdef MULTI_GASES
cvm(i) = (1.-(qv(i)+qd(i)))*cv_air*vicvqd(q(i,j,k,1:num_gas)) + qv(i)*cv_vap + ql(i)*c_liq + qs(i)*c_ice
#else
cvm(i) = (1.-(qv(i)+qd(i)))*cv_air + qv(i)*cv_vap + ql(i)*c_liq + qs(i)*c_ice
#endif
enddo
case default
Expand All @@ -3518,10 +3540,10 @@ end subroutine moist_cv
!>@brief The subroutine 'moist_cp' computes the FV3-consistent moist heat capacity under constant pressure,
!! including the heating capacity of water vapor and condensates.
subroutine moist_cp(is,ie, isd,ied, jsd,jed, km, j, k, nwat, sphum, liq_wat, rainwat, &
ice_wat, snowwat, graupel, q, qd, cpm, t1)
ice_wat, snowwat, graupel, hailwat, q, qd, cpm, t1)

integer, intent(in):: is, ie, isd,ied, jsd,jed, km, nwat, j, k
integer, intent(in):: sphum, liq_wat, rainwat, ice_wat, snowwat, graupel
integer, intent(in):: sphum, liq_wat, rainwat, ice_wat, snowwat, graupel, hailwat
#ifdef MULTI_GASES
real, intent(in), dimension(isd:ied,jsd:jed,km,num_gas):: q
#else
Expand Down Expand Up @@ -3617,6 +3639,20 @@ subroutine moist_cp(is,ie, isd,ied, jsd,jed, km, j, k, nwat, sphum, liq_wat, rai
cpm(i) = (1.-(qv(i)+qd(i)))*cp_air + qv(i)*cp_vapor + ql(i)*c_liq + qs(i)*c_ice
#endif
enddo

case(7)
do i=is,ie
qv(i) = q(i,j,k,sphum)
ql(i) = q(i,j,k,liq_wat) + q(i,j,k,rainwat)
qs(i) = q(i,j,k,ice_wat) + q(i,j,k,snowwat) + q(i,j,k,graupel) + q(i,j,k,hailwat)
qd(i) = ql(i) + qs(i)
#ifdef MULTI_GASES
cpm(i) = (1.-(qv(i)+qd(i)))*cp_air*vicpqd(q(i,j,k,:)) + qv(i)*cp_vapor + ql(i)*c_liq + qs(i)*c_ice
#else
cpm(i) = (1.-(qv(i)+qd(i)))*cp_air + qv(i)*cp_vapor + ql(i)*c_liq + qs(i)*c_ice
#endif
enddo

case default
!call mpp_error (NOTE, 'fv_mapz::moist_cp - using default cp_air')
do i=is,ie
Expand Down
Loading