diff --git a/.gitmodules b/.gitmodules
index c57c1a8235..dc3866ecd8 100644
--- a/.gitmodules
+++ b/.gitmodules
@@ -36,7 +36,7 @@
[submodule "atmos_phys"]
path = src/atmos_phys
url = https://github.com/ESCOMP/atmospheric_physics
- fxtag = atmos_phys0_07_001
+ fxtag = atmos_phys0_08_000
fxrequired = AlwaysRequired
fxDONOTUSEurl = https://github.com/ESCOMP/atmospheric_physics
diff --git a/bld/build-namelist b/bld/build-namelist
index 333d1eab40..3b2182bdc7 100755
--- a/bld/build-namelist
+++ b/bld/build-namelist
@@ -3658,7 +3658,6 @@ if (!$simple_phys) {
add_default($nl, 'zmconv_c0_ocn');
add_default($nl, 'zmconv_ke');
add_default($nl, 'zmconv_ke_lnd');
- add_default($nl, 'zmconv_org');
add_default($nl, 'zmconv_num_cin');
add_default($nl, 'zmconv_dmpdz');
add_default($nl, 'zmconv_tiedke_add');
diff --git a/bld/config_files/definition.xml b/bld/config_files/definition.xml
index 0b7b6bca45..9feec64a03 100644
--- a/bld/config_files/definition.xml
+++ b/bld/config_files/definition.xml
@@ -83,10 +83,6 @@ Switch to turn on UNICON package: 0 => off, 1 => on
Switch to turn on/off advecting CLUBB moments: 0 => no, 1 => yes
-
-Switch to turn on/off parameterization for sub-grid scale convective organization for the ZM deep convective scheme based
-on Mapes and Neale (2011): 0 => no, 1 => yes
-
PBL package: uw (University of Washington), hb (Holtslag and Boville), hbr
(Holtslag, Boville, and Rasch), clubb_sgs, spcam_sam1om, spcam_m2005, none.
diff --git a/bld/configure b/bld/configure
index 8f9f435c8c..82bae67421 100755
--- a/bld/configure
+++ b/bld/configure
@@ -116,8 +116,6 @@ OPTIONS
The user does not need to set this if one of the waccm chemistry options
is chosen.
-waccmx Build CAM/WACCM with WACCM upper Thermosphere/Ionosphere extended package
- -zmconv_org Include parameterization for sub-grid scale convective organization for the ZM deep convective scheme based
- on Mapes and Neale (2011)
Options relevent to SCAM mode:
@@ -317,7 +315,6 @@ GetOptions(
"version" => \$opts{'version'},
"waccm_phys" => \$opts{'waccm_phys'},
"waccmx" => \$opts{'waccmx'},
- "zmconv_org" => \$opts{'zmconv_org'},
) or usage();
# Give usage message.
@@ -939,16 +936,6 @@ if (defined $opts{'clubb_opts'}) {
my $clubb_do_adv = $cfg_ref->get('clubb_do_adv');
if ($print>=2) { print "clubb_do_adv: $clubb_do_adv$eol"; }
-#-----------------------------------------------------------------------------------------------
-# ZM convective organization
-
-if (defined $opts{'zmconv_org'}) {
- $cfg_ref->set('zmconv_org', $opts{'zmconv_org'});
-}
-
-my $zmconv_org = $cfg_ref->get('zmconv_org');
-if ($print>=2) { print "zmconv_org: $zmconv_org$eol"; }
-
#-----------------------------------------------------------------------------------------------
# Macro-physics package
@@ -1569,11 +1556,6 @@ else {
if ($print>=2) { print "Advected constituents added by $microphys_pkg microphysics: 10$eol"; }
}
- if ($zmconv_org == 1 ) {
- $nadv += 1;
- if ($print>=2) { print "Advected constituents added by $microphys_pkg microphysics: 8$eol"; }
- }
-
if ($clubb_do_adv) {
$nadv += 9;
if ($print>=2) { print "Advected constituents added by $microphys_pkg microphysics: 8$eol"; }
@@ -2344,6 +2326,8 @@ sub write_filepath
print $fh "$camsrcdir/src/atmos_phys/schemes/check_energy\n";
print $fh "$camsrcdir/src/atmos_phys/schemes/utilities\n";
+ print $fh "$camsrcdir/src/atmos_phys/schemes/cloud_fraction\n";
+
# Dynamics package and test utilities
print $fh "$camsrcdir/src/dynamics/$dyn\n";
if($dyn eq 'se') {
diff --git a/bld/namelist_files/namelist_defaults_cam.xml b/bld/namelist_files/namelist_defaults_cam.xml
index c7239134ba..53af1a16cf 100644
--- a/bld/namelist_files/namelist_defaults_cam.xml
+++ b/bld/namelist_files/namelist_defaults_cam.xml
@@ -2936,9 +2936,6 @@
5.0E-6
5.0E-6
- .false.
- .true.
-
5
1
1
diff --git a/bld/namelist_files/namelist_definition.xml b/bld/namelist_files/namelist_definition.xml
index 548891586c..d07e4c403f 100644
--- a/bld/namelist_files/namelist_definition.xml
+++ b/bld/namelist_files/namelist_definition.xml
@@ -3294,13 +3294,6 @@ Tunable evaporation efficiency in ZM deep convection scheme.
Default: set by build-namelist
-
-Include organization parameterization in ZM. This value is set to true automatically
-if -zmconv_org is set in configure.
-Default: .false., unless -zmconv_org set in configure
-
-
The number of negative buoyancy regions that are allowed before the convection top and CAPE calculations are completed.
diff --git a/doc/ChangeLog b/doc/ChangeLog
index 0aaa217c2a..1ad2135740 100644
--- a/doc/ChangeLog
+++ b/doc/ChangeLog
@@ -1,3 +1,91 @@
+
+===============================================================
+
+Tag name: cam6_4_065
+Originator(s): cacraig
+Date: Feb 11, 2025
+One-line Summary: ZM CCPP'ization round 4 (completes CCPP conversion of ZM)
+
+Github PR URL: https://github.com/ESCOMP/CAM/pull/1218
+
+Purpose of changes (include the issue number and title text for each relevant GitHub issue):
+ - Convert ZM to CCPP and move into atmospheric_physics github repo: https://github.com/ESCOMP/CAM/issues/873
+
+Describe any changes made to build system: N/A
+
+Describe any changes made to the namelist:
+ - Removed zmconv_org namelist as that partially implemented capability has been removed
+
+List any changes to the defaults for the boundary datasets: N/A
+
+Describe any substantial timing or memory changes: N/A
+
+Code reviewed by: nusbaume, jimmielin
+
+List all files eliminated: N/A
+D src/physics/cam/wv_sat_methods.F90
+D src/physics/cam/wv_saturation.F90
+D src/utils/error_messages.F90
+D src/utils/namelist_utils.F90
+ - Moved to atmospheric_physics (and currently reside in the to_be_ccppized directory)
+
+List all files added and what they do: N/A
+
+List all existing files that have been modified, and describe the changes:
+M .gitmodules
+ - update atmospheric_physics to bring in ZM changes
+
+M bld/build-namelist
+M bld/config_files/definition.xml
+M bld/configure
+M bld/namelist_files/namelist_defaults_cam.xml
+M bld/namelist_files/namelist_definition.xml
+ - remove zmconv_org namelist
+
+M src/physics/cam/cloud_fraction.F90
+ - moved cldfrc_fice to atmospheric_physics and ccppized
+
+M src/physics/cam/clubb_intr.F90
+ - removed difzm declarations as no longer needed
+
+M src/physics/cam/convect_shallow.F90
+M src/physics/cam/macrop_driver.F90
+M src/physics/cam/physpkg.F90
+M src/physics/cam/rk_stratiform.F90
+M src/physics/cam/zm_conv_intr.F90
+M src/physics/cam7/physpkg.F90
+M src/physics/simple/physpkg.F90
+ - various mods to get this to work with the routines that are ccppized
+
+M src/utils/cam_ccpp/ccpp_constituent_prop_mod.F90
+ - Add routine:
+ ccp_set_standard_name to set constituent's standard name
+ ccp_is_dry to return whether species is dry
+ ccp_set_dry to set constituent's dry property based on what is passed in
+
+If there were any failures reported from running test_driver.sh on any test
+platform, and checkin with these failures has been OK'd by the gatekeeper,
+then copy the lines from the td.*.status files for the failed tests to the
+appropriate machine below. All failed tests must be justified.
+
+derecho/intel/aux_cam: all BFB, except:
+ ERP_Ln9.f09_f09_mg17.FCSD_HCO.derecho_intel.cam-outfrq9s (Overall: FAIL)
+ SMS_Ld1.f09_f09_mg17.FCHIST_GC.derecho_intel.cam-outfrq1d (Overall: DIFF)
+ - pre-existing failure due to HEMCO not having reproducible results (issues #1018 and #856)
+
+ SMS_D_Ln9.f19_f19_mg17.FXHIST.derecho_intel.cam-outfrq9s_amie (Overall: FAIL)
+ SMS_D_Ln9_P1280x1.ne0CONUSne30x8_ne0CONUSne30x8_mt12.FCHIST.derecho_intel.cam-outfrq9s (Overall: FAIL)
+ - pre-existing failures due to build-namelist error requiring CLM/CTSM external update
+
+derecho/nvhpc/aux_cam: BFB
+
+izumi/nag/aux_cam: BFB
+
+izumi/gnu/aux_cam:
+ ERC_D_Ln9.f10_f10_mg37.QPSPCAMM.izumi_gnu.cam-outfrq3s (Overall: FAIL) details:
+ - New failure, but since SPCAM is being removed(PR 1217) will create this new pre-existing failure
+
+===============================================================
===============================================================
Tag name: cam6_4_064
diff --git a/src/physics/cam/cloud_fraction.F90 b/src/physics/cam/cloud_fraction.F90
index 3285862fae..821b84695d 100644
--- a/src/physics/cam/cloud_fraction.F90
+++ b/src/physics/cam/cloud_fraction.F90
@@ -5,7 +5,7 @@ module cloud_fraction
use shr_kind_mod, only: r8 => shr_kind_r8
use ppgrid, only: pcols, pver, pverp
- use ref_pres, only: pref_mid
+ use ref_pres, only: pref_mid
use spmd_utils, only: masterproc
use cam_logfile, only: iulog
use cam_abortutils, only: endrun
@@ -22,7 +22,6 @@ module cloud_fraction
cldfrc_init, &! Inititialization of cloud_fraction run-time parameters
cldfrc_getparams, &! public access of tuning parameters
cldfrc, &! Computation of cloud fraction
- cldfrc_fice, &! Calculate fraction of condensate in ice phase (radiation partitioning)
dp1, &! parameter for deep convection cloud fraction needed in clubb_intr
dp2 ! parameter for deep convection cloud fraction needed in clubb_intr
@@ -32,9 +31,9 @@ module cloud_fraction
! Top level
integer :: top_lev = 1
- ! Physics buffer indices
- integer :: sh_frac_idx = 0
- integer :: dp_frac_idx = 0
+ ! Physics buffer indices
+ integer :: sh_frac_idx = 0
+ integer :: dp_frac_idx = 0
! Namelist variables
logical :: cldfrc_freeze_dry ! switch for Vavrus correction
@@ -154,8 +153,8 @@ subroutine cldfrc_register
!-----------------------------------------------------------------------
- call pbuf_add_field('SH_FRAC', 'physpkg', dtype_r8, (/pcols,pver/), sh_frac_idx)
- call pbuf_add_field('DP_FRAC', 'physpkg', dtype_r8, (/pcols,pver/), dp_frac_idx)
+ call pbuf_add_field('SH_FRAC', 'physpkg', dtype_r8, (/pcols,pver/), sh_frac_idx)
+ call pbuf_add_field('DP_FRAC', 'physpkg', dtype_r8, (/pcols,pver/), dp_frac_idx)
end subroutine cldfrc_register
@@ -215,7 +214,7 @@ subroutine cldfrc_init
inversion_cld_off = .false.
endif
- if ( masterproc ) then
+ if ( masterproc ) then
write(iulog,*)'tuning parameters cldfrc_init: inversion_cld_off',inversion_cld_off
write(iulog,*)'tuning parameters cldfrc_init: dp1',dp1,'dp2',dp2,'sh1',sh1,'sh2',sh2
if (shallow_scheme .ne. 'UW' .or. shallow_scheme .eq. 'SPCAM_m2005' ) then
@@ -249,38 +248,38 @@ subroutine cldfrc(lchnk ,ncol , pbuf, &
cmfmc ,cmfmc2 ,landfrac,snowh ,concld ,cldst , &
ts ,sst ,ps ,zdu ,ocnfrac ,&
rhu00 ,cldice ,icecldf ,liqcldf ,relhum ,dindex )
- !-----------------------------------------------------------------------
- !
- ! Purpose:
- ! Compute cloud fraction
- !
- !
- ! Method:
+ !-----------------------------------------------------------------------
+ !
+ ! Purpose:
+ ! Compute cloud fraction
+ !
+ !
+ ! Method:
! This calculate cloud fraction using a relative humidity threshold
- ! The threshold depends upon pressure, and upon the presence or absence
- ! of convection as defined by a reasonably large vertical mass flux
+ ! The threshold depends upon pressure, and upon the presence or absence
+ ! of convection as defined by a reasonably large vertical mass flux
! entering that layer from below.
- !
+ !
! Author: Many. Last modified by Jim McCaa
- !
+ !
!-----------------------------------------------------------------------
use cam_history, only: outfld
use physconst, only: cappa, gravit, rair, tmelt
use wv_saturation, only: qsat, qsat_water, svp_ice_vect
use phys_grid, only: get_rlat_all_p, get_rlon_all_p
-
+
!RBN - Need this to write shallow,deep fraction to phys buffer.
!PJR - we should probably make seperate modules for determining convective
! clouds and make this one just responsible for relative humidity clouds
-
+
use physics_buffer, only: physics_buffer_desc, pbuf_get_field
! Arguments
integer, intent(in) :: lchnk ! chunk identifier
integer, intent(in) :: ncol ! number of atmospheric columns
integer, intent(in) :: dindex ! 0 or 1 to perturb rh
-
+
type(physics_buffer_desc), pointer :: pbuf(:)
real(r8), intent(in) :: pmid(pcols,pver) ! midpoint pressures
real(r8), intent(in) :: temp(pcols,pver) ! temperature
@@ -307,7 +306,7 @@ subroutine cldfrc(lchnk ,ncol , pbuf, &
real(r8), intent(out) :: clc(pcols) ! column convective cloud amount
real(r8), intent(out) :: cldst(pcols,pver) ! cloud fraction
real(r8), intent(out) :: rhu00(pcols,pver) ! RH threshold for cloud
- real(r8), intent(out) :: relhum(pcols,pver) ! RH
+ real(r8), intent(out) :: relhum(pcols,pver) ! RH
real(r8), intent(out) :: icecldf(pcols,pver) ! ice cloud fraction
real(r8), intent(out) :: liqcldf(pcols,pver) ! liquid cloud fraction (combined into cloud)
@@ -376,7 +375,7 @@ subroutine cldfrc(lchnk ,ncol , pbuf, &
! The idea is that the RH limits for condensation are strict only for
! water saturation
!
- ! Ice clouds are formed by explicit parameterization of ice nucleation.
+ ! Ice clouds are formed by explicit parameterization of ice nucleation.
! Closure for ice cloud fraction is done on available cloud ice, such that
! the in-cloud ice content matches an empirical fit
! thus, icecldf = min(cldice/icicval,1) where icicval = f(temp,cldice,numice)
@@ -385,17 +384,17 @@ subroutine cldfrc(lchnk ,ncol , pbuf, &
! No dA/dt term for ice?
!
! There are three co-existing cloud types: convective, inversion related low-level
- ! stratocumulus, and layered cloud (based on relative humidity). Layered and
- ! stratocumulus clouds do not compete with convective cloud for which one creates
- ! the most cloud. They contribute collectively to the total grid-box average cloud
- ! amount. This is reflected in the way in which the total cloud amount is evaluated
+ ! stratocumulus, and layered cloud (based on relative humidity). Layered and
+ ! stratocumulus clouds do not compete with convective cloud for which one creates
+ ! the most cloud. They contribute collectively to the total grid-box average cloud
+ ! amount. This is reflected in the way in which the total cloud amount is evaluated
! (a sum as opposed to a logical "or" operation)
!
!==================================================================================
! set defaults for rhu00
rhu00(:,:) = 2.0_r8
! define rh perturbation in order to estimate rhdfda
- rhpert = 0.01_r8
+ rhpert = 0.01_r8
!set Wang and Sassen IWC paramters
a=26.87_r8
@@ -460,7 +459,7 @@ subroutine cldfrc(lchnk ,ncol , pbuf, &
!
! Estimate of local convective cloud cover based on convective mass flux
- ! Modify local large-scale relative humidity to account for presence of
+ ! Modify local large-scale relative humidity to account for presence of
! convective cloud when evaluating relative humidity based layered cloud amount
!
concld(:ncol,top_lev:pver) = 0.0_r8
@@ -468,7 +467,7 @@ subroutine cldfrc(lchnk ,ncol , pbuf, &
! cloud mass flux in SI units of kg/m2/s; should produce typical numbers of 20%
! shallow and deep convective cloudiness are evaluated separately (since processes
! are evaluated separately) and summed
- !
+ !
#ifndef PERGRO
do k=top_lev,pver
do i=1,ncol
@@ -488,7 +487,7 @@ subroutine cldfrc(lchnk ,ncol , pbuf, &
! ****** Compute layer cloudiness ******
!
!====================================================================
- ! Begin the evaluation of layered cloud amount based on (modified) RH
+ ! Begin the evaluation of layered cloud amount based on (modified) RH
!====================================================================
!
numkcld = pver
@@ -517,7 +516,7 @@ subroutine cldfrc(lchnk ,ncol , pbuf, &
! SJV: decrease cloud amount if very low water vapor content
! (thus very cold): "freeze dry"
if (cldfrc_freeze_dry) then
- rhcloud(i,k) = rhcloud(i,k)*max(0.15_r8,min(1.0_r8,q(i,k)/0.0030_r8))
+ rhcloud(i,k) = rhcloud(i,k)*max(0.15_r8,min(1.0_r8,q(i,k)/0.0030_r8))
endif
else if ( pmid(i,k).lt.premit ) then
@@ -537,7 +536,7 @@ subroutine cldfrc(lchnk ,ncol , pbuf, &
! linear rh threshold transition between thresholds for low & high cloud
!
rhwght = (premib-(max(pmid(i,k),premit)))/(premib-premit)
-
+
if (land(i) .and. (snowh(i) <= 0.000001_r8)) then
rhlim = rhminh*rhwght + (rhminl - rhminl_adj_land)*(1.0_r8-rhwght)
else
@@ -591,7 +590,7 @@ subroutine cldfrc(lchnk ,ncol , pbuf, &
!--------ICE CLOUD OPTION 3--------Wood & Field 2000 (JAS)
! eq 6: cloud fraction = 1 - exp (-K * qc/qsati)
-
+
icecldf(i,k)=1._r8 - exp(-Kc*cldice(i,k)/(qs(i,k)*(esi(i,k)/esl(i,k))))
icecldf(i,k)=max(0._r8,min(icecldf(i,k),1._r8))
else
@@ -634,7 +633,7 @@ subroutine cldfrc(lchnk ,ncol , pbuf, &
cloud(i,k) = rhcloud(i,k)
end if
end do
- end do
+ end do
!
! Add in the marine strat
! MARINE STRATUS SHOULD BE A SPECIAL CASE OF LAYERED CLOUD
@@ -644,20 +643,20 @@ subroutine cldfrc(lchnk ,ncol , pbuf, &
!===================================================================================
!
! SOME OBSERVATIONS ABOUT THE FOLLOWING SECTION OF CODE (missed in earlier look)
- ! K700 IS SET AS A CONSTANT BASED ON HYBRID COORDINATE: IT DOES NOT DEPEND ON
- ! LOCAL PRESSURE; THERE IS NO PRESSURE RAMP => LOOKS LEVEL DEPENDENT AND
+ ! K700 IS SET AS A CONSTANT BASED ON HYBRID COORDINATE: IT DOES NOT DEPEND ON
+ ! LOCAL PRESSURE; THERE IS NO PRESSURE RAMP => LOOKS LEVEL DEPENDENT AND
! DISCONTINUOUS IN SPACE (I.E., STRATUS WILL END SUDDENLY WITH NO TRANSITION)
!
! IT APPEARS THAT STRAT IS EVALUATED ACCORDING TO KLEIN AND HARTMANN; HOWEVER,
! THE ACTUAL STRATUS AMOUNT (CLDST) APPEARS TO DEPEND DIRECTLY ON THE RH BELOW
- ! THE STRONGEST PART OF THE LOW LEVEL INVERSION.
+ ! THE STRONGEST PART OF THE LOW LEVEL INVERSION.
!PJR answers: 1) the rh limitation is a physical/mathematical limitation
! cant have more cloud than there is RH
! allowed the cloud to exist two layers below the inversion
! because the numerics frequently make 50% relative humidity
! in level below the inversion which would allow no cloud
! 2) since the cloud is only allowed over ocean, it should
- ! be very insensitive to surface pressure (except due to
+ ! be very insensitive to surface pressure (except due to
! spectral ringing, which also causes so many other problems
! I didnt worry about it.
!
@@ -738,77 +737,4 @@ end subroutine cldfrc
!================================================================================================
- subroutine cldfrc_fice(ncol, t, fice, fsnow)
-!
-! Compute the fraction of the total cloud water which is in ice phase.
-! The fraction depends on temperature only.
-! This is the form that was used for radiation, the code came from cldefr originally
-!
-! Author: B. A. Boville Sept 10, 2002
-! modified: PJR 3/13/03 (added fsnow to ascribe snow production for convection )
-!-----------------------------------------------------------------------
- use physconst, only: tmelt
-
-! Arguments
- integer, intent(in) :: ncol ! number of active columns
- real(r8), intent(in) :: t(:,:) ! temperature
-
- real(r8), intent(out) :: fice(:,:) ! Fractional ice content within cloud
- real(r8), intent(out) :: fsnow(:,:) ! Fractional snow content for convection
-
-! Local variables
- real(r8) :: tmax_fice ! max temperature for cloud ice formation
- real(r8) :: tmin_fice ! min temperature for cloud ice formation
- real(r8) :: tmax_fsnow ! max temperature for transition to convective snow
- real(r8) :: tmin_fsnow ! min temperature for transition to convective snow
-
- integer :: i,k ! loop indexes
-
-!-----------------------------------------------------------------------
-
- tmax_fice = tmelt - 10._r8 ! max temperature for cloud ice formation
- tmin_fice = tmax_fice - 30._r8 ! min temperature for cloud ice formation
- tmax_fsnow = tmelt ! max temperature for transition to convective snow
- tmin_fsnow = tmelt - 5._r8 ! min temperature for transition to convective snow
-
- fice(:,:top_lev-1) = 0._r8
- fsnow(:,:top_lev-1) = 0._r8
-
-! Define fractional amount of cloud that is ice
- do k=top_lev,pver
- do i=1,ncol
-
-! If warmer than tmax then water phase
- if (t(i,k) > tmax_fice) then
- fice(i,k) = 0.0_r8
-
-! If colder than tmin then ice phase
- else if (t(i,k) < tmin_fice) then
- fice(i,k) = 1.0_r8
-
-! Otherwise mixed phase, with ice fraction decreasing linearly from tmin to tmax
- else
- fice(i,k) =(tmax_fice - t(i,k)) / (tmax_fice - tmin_fice)
- end if
-
-! snow fraction partitioning
-
-! If warmer than tmax then water phase
- if (t(i,k) > tmax_fsnow) then
- fsnow(i,k) = 0.0_r8
-
-! If colder than tmin then ice phase
- else if (t(i,k) < tmin_fsnow) then
- fsnow(i,k) = 1.0_r8
-
-! Otherwise mixed phase, with ice fraction decreasing linearly from tmin to tmax
- else
- fsnow(i,k) =(tmax_fsnow - t(i,k)) / (tmax_fsnow - tmin_fsnow)
- end if
-
- end do
- end do
-
- end subroutine cldfrc_fice
-
end module cloud_fraction
diff --git a/src/physics/cam/clubb_intr.F90 b/src/physics/cam/clubb_intr.F90
index 27acf0f14f..9bbed56277 100644
--- a/src/physics/cam/clubb_intr.F90
+++ b/src/physics/cam/clubb_intr.F90
@@ -478,7 +478,6 @@ module clubb_intr
integer :: &
dlfzm_idx = -1, & ! ZM detrained convective cloud water mixing ratio.
- difzm_idx = -1, & ! ZM detrained convective cloud ice mixing ratio.
dnlfzm_idx = -1, & ! ZM detrained convective cloud water num concen.
dnifzm_idx = -1 ! ZM detrained convective cloud ice num concen.
@@ -2525,7 +2524,6 @@ subroutine clubb_tend_cam( state, ptend_all, pbuf, hdtime, &
! ZM microphysics
real(r8), pointer :: dlfzm(:,:) ! ZM detrained convective cloud water mixing ratio.
- real(r8), pointer :: difzm(:,:) ! ZM detrained convective cloud ice mixing ratio.
real(r8), pointer :: dnlfzm(:,:) ! ZM detrained convective cloud water num concen.
real(r8), pointer :: dnifzm(:,:) ! ZM detrained convective cloud ice num concen.
diff --git a/src/physics/cam/convect_shallow.F90 b/src/physics/cam/convect_shallow.F90
index 902187eb24..e757e1eef4 100644
--- a/src/physics/cam/convect_shallow.F90
+++ b/src/physics/cam/convect_shallow.F90
@@ -14,10 +14,12 @@ module convect_shallow
use physconst, only : cpair, zvir
use ppgrid, only : pver, pcols, pverp
use zm_conv_evap, only : zm_conv_evap_run
- use zm_conv_intr, only : zmconv_ke, zmconv_ke_lnd, zmconv_org
+ use zm_conv_intr, only : zmconv_ke, zmconv_ke_lnd
use cam_history, only : outfld, addfld, horiz_only
use cam_logfile, only : iulog
use phys_control, only : phys_getopts
+ use cloud_fraction_fice, only: cloud_fraction_fice_run
+ use ref_pres, only: trop_cloud_top_lev
implicit none
private
@@ -473,10 +475,21 @@ subroutine convect_shallow_tend( ztodt , cmfmc , &
real(r8), pointer, dimension(:,:) :: cmfmc2 ! (pcols,pverp) Updraft mass flux by shallow convection [ kg/s/m2 ]
real(r8), pointer, dimension(:,:) :: sh_e_ed_ratio ! (pcols,pver) fer/(fer+fdr) from uwschu
+ real(r8), dimension(pcols,pver) :: fsnow_conv
+ real(r8), dimension(pcols,pver) :: fice
+
logical :: lq(pcnst)
type(unicon_out_t) :: unicon_out
+ character(len=40) :: scheme_name
+ character(len=16) :: macrop_scheme
+ character(len=512):: errmsg
+ integer :: errflg
+ integer :: top_lev
+
+
+
! ----------------------- !
! Main Computation Begins !
! ----------------------- !
@@ -872,16 +885,25 @@ subroutine convect_shallow_tend( ztodt , cmfmc , &
tend_s_snwprd(:,:) = 0._r8
tend_s_snwevmlt(:,:) = 0._r8
snow(:) = 0._r8
+ fice(:,:) = 0._r8
+ fsnow_conv(:,:) = 0._r8
!REMOVECAM_END
+ top_lev = 1
+ call phys_getopts (macrop_scheme_out = macrop_scheme)
+ if ( .not. (macrop_scheme == "rk" .or. macrop_scheme == "SPCAM_sam1mom")) top_lev = trop_cloud_top_lev
+
+ call cloud_fraction_fice_run(ncol, state1%t(1:ncol,:), tmelt, top_lev, pver, fice(1:ncol,:), fsnow_conv(1:ncol,:))
+
call zm_conv_evap_run(state1%ncol, pver, pverp, &
gravit, latice, latvap, tmelt, &
- cpair, zmconv_ke, zmconv_ke_lnd, zmconv_org, &
+ cpair, zmconv_ke, zmconv_ke_lnd, &
state1%t(:ncol,:),state1%pmid(:ncol,:),state1%pdel(:ncol,:),state1%q(:ncol,:pver,1), &
landfracdum(:ncol), &
ptend_loc%s(:ncol,:), tend_s_snwprd(:ncol,:), tend_s_snwevmlt(:ncol,:), ptend_loc%q(:ncol,:pver,1), &
rprdsh(:ncol,:), cld(:ncol,:), ztodt, &
- precc(:ncol), snow(:ncol), ntprprd(:ncol,:), ntsnprd(:ncol,:), flxprec(:ncol,:), flxsnow(:ncol,:) )
+ precc(:ncol), snow(:ncol), ntprprd(:ncol,:), ntsnprd(:ncol,:), fsnow_conv(:ncol,:), flxprec(:ncol,:), flxsnow(:ncol,:),&
+ scheme_name, errmsg, errflg)
! ---------------------------------------------- !
! record history variables from zm_conv_evap_run !
diff --git a/src/physics/cam/macrop_driver.F90 b/src/physics/cam/macrop_driver.F90
index d381387bfc..26217c2a8c 100644
--- a/src/physics/cam/macrop_driver.F90
+++ b/src/physics/cam/macrop_driver.F90
@@ -87,7 +87,6 @@ module macrop_driver
integer :: &
dlfzm_idx = -1, & ! ZM detrained convective cloud water mixing ratio.
- difzm_idx = -1, & ! ZM detrained convective cloud ice mixing ratio.
dnlfzm_idx = -1, & ! ZM detrained convective cloud water num concen.
dnifzm_idx = -1 ! ZM detrained convective cloud ice num concen.
@@ -395,7 +394,8 @@ subroutine macrop_driver_tend( &
! !
!-------------------------------------------------------- !
- use cloud_fraction, only: cldfrc, cldfrc_fice
+ use cloud_fraction, only: cldfrc
+ use cloud_fraction_fice, only: cloud_fraction_fice_run
use physics_types, only: physics_state, physics_ptend
use physics_types, only: physics_ptend_init, physics_update
use physics_types, only: physics_ptend_sum, physics_state_copy
@@ -486,7 +486,6 @@ subroutine macrop_driver_tend( &
! ZM microphysics
real(r8), pointer :: dlfzm(:,:) ! ZM detrained convective cloud water mixing ratio.
- real(r8), pointer :: difzm(:,:) ! ZM detrained convective cloud ice mixing ratio.
real(r8), pointer :: dnlfzm(:,:) ! ZM detrained convective cloud water num concen.
real(r8), pointer :: dnifzm(:,:) ! ZM detrained convective cloud ice num concen.
@@ -872,8 +871,8 @@ subroutine macrop_driver_tend( &
fice(:,:) = 0._r8
fsnow(:,:) = 0._r8
!REMOVECAM_END
- call cldfrc_fice( ncol, state_loc%t(:ncol,:), fice(:ncol,:), fsnow(:ncol,:) )
+ call cloud_fraction_fice_run(ncol, state_loc%t(:ncol,:), tmelt, top_lev, pver, fice(:ncol,:), fsnow(:ncol,:))
lq(:) = .FALSE.
diff --git a/src/physics/cam/physpkg.F90 b/src/physics/cam/physpkg.F90
index 40bcfc3482..c9822b18fc 100644
--- a/src/physics/cam/physpkg.F90
+++ b/src/physics/cam/physpkg.F90
@@ -731,6 +731,7 @@ subroutine phys_init( phys_state, phys_tend, pbuf2d, cam_in, cam_out )
use co2_cycle, only: co2_init, co2_transport
use convect_deep, only: convect_deep_init
use convect_shallow, only: convect_shallow_init
+ use constituents, only: cnst_get_ind
use cam_diagnostics, only: diag_init
use gw_drag, only: gw_init
use radheat, only: radheat_init
@@ -787,7 +788,7 @@ subroutine phys_init( phys_state, phys_tend, pbuf2d, cam_in, cam_out )
! local variables
integer :: lchnk
- integer :: ierr
+ integer :: ierr, ixq
logical :: history_budget ! output tendencies and state variables for
! temperature, water vapor, cloud
@@ -966,7 +967,8 @@ subroutine phys_init( phys_state, phys_tend, pbuf2d, cam_in, cam_out )
! Initialize CAM CCPP constituent properties array
! for use in CCPP-ized physics schemes:
- call ccpp_const_props_init()
+ call cnst_get_ind('Q', ixq)
+ call ccpp_const_props_init(ixq)
! Initialize qneg3 and qneg4
call qneg_init()
diff --git a/src/physics/cam/rk_stratiform.F90 b/src/physics/cam/rk_stratiform.F90
index 002300dbfd..148ea3fd28 100644
--- a/src/physics/cam/rk_stratiform.F90
+++ b/src/physics/cam/rk_stratiform.F90
@@ -427,7 +427,8 @@ subroutine rk_stratiform_tend( &
! !
!-------------------------------------------------------- !
- use cloud_fraction, only: cldfrc, cldfrc_fice
+ use cloud_fraction, only: cldfrc
+ use cloud_fraction_fice, only: cloud_fraction_fice_run
use physics_types, only: physics_state, physics_ptend
use physics_types, only: physics_ptend_init, physics_update
use physics_types, only: physics_ptend_sum, physics_state_copy
@@ -440,7 +441,7 @@ subroutine rk_stratiform_tend( &
use phys_control, only: cam_physpkg_is
use tropopause, only: tropopause_find_cam
use phys_grid, only: get_rlat_all_p
- use physconst, only: pi
+ use physconst, only: pi, tmelt
! Arguments
type(physics_state), intent(in) :: state ! State variables
@@ -577,6 +578,9 @@ subroutine rk_stratiform_tend( &
real(r8) :: dlat(pcols)
real(r8), parameter :: rad2deg = 180._r8/pi
+ integer :: top_lev
+
+
! ======================================================================
lchnk = state%lchnk
@@ -812,7 +816,9 @@ subroutine rk_stratiform_tend( &
fice(:,:) = 0._r8
fsnow(:,:) = 0._r8
!REMOVECAM_END
- call cldfrc_fice(ncol, state1%t(1:ncol,:), fice(1:ncol,:), fsnow(1:ncol,:))
+ top_lev = 1
+ call cloud_fraction_fice_run(ncol, state1%t(:ncol,:), tmelt, top_lev, pver, fice(:ncol,:), fsnow(:ncol,:))
+
! Perform repartitioning of stratiform condensate.
! Corresponding heating tendency will be added later.
diff --git a/src/physics/cam/wv_sat_methods.F90 b/src/physics/cam/wv_sat_methods.F90
deleted file mode 100644
index bb2ffeb45b..0000000000
--- a/src/physics/cam/wv_sat_methods.F90
+++ /dev/null
@@ -1,759 +0,0 @@
-module wv_sat_methods
-
-! This portable module contains all CAM methods for estimating
-! the saturation vapor pressure of water.
-!
-! wv_saturation provides CAM-specific interfaces and utilities
-! based on these formulae.
-!
-! Typical usage of this module:
-!
-! Init:
-! call wv_sat_methods_init(r8, , errstring)
-!
-! Get scheme index from a name string:
-! scheme_idx = wv_sat_get_scheme_idx(scheme_name)
-! if (.not. wv_sat_valid_idx(scheme_idx))
-!
-! Get pressures:
-! es = wv_sat_svp_water(t, scheme_idx)
-! es = wv_sat_svp_ice(t, scheme_idx)
-!
-! Use ice/water transition range:
-! es = wv_sat_svp_trice(t, ttrice, scheme_idx)
-!
-! Note that elemental functions cannot be pointed to, nor passed
-! as arguments. If you need to do either, it is recommended to
-! wrap the function so that it can be given an explicit (non-
-! elemental) interface.
-
-implicit none
-private
-save
-
-integer, parameter :: r8 = selected_real_kind(12) ! 8 byte real
-
-integer, parameter :: VLENS = 128 ! vector length for a GPU kernel
-
-real(r8) :: tmelt ! Melting point of water at 1 atm (K)
-real(r8) :: h2otrip ! Triple point temperature of water (K)
-real(r8) :: tboil ! Boiling point of water at 1 atm (K)
-
-real(r8) :: ttrice ! Ice-water transition range
-
-real(r8) :: epsilo ! Ice-water transition range
-real(r8) :: omeps ! 1._r8 - epsilo
-
-! Indices representing individual schemes
-integer, parameter :: Invalid_idx = -1
-integer, parameter :: GoffGratch_idx = 1
-integer, parameter :: MurphyKoop_idx = 2
-integer, parameter :: Bolton_idx = 3
-
-! Index representing the current default scheme.
-integer, parameter :: initial_default_idx = GoffGratch_idx
-integer :: default_idx = initial_default_idx
-
-!$acc declare create (epsilo, tmelt, tboil, omeps, h2otrip, ttrice)
-
-public wv_sat_methods_init
-public wv_sat_get_scheme_idx
-public wv_sat_valid_idx
-
-public wv_sat_set_default
-public wv_sat_reset_default
-
-public wv_sat_qsat_water, wv_sat_qsat_water_vect
-public wv_sat_qsat_ice, wv_sat_qsat_ice_vect
-
-public wv_sat_svp_trans, wv_sat_svp_trans_vect
-
-! pressure -> humidity conversion
-public wv_sat_svp_to_qsat, wv_sat_svp_to_qsat_vect
-
-! Combined qsat operations
-public wv_sat_qsat_trans
-
-public wv_sat_svp_water, wv_sat_svp_water_vect
-public wv_sat_svp_ice, wv_sat_svp_ice_vect
-
-contains
-
-!---------------------------------------------------------------------
-! ADMINISTRATIVE FUNCTIONS
-!---------------------------------------------------------------------
-
-! Get physical constants
-subroutine wv_sat_methods_init(kind, tmelt_in, h2otrip_in, tboil_in, &
- ttrice_in, epsilo_in, errstring)
- integer, intent(in) :: kind
- real(r8), intent(in) :: tmelt_in
- real(r8), intent(in) :: h2otrip_in
- real(r8), intent(in) :: tboil_in
- real(r8), intent(in) :: ttrice_in
- real(r8), intent(in) :: epsilo_in
- character(len=*), intent(out) :: errstring
-
- errstring = ' '
-
- if (kind /= r8) then
- write(errstring,*) 'wv_sat_methods_init: ERROR: ', &
- kind,' was input kind but ',r8,' is internal kind.'
- return
- end if
-
- if (ttrice_in < 0._r8) then
- write(errstring,*) 'wv_sat_methods_init: ERROR: ', &
- ttrice_in,' was input for ttrice, but negative range is invalid.'
- return
- end if
-
- tmelt = tmelt_in
- h2otrip = h2otrip_in
- tboil = tboil_in
- ttrice = ttrice_in
- epsilo = epsilo_in
-
- omeps = 1._r8 - epsilo
-
- !$acc update device (tmelt,h2otrip,tboil,ttrice,epsilo,omeps)
-
-end subroutine wv_sat_methods_init
-
-! Look up index by name.
-pure function wv_sat_get_scheme_idx(name) result(idx)
- character(len=*), intent(in) :: name
- integer :: idx
-
- select case (name)
- case("GoffGratch")
- idx = GoffGratch_idx
- case("MurphyKoop")
- idx = MurphyKoop_idx
- case("Bolton")
- idx = Bolton_idx
- case default
- idx = Invalid_idx
- end select
-
-end function wv_sat_get_scheme_idx
-
-! Check validity of an index from the above routine.
-pure function wv_sat_valid_idx(idx) result(status)
- integer, intent(in) :: idx
- logical :: status
-
- status = (idx /= Invalid_idx)
-
-end function wv_sat_valid_idx
-
-! Set default scheme (otherwise, Goff & Gratch is default)
-! Returns a logical representing success (.true.) or
-! failure (.false.).
-function wv_sat_set_default(name) result(status)
- character(len=*), intent(in) :: name
- logical :: status
-
- ! Don't want to overwrite valid default with invalid,
- ! so assign to temporary and check it first.
- integer :: tmp_idx
-
- tmp_idx = wv_sat_get_scheme_idx(name)
-
- status = wv_sat_valid_idx(tmp_idx)
-
- if (status) default_idx = tmp_idx
-
-end function wv_sat_set_default
-
-! Reset default scheme to initial value.
-! The same thing can be accomplished with wv_sat_set_default;
-! the real reason to provide this routine is to reset the
-! module for testing purposes.
-subroutine wv_sat_reset_default()
-
- default_idx = initial_default_idx
-
-end subroutine wv_sat_reset_default
-
-!---------------------------------------------------------------------
-! UTILITIES
-!---------------------------------------------------------------------
-
-! Get saturation specific humidity given pressure and SVP.
-! Specific humidity is limited to range 0-1.
-function wv_sat_svp_to_qsat(es, p) result(qs)
- real(r8), intent(in) :: es ! SVP
- real(r8), intent(in) :: p ! Current pressure.
- real(r8) :: qs
-
- ! If pressure is less than SVP, set qs to maximum of 1.
- if ( (p - es) <= 0._r8 ) then
- qs = 1.0_r8
- else
- qs = epsilo*es / (p - omeps*es)
- end if
-
-end function wv_sat_svp_to_qsat
-
-! Get saturation specific humidity given pressure and SVP.
-! Specific humidity is limited to range 0-1.
-subroutine wv_sat_svp_to_qsat_vect(es, p, qs, vlen)
-
- integer, intent(in) :: vlen
- real(r8), intent(in) :: es(vlen) ! SVP
- real(r8), intent(in) :: p(vlen) ! Current pressure.
- real(r8), intent(out) :: qs(vlen)
- integer :: i
-
- ! If pressure is less than SVP, set qs to maximum of 1.
-
- !$acc data present (es,p,qs)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i=1,vlen
- if ( (p(i) - es(i)) <= 0._r8 ) then
- qs(i) = 1.0_r8
- else
- qs(i) = epsilo*es(i) / (p(i) - omeps*es(i))
- end if
- end do
- !$acc end parallel
-
- !$acc end data
-end subroutine wv_sat_svp_to_qsat_vect
-
-subroutine wv_sat_qsat_water(t, p, es, qs, idx)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate SVP over water at a given temperature, and then !
- ! calculate and return saturation specific humidity. !
- !------------------------------------------------------------------!
-
- ! Inputs
- real(r8), intent(in) :: t ! Temperature
- real(r8), intent(in) :: p ! Pressure
- ! Outputs
- real(r8), intent(out) :: es ! Saturation vapor pressure
- real(r8), intent(out) :: qs ! Saturation specific humidity
-
- integer, intent(in), optional :: idx ! Scheme index
-
- es = wv_sat_svp_water(t, idx)
-
- qs = wv_sat_svp_to_qsat(es, p)
-
- ! Ensures returned es is consistent with limiters on qs.
- es = min(es, p)
-
-end subroutine wv_sat_qsat_water
-
-subroutine wv_sat_qsat_water_vect(t, p, es, qs, vlen, idx)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate SVP over water at a given temperature, and then !
- ! calculate and return saturation specific humidity. !
- !------------------------------------------------------------------!
- ! Inputs
-
- integer, intent(in) :: vlen
- real(r8), intent(in) :: t(vlen) ! Temperature
- real(r8), intent(in) :: p(vlen) ! Pressure
- ! Outputs
- real(r8), intent(out) :: es(vlen) ! Saturation vapor pressure
- real(r8), intent(out) :: qs(vlen) ! Saturation specific humidity
-
- integer, intent(in), optional :: idx ! Scheme index
- integer :: i
-
- !$acc data present (t,p,es,qs)
-
- call wv_sat_svp_water_vect(t, es, vlen, idx)
- call wv_sat_svp_to_qsat_vect(es, p, qs, vlen)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i=1,vlen
- ! Ensures returned es is consistent with limiters on qs.
- es(i) = min(es(i), p(i))
- enddo
- !$acc end parallel
-
- !$acc end data
-end subroutine wv_sat_qsat_water_vect
-
-subroutine wv_sat_qsat_ice(t, p, es, qs, idx)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate SVP over ice at a given temperature, and then !
- ! calculate and return saturation specific humidity. !
- !------------------------------------------------------------------!
-
- ! Inputs
- real(r8), intent(in) :: t ! Temperature
- real(r8), intent(in) :: p ! Pressure
- ! Outputs
- real(r8), intent(out) :: es ! Saturation vapor pressure
- real(r8), intent(out) :: qs ! Saturation specific humidity
-
- integer, intent(in), optional :: idx ! Scheme index
-
- es = wv_sat_svp_ice(t, idx)
-
- qs = wv_sat_svp_to_qsat(es, p)
-
- ! Ensures returned es is consistent with limiters on qs.
- es = min(es, p)
-
-end subroutine wv_sat_qsat_ice
-
-subroutine wv_sat_qsat_ice_vect(t, p, es, qs, vlen, idx)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate SVP over ice at a given temperature, and then !
- ! calculate and return saturation specific humidity. !
- !------------------------------------------------------------------!
- ! Inputs
-
- integer, intent(in) :: vlen
- real(r8), intent(in) :: t(vlen) ! Temperature
- real(r8), intent(in) :: p(vlen) ! Pressure
- ! Outputs
- real(r8), intent(out) :: es(vlen) ! Saturation vapor pressure
- real(r8), intent(out) :: qs(vlen) ! Saturation specific humidity
-
- integer, intent(in), optional :: idx ! Scheme index
- integer :: i
-
- !$acc data present (t,p,es,qs)
-
- call wv_sat_svp_ice_vect(t, es, vlen, idx)
- call wv_sat_svp_to_qsat_vect(es, p, qs, vlen)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i=1,vlen
- ! Ensures returned es is consistent with limiters on qs.
- es(i) = min(es(i), p(i))
- enddo
- !$acc end parallel
-
- !$acc end data
-end subroutine wv_sat_qsat_ice_vect
-
-subroutine wv_sat_qsat_trans(t, p, es, qs, idx)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate SVP over ice at a given temperature, and then !
- ! calculate and return saturation specific humidity. !
- !------------------------------------------------------------------!
-
- ! Inputs
- real(r8), intent(in) :: t ! Temperature
- real(r8), intent(in) :: p ! Pressure
- ! Outputs
- real(r8), intent(out) :: es ! Saturation vapor pressure
- real(r8), intent(out) :: qs ! Saturation specific humidity
-
- integer, intent(in), optional :: idx ! Scheme index
-
- es = wv_sat_svp_trans(t, idx)
-
- qs = wv_sat_svp_to_qsat(es, p)
-
- ! Ensures returned es is consistent with limiters on qs.
- es = min(es, p)
-
-end subroutine wv_sat_qsat_trans
-
-!---------------------------------------------------------------------
-! SVP INTERFACE FUNCTIONS
-!---------------------------------------------------------------------
-
-function wv_sat_svp_water(t, idx) result(es)
- real(r8), intent(in) :: t
- integer, intent(in), optional :: idx
- real(r8) :: es
-
- integer :: use_idx
-
- if (present(idx)) then
- use_idx = idx
- else
- use_idx = default_idx
- end if
-
- select case (use_idx)
- case(GoffGratch_idx)
- es = GoffGratch_svp_water(t)
- case(MurphyKoop_idx)
- es = MurphyKoop_svp_water(t)
- case(Bolton_idx)
- es = Bolton_svp_water(t)
- end select
-
-end function wv_sat_svp_water
-
-subroutine wv_sat_svp_water_vect(t, es, vlen, idx)
- integer, intent(in) :: vlen
- real(r8), intent(in) :: t(vlen)
- integer, intent(in), optional :: idx
- real(r8), intent(out) :: es(vlen)
- integer :: i
- integer :: use_idx
-
- !$acc data present (t,es)
-
- if (present(idx)) then
- use_idx = idx
- else
- use_idx = default_idx
- end if
-
- select case (use_idx)
- case(GoffGratch_idx)
- call GoffGratch_svp_water_vect(t,es,vlen)
- case(MurphyKoop_idx)
- call MurphyKoop_svp_water_vect(t,es,vlen)
- case(Bolton_idx)
- call Bolton_svp_water_vect(t,es,vlen)
- end select
-
- !$acc end data
-end subroutine wv_sat_svp_water_vect
-
-function wv_sat_svp_ice(t, idx) result(es)
- real(r8), intent(in) :: t
- integer, intent(in), optional :: idx
- real(r8) :: es
-
- integer :: use_idx
-
- if (present(idx)) then
- use_idx = idx
- else
- use_idx = default_idx
- end if
-
- select case (use_idx)
- case(GoffGratch_idx)
- es = GoffGratch_svp_ice(t)
- case(MurphyKoop_idx)
- es = MurphyKoop_svp_ice(t)
- case(Bolton_idx)
- es = Bolton_svp_water(t)
- end select
-
-end function wv_sat_svp_ice
-
-subroutine wv_sat_svp_ice_vect(t, es, vlen, idx)
- integer, intent(in) :: vlen
- real(r8), intent(in) :: t(vlen)
- integer, intent(in), optional :: idx
- real(r8), intent(out) :: es(vlen)
- integer :: i
-
- integer :: use_idx
-
- !$acc data present (t,es)
-
- if (present(idx)) then
- use_idx = idx
- else
- use_idx = default_idx
- end if
-
- select case (use_idx)
- case(GoffGratch_idx)
- call GoffGratch_svp_ice_vect(t,es,vlen)
- case(MurphyKoop_idx)
- call MurphyKoop_svp_ice_vect(t,es,vlen)
- case(Bolton_idx)
- call Bolton_svp_water_vect(t,es,vlen)
- end select
-
- !$acc end data
-end subroutine wv_sat_svp_ice_vect
-
-function wv_sat_svp_trans(t, idx) result(es)
-
- real(r8), intent(in) :: t
- integer, intent(in), optional :: idx
- real(r8) :: es
-
- real(r8) :: esice ! Saturation vapor pressure over ice
- real(r8) :: weight ! Intermediate scratch variable for es transition
-
-!
-! Water
-!
- if (t >= (tmelt - ttrice)) then
- es = wv_sat_svp_water(t,idx)
- else
- es = 0.0_r8
- end if
-
-!
-! Ice
-!
- if (t < tmelt) then
-
- esice = wv_sat_svp_ice(t,idx)
-
- if ( (tmelt - t) > ttrice ) then
- weight = 1.0_r8
- else
- weight = (tmelt - t)/ttrice
- end if
-
- es = weight*esice + (1.0_r8 - weight)*es
- end if
-
-end function wv_sat_svp_trans
-
-subroutine wv_sat_svp_trans_vect(t, es, vlen, idx)
-
- integer, intent(in) :: vlen
- real(r8), intent(in) :: t(vlen)
- integer, intent(in), optional :: idx
- real(r8), intent(out) :: es(vlen)
-
- real(r8) :: esice(vlen) ! Saturation vapor pressure over ice
- real(r8) :: weight ! Intermediate scratch variable for es transition
- integer :: i
-
- !$acc data present (t,es) &
- !$acc create (esice)
-
-!
-! Water
-!
- call wv_sat_svp_water_vect(t,es,vlen,idx)
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i = 1, vlen
- if (t(i) < (tmelt - ttrice)) then
- es(i) = 0.0_r8
- end if
- end do
- !$acc end parallel
-!
-! Ice
-!
- call wv_sat_svp_ice_vect(t,esice,vlen,idx)
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i = 1, vlen
- if (t(i) < tmelt) then
- if ( (tmelt - t(i)) > ttrice ) then
- weight = 1.0_r8
- else
- weight = (tmelt - t(i))/ttrice
- end if
-
- es(i) = weight*esice(i) + (1.0_r8 - weight)*es(i)
- end if
- end do
- !$acc end parallel
-
- !$acc end data
-end subroutine wv_sat_svp_trans_vect
-
-!---------------------------------------------------------------------
-! SVP METHODS
-!---------------------------------------------------------------------
-
-! Goff & Gratch (1946)
-
-function GoffGratch_svp_water(t) result(es)
- real(r8), intent(in) :: t ! Temperature in Kelvin
- real(r8) :: es ! SVP in Pa
-
- ! uncertain below -70 C
- es = 10._r8**(-7.90298_r8*(tboil/t-1._r8)+ &
- 5.02808_r8*log10(tboil/t)- &
- 1.3816e-7_r8*(10._r8**(11.344_r8*(1._r8-t/tboil))-1._r8)+ &
- 8.1328e-3_r8*(10._r8**(-3.49149_r8*(tboil/t-1._r8))-1._r8)+ &
- log10(1013.246_r8))*100._r8
-
-end function GoffGratch_svp_water
-
-subroutine GoffGratch_svp_water_vect(t, es, vlen)
- integer, intent(in) :: vlen
- real(r8), intent(in) :: t(vlen) ! Temperature in Kelvin
- real(r8), intent(out) :: es(vlen) ! SVP in Pa
- real(r8) :: log_tboil
- integer :: i
-
- !$acc data present (t,es)
-
- ! Goff, J. A., and S. Gratch. “Low-Pressure Properties of Water from -160F
- ! to 212F.” Trans. Am. Soc. Heat. Vent. Eng. 52 (1946): 95–121.
- ! uncertain below -70 C
-
- log_tboil = log10(tboil)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i=1,vlen
- es(i) = 10._r8**(-7.90298_r8*(tboil/t(i)-1._r8)+ &
- 5.02808_r8*(log_tboil-log10(t(i)))- &
- 1.3816e-7_r8*(10._r8**(11.344_r8*(1._r8-t(i)/tboil))-1._r8)+ &
- 8.1328e-3_r8*(10._r8**(-3.49149_r8*(tboil/t(i)-1._r8))-1._r8)+ &
- log10(1013.246_r8))*100._r8
- enddo
- !$acc end parallel
-
- !$acc end data
-end subroutine GoffGratch_svp_water_vect
-
-function GoffGratch_svp_ice(t) result(es)
- real(r8), intent(in) :: t ! Temperature in Kelvin
- real(r8) :: es ! SVP in Pa
-
- ! good down to -100 C
- es = 10._r8**(-9.09718_r8*(h2otrip/t-1._r8)-3.56654_r8* &
- log10(h2otrip/t)+0.876793_r8*(1._r8-t/h2otrip)+ &
- log10(6.1071_r8))*100._r8
-
-end function GoffGratch_svp_ice
-
-subroutine GoffGratch_svp_ice_vect(t, es, vlen)
- integer, intent(in) :: vlen
- real(r8), intent(in) :: t(vlen) ! Temperature in Kelvin
- real(r8), intent(out) :: es(vlen) ! SVP in Pa
- real(r8), parameter :: log_param = log10(6.1071_r8)
- integer :: i
- ! good down to -100 C
-
- !$acc data present (t,es)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i=1,vlen
- es(i) = 10._r8**(-9.09718_r8*(h2otrip/t(i)-1._r8)-3.56654_r8* &
- log10(h2otrip/t(i))+0.876793_r8*(1._r8-t(i)/h2otrip)+ &
- log_param)*100._r8
- enddo
- !$acc end parallel
-
- !$acc end data
-end subroutine GoffGratch_svp_ice_vect
-
-! Murphy & Koop (2005)
-
-function MurphyKoop_svp_water(t) result(es)
- real(r8), intent(in) :: t ! Temperature in Kelvin
- real(r8) :: es ! SVP in Pa
-
- ! (good for 123 < T < 332 K)
- es = exp(54.842763_r8 - (6763.22_r8 / t) - (4.210_r8 * log(t)) + &
- (0.000367_r8 * t) + (tanh(0.0415_r8 * (t - 218.8_r8)) * &
- (53.878_r8 - (1331.22_r8 / t) - (9.44523_r8 * log(t)) + &
- 0.014025_r8 * t)))
-
-end function MurphyKoop_svp_water
-
-subroutine MurphyKoop_svp_water_vect(t, es, vlen)
- integer, intent(in) :: vlen
- real(r8), intent(in) :: t(vlen) ! Temperature in Kelvin
- real(r8), intent(out) :: es(vlen) ! SVP in Pa
-
- integer :: i
- ! Murphy, D. M., and T. Koop. “Review of the Vapour Pressure of Ice and
- ! Supercooled Water for Atmospheric Applications.” Q. J. R. Meteorol.
- ! Soc. 131, no. 608 (2005): 1539–65. 10.1256/qj.04.94
- ! (good for 123 < T < 332 K)
-
- !$acc data present (t,es)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i = 1, vlen
- es(i) = exp(54.842763_r8 - (6763.22_r8 / t(i)) - (4.210_r8 * log(t(i))) + &
- (0.000367_r8 * t(i)) + (tanh(0.0415_r8 * (t(i) - 218.8_r8)) * &
- (53.878_r8 - (1331.22_r8 / t(i)) - (9.44523_r8 * log(t(i))) + &
- 0.014025_r8 * t(i))))
- end do
- !$acc end parallel
-
- !$acc end data
-end subroutine MurphyKoop_svp_water_vect
-
-function MurphyKoop_svp_ice(t) result(es)
- real(r8), intent(in) :: t ! Temperature in Kelvin
- real(r8) :: es ! SVP in Pa
-
- ! (good down to 110 K)
- es = exp(9.550426_r8 - (5723.265_r8 / t) + (3.53068_r8 * log(t)) &
- - (0.00728332_r8 * t))
-
-end function MurphyKoop_svp_ice
-
-subroutine MurphyKoop_svp_ice_vect(t, es, vlen)
- integer, intent(in) :: vlen
- real(r8), intent(in) :: t(vlen) ! Temperature in Kelvin
- real(r8), intent(out) :: es(vlen) ! SVP in Pa
-
- integer :: i
- ! (good down to 110 K)
-
- !$acc data present (t,es)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i = 1, vlen
- es(i) = exp(9.550426_r8 - (5723.265_r8 / t(i)) + (3.53068_r8 * log(t(i))) &
- - (0.00728332_r8 * t(i)))
- end do
- !$acc end parallel
-
- !$acc end data
-end subroutine MurphyKoop_svp_ice_vect
-
-! Bolton (1980)
-! zm_conv deep convection scheme contained this SVP calculation.
-! It appears to be from D. Bolton, 1980, Monthly Weather Review.
-! Unlike the other schemes, no distinct ice formula is associated
-! with it. (However, a Bolton ice formula exists in CLUBB.)
-
-! The original formula used degrees C, but this function
-! takes Kelvin and internally converts.
-
-function Bolton_svp_water(t) result(es)
- real(r8),parameter :: c1 = 611.2_r8
- real(r8),parameter :: c2 = 17.67_r8
- real(r8),parameter :: c3 = 243.5_r8
-
- real(r8), intent(in) :: t ! Temperature in Kelvin
- real(r8) :: es ! SVP in Pa
-
- es = c1*exp( (c2*(t - tmelt))/((t - tmelt)+c3) )
-
-end function Bolton_svp_water
-
-subroutine Bolton_svp_water_vect(t, es,vlen)
- real(r8),parameter :: c1 = 611.2_r8
- real(r8),parameter :: c2 = 17.67_r8
- real(r8),parameter :: c3 = 243.5_r8
-
- integer, intent(in) :: vlen
- real(r8), intent(in) :: t(vlen) ! Temperature in Kelvin
- real(r8), intent(out) :: es(vlen) ! SVP in Pa
-
- integer :: i
-
- !$acc data present (t,es)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i = 1, vlen
- es(i) = c1*exp( (c2*(t(i) - tmelt))/((t(i) - tmelt)+c3) )
- end do
- !$acc end parallel
-
- !$acc end data
-end subroutine Bolton_svp_water_vect
-
-end module wv_sat_methods
diff --git a/src/physics/cam/wv_saturation.F90 b/src/physics/cam/wv_saturation.F90
deleted file mode 100644
index ac94482e20..0000000000
--- a/src/physics/cam/wv_saturation.F90
+++ /dev/null
@@ -1,1484 +0,0 @@
-module wv_saturation
-
-!--------------------------------------------------------------------!
-! Module Overview: !
-! !
-! This module provides an interface to wv_sat_methods, providing !
-! saturation vapor pressure and related calculations to CAM. !
-! !
-! The original wv_saturation codes were introduced by J. J. Hack, !
-! February 1990. The code has been extensively rewritten since then, !
-! including a total refactoring in Summer 2012. !
-! !
-!--------------------------------------------------------------------!
-! Methods: !
-! !
-! Pure water/ice saturation vapor pressures are calculated on the !
-! fly, with the specific method determined by a runtime option. !
-! Mixed phase SVP is interpolated from the internal table, estbl, !
-! which is created during initialization. !
-! !
-! The default method for calculating SVP is determined by a namelist !
-! option, and used whenever svp_water/ice or qsat are called. !
-! !
-!--------------------------------------------------------------------!
-
-use shr_kind_mod, only: r8 => shr_kind_r8
-use physconst, only: epsilo, &
- latvap, &
- latice, &
- rh2o, &
- cpair, &
- tmelt, &
- h2otrip
-
-use wv_sat_methods, only: &
- svp_to_qsat => wv_sat_svp_to_qsat, &
- svp_to_qsat_vect => wv_sat_svp_to_qsat_vect
-
-implicit none
-private
-save
-
-! Public interfaces
-! Namelist, initialization, finalization
-public wv_sat_readnl
-public wv_sat_init
-public wv_sat_final
-
-! Saturation vapor pressure calculations
-public svp_water, svp_water_vect
-public svp_ice, svp_ice_vect
-
-! Mixed phase (water + ice) saturation vapor pressure table lookup
-public estblf
-
-public svp_to_qsat
-
-! Subroutines that return both SVP and humidity
-! Optional arguments do temperature derivatives
-interface qsat
- module procedure qsat_line
- module procedure qsat_vect
- module procedure qsat_2D
-end interface
-public qsat ! Mixed phase
-interface qsat_water
- module procedure qsat_water_line
- module procedure qsat_water_vect
- module procedure qsat_water_2D
-end interface
-public qsat_water ! SVP over water only
-interface qsat_ice
- module procedure qsat_ice_line
- module procedure qsat_ice_vect
- module procedure qsat_ice_2D
-end interface
-public qsat_ice ! SVP over ice only
-
-! Wet bulb temperature solver
-public :: findsp_vc, findsp
-
-! Data
-
-! This value is slightly high, but it seems to be the value for the
-! steam point of water originally (and most frequently) used in the
-! Goff & Gratch scheme.
-real(r8), parameter :: tboil = 373.16_r8
-
-! Table of saturation vapor pressure values (estbl) from tmin to
-! tmax+1 Kelvin, in one degree increments. ttrice defines the
-! transition region, estbl contains a combination of ice & water
-! values.
-! Make these public parameters in case another module wants to see the
-! extent of the table.
- real(r8), public, parameter :: tmin = 127.16_r8
- real(r8), public, parameter :: tmax = 375.16_r8
-
- real(r8), parameter :: ttrice = 20.00_r8 ! transition range from es over H2O to es over ice
-
- integer :: plenest ! length of estbl
- real(r8), allocatable :: estbl(:) ! table values of saturation vapor pressure
-
- real(r8) :: omeps ! 1.0_r8 - epsilo
-
- real(r8) :: c3 ! parameter used by findsp
-
- ! Set coefficients for polynomial approximation of difference
- ! between saturation vapor press over water and saturation pressure
- ! over ice for -ttrice < t < 0 (degrees C). NOTE: polynomial is
- ! valid in the range -40 < t < 0 (degrees C).
- real(r8) :: pcf(5) = (/ &
- 5.04469588506e-01_r8, &
- -5.47288442819e+00_r8, &
- -3.67471858735e-01_r8, &
- -8.95963532403e-03_r8, &
- -7.78053686625e-05_r8 /)
-
-! --- Degree 6 approximation ---
-! real(r8) :: pcf(6) = (/ &
-! 7.63285250063e-02, &
-! 5.86048427932e+00, &
-! 4.38660831780e-01, &
-! 1.37898276415e-02, &
-! 2.14444472424e-04, &
-! 1.36639103771e-06 /)
-
- integer, parameter :: VLENS = 128 ! vector length for a GPU kernel
-
- !$acc declare create (plenest,estbl,omeps,c3,pcf)
-
-contains
-
-!---------------------------------------------------------------------
-! ADMINISTRATIVE FUNCTIONS
-!---------------------------------------------------------------------
-
-subroutine wv_sat_readnl(nlfile)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Get runtime options for wv_saturation. !
- !------------------------------------------------------------------!
-
- use wv_sat_methods, only: wv_sat_get_scheme_idx, &
- wv_sat_valid_idx, &
- wv_sat_set_default
-
- use spmd_utils, only: masterproc
- use namelist_utils, only: find_group_name
- use units, only: getunit, freeunit
- use mpishorthand
- use cam_abortutils, only: endrun
-
- character(len=*), intent(in) :: nlfile ! filepath for file containing namelist input
-
- ! Local variables
- integer :: unitn, ierr
-
- character(len=32) :: wv_sat_scheme = "GoffGratch"
-
- character(len=*), parameter :: subname = 'wv_sat_readnl'
-
- namelist /wv_sat_nl/ wv_sat_scheme
- !-----------------------------------------------------------------------------
-
- if (masterproc) then
- unitn = getunit()
- open( unitn, file=trim(nlfile), status='old' )
- call find_group_name(unitn, 'wv_sat_nl', status=ierr)
- if (ierr == 0) then
- read(unitn, wv_sat_nl, iostat=ierr)
- if (ierr /= 0) then
- call endrun(subname // ':: ERROR reading namelist')
- return
- end if
- end if
- close(unitn)
- call freeunit(unitn)
-
- end if
-
-#ifdef SPMD
- call mpibcast(wv_sat_scheme, len(wv_sat_scheme) , mpichar, 0, mpicom)
-#endif
-
- if (.not. wv_sat_set_default(wv_sat_scheme)) then
- call endrun('wv_sat_readnl :: Invalid wv_sat_scheme.')
- return
- end if
-
-end subroutine wv_sat_readnl
-
-subroutine wv_sat_init
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Initialize module (e.g. setting parameters, initializing the !
- ! SVP lookup table). !
- !------------------------------------------------------------------!
-
- use wv_sat_methods, only: wv_sat_methods_init, &
- wv_sat_get_scheme_idx, &
- wv_sat_valid_idx
- use spmd_utils, only: masterproc
- use cam_logfile, only: iulog
- use cam_abortutils, only: endrun
- use shr_assert_mod, only: shr_assert_in_domain
- use error_messages, only: handle_errmsg
-
- integer :: status
-
- ! For wv_sat_methods error reporting.
- character(len=256) :: errstring
-
- ! For generating internal SVP table.
- real(r8) :: t ! Temperature
- integer :: i ! Increment counter
-
- ! Precalculated because so frequently used.
- omeps = 1.0_r8 - epsilo
-
- ! Transition range method is only valid for transition temperatures at:
- ! -40 deg C < T < 0 deg C
- call shr_assert_in_domain(ttrice, ge=0._r8, le=40._r8, varname="ttrice",&
- msg="wv_sat_init: Invalid transition temperature range.")
-
-! This parameter uses a hardcoded 287.04_r8?
- c3 = 287.04_r8*(7.5_r8*log(10._r8))/cpair
-
-! Init "methods" module containing actual SVP formulae.
-
- call wv_sat_methods_init(r8, tmelt, h2otrip, tboil, ttrice, &
- epsilo, errstring)
-
- call handle_errmsg(errstring, subname="wv_sat_methods_init")
-
- ! Add two to make the table slightly too big, just in case.
- plenest = ceiling(tmax-tmin) + 2
-
- ! Allocate SVP table.
- allocate(estbl(plenest), stat=status)
- if (status /= 0) then
- call endrun('wv_sat_init :: ERROR allocating saturation vapor pressure table')
- return
- end if
-
- do i = 1, plenest
- estbl(i) = svp_trans(tmin + real(i-1,r8))
- end do
-
- !$acc update device (plenest,estbl,omeps,c3,pcf)
-
- if (masterproc) then
- write(iulog,*)' *** SATURATION VAPOR PRESSURE TABLE COMPLETED ***'
- end if
-
-end subroutine wv_sat_init
-
-subroutine wv_sat_final
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Deallocate global variables in module. !
- !------------------------------------------------------------------!
- use cam_abortutils, only: endrun
-
- integer :: status
-
- if (allocated(estbl)) then
-
- deallocate(estbl, stat=status)
-
- if (status /= 0) then
- call endrun('wv_sat_final :: ERROR deallocating table')
- return
- end if
-
- end if
-
-end subroutine wv_sat_final
-
-!---------------------------------------------------------------------
-! DEFAULT SVP FUNCTIONS
-!---------------------------------------------------------------------
-
-! Compute saturation vapor pressure over water
-function svp_water(t) result(es)
-
- use wv_sat_methods, only: &
- wv_sat_svp_water
-
- real(r8), intent(in) :: t ! Temperature (K)
- real(r8) :: es ! SVP (Pa)
-
- es = wv_sat_svp_water(t)
-
-end function svp_water
-
-! Compute saturation vapor pressure over water
-subroutine svp_water_vect(t, es, vlen)
-
- use wv_sat_methods, only: &
- wv_sat_svp_water_vect
-
- integer, intent(in) :: vlen
- real(r8), intent(in) :: t(vlen) ! Temperature (K)
- real(r8), intent(out) :: es(vlen) ! SVP (Pa)
-
- !$acc data copyin (t) copyout (es)
-
- call wv_sat_svp_water_vect(t, es, vlen)
-
- !$acc end data
-end subroutine svp_water_vect
-
-! Compute saturation vapor pressure over ice
-function svp_ice(t) result(es)
-
- use wv_sat_methods, only: &
- wv_sat_svp_ice
-
- real(r8), intent(in) :: t ! Temperature (K)
- real(r8) :: es ! SVP (Pa)
-
- es = wv_sat_svp_ice(t)
-
-end function svp_ice
-
-! Compute saturation vapor pressure over ice
-subroutine svp_ice_vect(t, es, vlen)
-
- use wv_sat_methods, only: &
- wv_sat_svp_ice_vect
-
- integer, intent(in) :: vlen
- real(r8), intent(in) :: t(vlen) ! Temperature (K)
- real(r8), intent(out) :: es(vlen) ! SVP (Pa)
-
- !$acc data copyin(t) copyout(es)
-
- call wv_sat_svp_ice_vect(t, es, vlen)
-
- !$acc end data
-end subroutine svp_ice_vect
-
-! Compute saturation vapor pressure with an ice-water transition
-function svp_trans(t) result(es)
-
- use wv_sat_methods, only: &
- wv_sat_svp_trans
-
- real(r8), intent(in) :: t ! Temperature (K)
- real(r8) :: es ! SVP (Pa)
-
- es = wv_sat_svp_trans(t)
-
-end function svp_trans
-
-! Compute saturation vapor pressure with an ice-water transition
-subroutine svp_trans_vect(t, es, vlen)
-
- use wv_sat_methods, only: &
- wv_sat_svp_trans_vect
-
- integer, intent(in) :: vlen
- real(r8), intent(in) :: t(vlen) ! Temperature (K)
- real(r8), intent(out) :: es(vlen) ! SVP (Pa)
-
- !$acc data copyin(t) copyout(es)
-
- call wv_sat_svp_trans_vect(t, es, vlen)
-
- !$acc end data
-end subroutine svp_trans_vect
-
-!---------------------------------------------------------------------
-! UTILITIES
-!---------------------------------------------------------------------
-
-! Does linear interpolation from nearest values found
-! in the table (estbl).
-elemental function estblf(t) result(es)
-
- real(r8), intent(in) :: t ! Temperature
- real(r8) :: es ! SVP (Pa)
-
- integer :: i ! Index for t in the table
- real(r8) :: t_tmp ! intermediate temperature for es look-up
-
- real(r8) :: weight ! Weight for interpolation
-
- t_tmp = max(min(t,tmax)-tmin, 0._r8) ! Number of table entries above tmin
- i = int(t_tmp) + 1 ! Corresponding index.
- weight = t_tmp - aint(t_tmp, r8) ! Fractional part of t_tmp (for interpolation).
- es = (1._r8 - weight)*estbl(i) + weight*estbl(i+1)
-
-end function estblf
-
-! Does linear interpolation from nearest values found
-! in the table (estbl).
-subroutine estblf_vect(t, es, vlen)
-
- integer, intent(in) :: vlen
- real(r8), dimension(vlen), intent(in) :: t ! Temperature
- real(r8), dimension(vlen), intent(out) :: es ! SVP (Pa)
-
- integer :: i ! Index for t in the table
- integer :: j
- real(r8) :: t_tmp ! intermediate temperature for es look-up
-
- real(r8) :: weight ! Weight for interpolation
-
- !$acc data present (t,es)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector private(t_tmp,weight,i)
- do j = 1, vlen
- t_tmp = max(min(t(j),tmax)-tmin, 0._r8) ! Number of table entries above tmin
- i = int(t_tmp) + 1 ! Corresponding index.
- weight = t_tmp - aint(t_tmp, r8) ! Fractional part of t_tmp (for interpolation).
- es(j) = (1._r8 - weight)*estbl(i) + weight*estbl(i+1)
- end do
- !$acc end parallel
-
- !$acc end data
-end subroutine estblf_vect
-
-! Get enthalpy based only on temperature
-! and specific humidity.
-elemental function tq_enthalpy(t, q, hltalt) result(enthalpy)
-
- real(r8), intent(in) :: t ! Temperature
- real(r8), intent(in) :: q ! Specific humidity
- real(r8), intent(in) :: hltalt ! Modified hlat for T derivatives
-
- real(r8) :: enthalpy
-
- enthalpy = cpair * t + hltalt * q
-
-end function tq_enthalpy
-
-! Get enthalpy based only on temperature
-! and specific humidity.
-subroutine tq_enthalpy_vect(t, q, hltalt, enthalpy, vlen)
-
- integer, intent(in) :: vlen
- real(r8), dimension(vlen), intent(in) :: t ! Temperature
- real(r8), dimension(vlen), intent(in) :: q ! Specific humidity
- real(r8), dimension(vlen), intent(in) :: hltalt ! Modified hlat for T derivatives
-
- real(r8), dimension(vlen), intent(out) :: enthalpy
-
- integer :: i
-
- !$acc data present(t,q,hltalt,enthalpy)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i = 1, vlen
- enthalpy(i) = cpair * t(i) + hltalt(i) * q(i)
- end do
- !$acc end parallel
-
- !$acc end data
-end subroutine tq_enthalpy_vect
-
-!---------------------------------------------------------------------
-! LATENT HEAT OF VAPORIZATION CORRECTIONS
-!---------------------------------------------------------------------
-
-elemental subroutine no_ip_hltalt(t, hltalt)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate latent heat of vaporization of pure liquid water at !
- ! a given temperature. !
- !------------------------------------------------------------------!
-
- ! Inputs
- real(r8), intent(in) :: t ! Temperature
- ! Outputs
- real(r8), intent(out) :: hltalt ! Appropriately modified hlat
-
- hltalt = latvap
-
- ! Account for change of latvap with t above freezing where
- ! constant slope is given by -2369 j/(kg c) = cpv - cw
- if (t >= tmelt) then
- hltalt = hltalt - 2369.0_r8*(t-tmelt)
- end if
-
-end subroutine no_ip_hltalt
-
-subroutine no_ip_hltalt_vect(t, hltalt, vlen)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate latent heat of vaporization of pure liquid water at !
- ! a given temperature. !
- !------------------------------------------------------------------!
-
- ! Inputs
- integer, intent(in) :: vlen
- real(r8), dimension(vlen), intent(in) :: t ! Temperature
- ! Outputs
- real(r8), dimension(vlen), intent(out) :: hltalt ! Appropriately modified hlat
-
- integer :: i
-
- !$acc data present(t,hltalt)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i = 1, vlen
- hltalt(i) = latvap
- ! Account for change of latvap with t above freezing where
- ! constant slope is given by -2369 j/(kg c) = cpv - cw
- if (t(i) >= tmelt) then
- hltalt(i) = hltalt(i) - 2369.0_r8*(t(i)-tmelt)
- end if
- end do
- !$acc end parallel
-
- !$acc end data
-end subroutine no_ip_hltalt_vect
-
-elemental subroutine calc_hltalt(t, hltalt, tterm)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate latent heat of vaporization of water at a given !
- ! temperature, taking into account the ice phase if temperature !
- ! is below freezing. !
- ! Optional argument also calculates a term used to calculate !
- ! d(es)/dT within the water-ice transition range. !
- !------------------------------------------------------------------!
-
- ! Inputs
- real(r8), intent(in) :: t ! Temperature
- ! Outputs
- real(r8), intent(out) :: hltalt ! Appropriately modified hlat
- ! Term to account for d(es)/dT in transition region.
- real(r8), intent(out), optional :: tterm
-
- ! Local variables
- real(r8) :: tc ! Temperature in degrees C
- real(r8) :: weight ! Weight for es transition from water to ice
- ! Loop iterator
- integer :: i
-
- if (present(tterm)) tterm = 0.0_r8
-
- call no_ip_hltalt(t,hltalt)
- if (t < tmelt) then
- ! Weighting of hlat accounts for transition from water to ice.
- tc = t - tmelt
-
- if (tc >= -ttrice) then
- weight = -tc/ttrice
-
- ! polynomial expression approximates difference between es
- ! over water and es over ice from 0 to -ttrice (C) (max of
- ! ttrice is 40): required for accurate estimate of es
- ! derivative in transition range from ice to water
- if (present(tterm)) then
- do i = size(pcf), 1, -1
- tterm = pcf(i) + tc*tterm
- end do
- tterm = tterm/ttrice
- end if
-
- else
- weight = 1.0_r8
- end if
-
- hltalt = hltalt + weight*latice
-
- end if
-
-end subroutine calc_hltalt
-
-subroutine calc_hltalt_vect(t, hltalt, vlen, tterm)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate latent heat of vaporization of water at a given !
- ! temperature, taking into account the ice phase if temperature !
- ! is below freezing. !
- ! Optional argument also calculates a term used to calculate !
- ! d(es)/dT within the water-ice transition range. !
- !------------------------------------------------------------------!
-
- ! Inputs
- integer, intent(in) :: vlen
- real(r8), dimension(vlen), intent(in) :: t ! Temperature
- ! Outputs
- real(r8), dimension(vlen), intent(out) :: hltalt ! Appropriately modified hlat
- ! Term to account for d(es)/dT in transition region.
- real(r8), dimension(vlen), intent(out), optional :: tterm
-
- ! Local variables
- real(r8) :: tc ! Temperature in degrees C
- real(r8) :: weight ! Weight for es transition from water to ice
- logical :: present_tterm
- ! Loop iterator
- integer :: i, j, size_pcf
-
- present_tterm = present(tterm)
- size_pcf = size(pcf)
-
- !$acc data present(t,hltalt,tterm)
-
- if (present_tterm) then
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i = 1, vlen
- tterm(i) = 0.0_r8
- end do
- !$acc end parallel
- end if
-
- call no_ip_hltalt_vect(t,hltalt,vlen)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector private(tc,weight)
- do j = 1, vlen
- if (t(j) < tmelt) then
- ! Weighting of hlat accounts for transition from water to ice.
- tc = t(j) - tmelt
-
- if (tc >= -ttrice) then
- weight = -tc/ttrice
-
- ! polynomial expression approximates difference between es
- ! over water and es over ice from 0 to -ttrice (C) (max of
- ! ttrice is 40): required for accurate estimate of es
- ! derivative in transition range from ice to water
- if (present_tterm) then
- !$acc loop seq
- do i = size_pcf, 1, -1
- tterm(j) = pcf(i) + tc*tterm(j)
- end do
- tterm(j) = tterm(j)/ttrice
- end if
-
- else
- weight = 1.0_r8
- end if
-
- hltalt(j) = hltalt(j) + weight*latice
-
- end if
- end do
- !$acc end parallel
-
- !$acc end data
-end subroutine calc_hltalt_vect
-
-!---------------------------------------------------------------------
-! OPTIONAL OUTPUTS
-!---------------------------------------------------------------------
-
-! Temperature derivative outputs, for qsat_*
-subroutine deriv_outputs_line(t, p, es, qs, hltalt, tterm, &
- gam, dqsdt)
-
- ! Inputs
- real(r8), intent(in) :: t ! Temperature
- real(r8), intent(in) :: p ! Pressure
- real(r8), intent(in) :: es ! Saturation vapor pressure
- real(r8), intent(in) :: qs ! Saturation specific humidity
- real(r8), intent(in) :: hltalt ! Modified latent heat
- real(r8), intent(in) :: tterm ! Extra term for d(es)/dT in
- ! transition region.
-
- ! Outputs
- real(r8), intent(out), optional :: gam ! (hltalt/cpair)*(d(qs)/dt)
- real(r8), intent(out), optional :: dqsdt ! (d(qs)/dt)
-
- ! Local variables
- real(r8) :: desdt ! d(es)/dt
- real(r8) :: dqsdt_loc ! local copy of dqsdt
-
- if (qs == 1.0_r8) then
- dqsdt_loc = 0._r8
- else
- desdt = hltalt*es/(rh2o*t*t) + tterm
- dqsdt_loc = qs*p*desdt/(es*(p-omeps*es))
- end if
-
- if (present(dqsdt)) dqsdt = dqsdt_loc
- if (present(gam)) gam = dqsdt_loc * (hltalt/cpair)
-
-end subroutine deriv_outputs_line
-
-! Temperature derivative outputs, for qsat_*
-subroutine deriv_outputs_vect(t, p, es, qs, hltalt, tterm, vlen, &
- gam, dqsdt)
-
- ! Inputs
- integer, intent(in) :: vlen
- real(r8), dimension(vlen), intent(in) :: t ! Temperature
- real(r8), dimension(vlen), intent(in) :: p ! Pressure
- real(r8), dimension(vlen), intent(in) :: es ! Saturation vapor pressure
- real(r8), dimension(vlen), intent(in) :: qs ! Saturation specific humidity
- real(r8), dimension(vlen), intent(in) :: hltalt ! Modified latent heat
- real(r8), dimension(vlen), intent(in) :: tterm ! Extra term for d(es)/dT in
- ! transition region.
-
- ! Outputs
- real(r8), dimension(vlen), intent(out), optional :: gam ! (hltalt/cpair)*(d(qs)/dt)
- real(r8), dimension(vlen), intent(out), optional :: dqsdt ! (d(qs)/dt)
-
- ! Local variables
- real(r8) :: desdt ! d(es)/dt
- real(r8) :: dqsdt_loc ! local copy of dqsdt
- logical :: present_dqsdt, present_gam
- integer :: i
-
- present_dqsdt = present(dqsdt)
- present_gam = present(gam)
-
- !$acc data present(t,p,es,qs,hltalt,tterm,gam,dqsdt)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector private(dqsdt_loc,desdt)
- do i = 1, vlen
- if (qs(i) == 1.0_r8) then
- dqsdt_loc = 0._r8
- else
- desdt = hltalt(i)*es(i)/(rh2o*t(i)*t(i)) + tterm(i)
- dqsdt_loc = qs(i)*p(i)*desdt/(es(i)*(p(i)-omeps*es(i)))
- end if
-
- if (present_dqsdt) dqsdt(i) = dqsdt_loc
- if (present_gam) gam(i) = dqsdt_loc * (hltalt(i)/cpair)
- end do
- !$acc end parallel
-
- !$acc end data
-end subroutine deriv_outputs_vect
-
-!---------------------------------------------------------------------
-! QSAT (SPECIFIC HUMIDITY) PROCEDURES
-!---------------------------------------------------------------------
-
-subroutine qsat_line(t, p, es, qs, gam, dqsdt, enthalpy)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Look up and return saturation vapor pressure from precomputed !
- ! table, then calculate and return saturation specific humidity. !
- ! Optionally return various temperature derivatives or enthalpy !
- ! at saturation. !
- !------------------------------------------------------------------!
-
- ! Inputs
- real(r8), intent(in) :: t ! Temperature
- real(r8), intent(in) :: p ! Pressure
- ! Outputs
- real(r8), intent(out) :: es ! Saturation vapor pressure
- real(r8), intent(out) :: qs ! Saturation specific humidity
-
- real(r8), intent(out), optional :: gam ! (l/cpair)*(d(qs)/dt)
- real(r8), intent(out), optional :: dqsdt ! (d(qs)/dt)
- real(r8), intent(out), optional :: enthalpy ! cpair*t + hltalt*q
-
- ! Local variables
- real(r8) :: hltalt ! Modified latent heat for T derivatives
- real(r8) :: tterm ! Account for d(es)/dT in transition region
-
- es = estblf(t)
-
- qs = svp_to_qsat(es, p)
-
- ! Ensures returned es is consistent with limiters on qs.
- es = min(es, p)
-
- ! Calculate optional arguments.
- if (present(gam) .or. present(dqsdt) .or. present(enthalpy)) then
-
- ! "generalized" analytic expression for t derivative of es
- ! accurate to within 1 percent for 173.16 < t < 373.16
- call calc_hltalt(t, hltalt, tterm)
-
- if (present(enthalpy)) enthalpy = tq_enthalpy(t, qs, hltalt)
-
- call deriv_outputs_line(t, p, es, qs, hltalt, tterm, &
- gam=gam, dqsdt=dqsdt)
-
- end if
-
-end subroutine qsat_line
-
-subroutine qsat_vect(t, p, es, qs, vlen, gam, dqsdt, enthalpy)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Look up and return saturation vapor pressure from precomputed !
- ! table, then calculate and return saturation specific humidity. !
- ! Optionally return various temperature derivatives or enthalpy !
- ! at saturation. !
- !------------------------------------------------------------------!
-
- ! Inputs
- integer, intent(in) :: vlen
- real(r8), dimension(vlen), intent(in) :: t ! Temperature
- real(r8), dimension(vlen), intent(in) :: p ! Pressure
- ! Outputs
- real(r8), dimension(vlen), intent(out) :: es ! Saturation vapor pressure
- real(r8), dimension(vlen), intent(out) :: qs ! Saturation specific humidity
-
- real(r8), dimension(vlen), intent(out), optional :: gam ! (l/cpair)*(d(qs)/dt)
- real(r8), dimension(vlen), intent(out), optional :: dqsdt ! (d(qs)/dt)
- real(r8), dimension(vlen), intent(out), optional :: enthalpy ! cpair*t + hltalt*q
-
- ! Local variables
- real(r8), dimension(vlen) :: hltalt ! Modified latent heat for T derivatives
- real(r8), dimension(vlen) :: tterm ! Account for d(es)/dT in transition region
- integer :: i
- logical :: present_gam, present_dqsdt, present_enthalpy
-
- present_gam = present(gam)
- present_dqsdt = present(dqsdt)
- present_enthalpy = present(enthalpy)
-
- !$acc data copyin (t,p) &
- !$acc copyout (es,qs,gam,dqsdt,enthalpy) &
- !$acc create (hltalt,tterm)
-
- call estblf_vect(t, es, vlen)
-
- call svp_to_qsat_vect(es, p, qs, vlen)
-
- ! Ensures returned es is consistent with limiters on qs.
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i = 1, vlen
- es(i) = min(es(i), p(i))
- end do
- !$acc end parallel
-
- ! Calculate optional arguments.
- if (present_gam .or. present_dqsdt .or. present_enthalpy) then
-
- ! "generalized" analytic expression for t derivative of es
- ! accurate to within 1 percent for 173.16 < t < 373.16
- call calc_hltalt_vect(t, hltalt, vlen, tterm)
-
- if (present_enthalpy) call tq_enthalpy_vect(t, qs, hltalt, enthalpy, vlen)
-
- call deriv_outputs_vect(t, p, es, qs, hltalt, tterm, vlen, &
- gam=gam, dqsdt=dqsdt)
-
- end if
-
- !$acc end data
-end subroutine qsat_vect
-
-subroutine qsat_2D(t, p, es, qs, dim1, dim2, gam, dqsdt, enthalpy)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Look up and return saturation vapor pressure from precomputed !
- ! table, then calculate and return saturation specific humidity. !
- ! Optionally return various temperature derivatives or enthalpy !
- ! at saturation. !
- !------------------------------------------------------------------!
-
- ! Inputs
- integer, intent(in) :: dim1, dim2
- real(r8), dimension(dim1,dim2), intent(in) :: t ! Temperature
- real(r8), dimension(dim1,dim2), intent(in) :: p ! Pressure
- ! Outputs
- real(r8), dimension(dim1,dim2), intent(out) :: es ! Saturation vapor pressure
- real(r8), dimension(dim1,dim2), intent(out) :: qs ! Saturation specific humidity
-
- real(r8), dimension(dim1,dim2), intent(out), optional :: gam ! (l/cpair)*(d(qs)/dt)
- real(r8), dimension(dim1,dim2), intent(out), optional :: dqsdt ! (d(qs)/dt)
- real(r8), dimension(dim1,dim2), intent(out), optional :: enthalpy ! cpair*t + hltalt*q
-
- ! Local variables
- real(r8), dimension(dim1,dim2) :: hltalt ! Modified latent heat for T derivatives
- real(r8), dimension(dim1,dim2) :: tterm ! Account for d(es)/dT in transition region
- integer :: i, k, vlen
- logical :: present_gam, present_dqsdt, present_enthalpy
-
- vlen = dim1 * dim2
- present_gam = present(gam)
- present_dqsdt = present(dqsdt)
- present_enthalpy = present(enthalpy)
-
- !$acc data copyin (t,p) &
- !$acc copyout (es,qs,gam,dqsdt,enthalpy) &
- !$acc create (hltalt,tterm)
-
- call estblf_vect(t, es, vlen)
-
- call svp_to_qsat_vect(es, p, qs, vlen)
-
- ! Ensures returned es is consistent with limiters on qs.
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector collapse(2)
- do k = 1, dim2
- do i = 1, dim1
- es(i,k) = min(es(i,k), p(i,k))
- end do
- end do
- !$acc end parallel
-
- ! Calculate optional arguments.
- if (present_gam .or. present_dqsdt .or. present_enthalpy) then
-
- ! "generalized" analytic expression for t derivative of es
- ! accurate to within 1 percent for 173.16 < t < 373.16
- call calc_hltalt_vect(t, hltalt, vlen, tterm)
-
- if (present_enthalpy) call tq_enthalpy_vect(t, qs, hltalt, enthalpy, vlen)
-
- call deriv_outputs_vect(t, p, es, qs, hltalt, tterm, vlen, &
- gam=gam, dqsdt=dqsdt)
-
- end if
-
- !$acc end data
-end subroutine qsat_2D
-
-subroutine qsat_water_line(t, p, es, qs, gam, dqsdt, enthalpy)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate SVP over water at a given temperature, and then !
- ! calculate and return saturation specific humidity. !
- ! Optionally return various temperature derivatives or enthalpy !
- ! at saturation. !
- !------------------------------------------------------------------!
-
- use wv_sat_methods, only: wv_sat_qsat_water
-
- ! Inputs
- real(r8), intent(in) :: t ! Temperature
- real(r8), intent(in) :: p ! Pressure
- ! Outputs
- real(r8), intent(out) :: es ! Saturation vapor pressure
- real(r8), intent(out) :: qs ! Saturation specific humidity
-
- real(r8), intent(out), optional :: gam ! (l/cpair)*(d(qs)/dt)
- real(r8), intent(out), optional :: dqsdt ! (d(qs)/dt)
- real(r8), intent(out), optional :: enthalpy ! cpair*t + hltalt*q
-
- ! Local variables
- real(r8) :: hltalt ! Modified latent heat for T derivatives
-
- call wv_sat_qsat_water(t, p, es, qs)
-
- if (present(gam) .or. present(dqsdt) .or. present(enthalpy)) then
-
- ! "generalized" analytic expression for t derivative of es
- ! accurate to within 1 percent for 173.16 < t < 373.16
- call no_ip_hltalt(t, hltalt)
-
- if (present(enthalpy)) enthalpy = tq_enthalpy(t, qs, hltalt)
-
- ! For pure water/ice transition term is 0.
- call deriv_outputs_line(t, p, es, qs, hltalt, 0._r8, &
- gam=gam, dqsdt=dqsdt)
-
- end if
-
-end subroutine qsat_water_line
-
-subroutine qsat_water_vect(t, p, es, qs, vlen, gam, dqsdt, enthalpy)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate SVP over water at a given temperature, and then !
- ! calculate and return saturation specific humidity. !
- ! Optionally return various temperature derivatives or enthalpy !
- ! at saturation. !
- !------------------------------------------------------------------!
-
- use wv_sat_methods, only: wv_sat_qsat_water_vect
-
- ! Inputs
- integer, intent(in) :: vlen
- real(r8), dimension(vlen), intent(in) :: t ! Temperature
- real(r8), dimension(vlen), intent(in) :: p ! Pressure
- ! Outputs
- real(r8), dimension(vlen), intent(out) :: es ! Saturation vapor pressure
- real(r8), dimension(vlen), intent(out) :: qs ! Saturation specific humidity
-
- real(r8), dimension(vlen), intent(out), optional :: gam ! (l/cpair)*(d(qs)/dt)
- real(r8), dimension(vlen), intent(out), optional :: dqsdt ! (d(qs)/dt)
- real(r8), dimension(vlen), intent(out), optional :: enthalpy ! cpair*t + hltalt*q
-
- ! Local variables
- real(r8), dimension(vlen) :: hltalt ! Modified latent heat for T derivatives
- real(r8), dimension(vlen) :: tterm
- integer :: i
- logical :: present_gam, present_dqsdt, present_enthalpy
-
- present_gam = present(gam)
- present_dqsdt = present(dqsdt)
- present_enthalpy = present(enthalpy)
-
- !$acc data copyin (t,p) &
- !$acc copyout (es,qs,gam,dqsdt,enthalpy) &
- !$acc create (tterm,hltalt)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i = 1, vlen
- tterm(i) = 0._r8
- end do
- !$acc end parallel
-
- call wv_sat_qsat_water_vect(t, p, es, qs, vlen)
-
- if (present_gam .or. present_dqsdt .or. present_enthalpy) then
-
- ! "generalized" analytic expression for t derivative of es
- ! accurate to within 1 percent for 173.16 < t < 373.16
- call no_ip_hltalt_vect(t, hltalt, vlen)
-
- if (present_enthalpy) call tq_enthalpy_vect(t, qs, hltalt, enthalpy, vlen)
-
- ! For pure water/ice transition term is 0.
- call deriv_outputs_vect(t, p, es, qs, hltalt, tterm, vlen, &
- gam=gam, dqsdt=dqsdt)
-
- end if
-
- !$acc end data
-end subroutine qsat_water_vect
-
-subroutine qsat_water_2D(t, p, es, qs, dim1, dim2, gam, dqsdt, enthalpy)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate SVP over water at a given temperature, and then !
- ! calculate and return saturation specific humidity. !
- ! Optionally return various temperature derivatives or enthalpy !
- ! at saturation. !
- !------------------------------------------------------------------!
-
- use wv_sat_methods, only: wv_sat_qsat_water_vect
-
- ! Inputs
- integer, intent(in) :: dim1, dim2
- real(r8), dimension(dim1,dim2), intent(in) :: t ! Temperature
- real(r8), dimension(dim1,dim2), intent(in) :: p ! Pressure
- ! Outputs
- real(r8), dimension(dim1,dim2), intent(out) :: es ! Saturation vapor pressure
- real(r8), dimension(dim1,dim2), intent(out) :: qs ! Saturation specific humidity
-
- real(r8), dimension(dim1,dim2), intent(out), optional :: gam ! (l/cpair)*(d(qs)/dt)
- real(r8), dimension(dim1,dim2), intent(out), optional :: dqsdt ! (d(qs)/dt)
- real(r8), dimension(dim1,dim2), intent(out), optional :: enthalpy ! cpair*t + hltalt*q
-
- ! Local variables
- real(r8), dimension(dim1,dim2) :: hltalt ! Modified latent heat for T derivatives
- real(r8), dimension(dim1,dim2) :: tterm
- integer :: i, k, vlen
- logical :: present_gam, present_dqsdt, present_enthalpy
-
- vlen = dim1 * dim2
- present_gam = present(gam)
- present_dqsdt = present(dqsdt)
- present_enthalpy = present(enthalpy)
-
- !$acc data copyin (t,p) &
- !$acc copyout (es,qs,gam,dqsdt,enthalpy) &
- !$acc create (hltalt,tterm)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector collapse(2)
- do k = 1, dim2
- do i = 1, dim1
- tterm(i,k) = 0._r8
- end do
- end do
- !$acc end parallel
-
- call wv_sat_qsat_water_vect(t, p, es, qs, vlen)
-
- if (present_gam .or. present_dqsdt .or. present_enthalpy) then
-
- ! "generalized" analytic expression for t derivative of es
- ! accurate to within 1 percent for 173.16 < t < 373.16
- call no_ip_hltalt_vect(t, hltalt, vlen)
-
- if (present_enthalpy) call tq_enthalpy_vect(t, qs, hltalt, enthalpy, vlen)
-
- ! For pure water/ice transition term is 0.
- call deriv_outputs_vect(t, p, es, qs, hltalt, tterm, vlen, &
- gam=gam, dqsdt=dqsdt)
-
- end if
-
- !$acc end data
-end subroutine qsat_water_2D
-
-subroutine qsat_ice_line(t, p, es, qs, gam, dqsdt, enthalpy)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate SVP over ice at a given temperature, and then !
- ! calculate and return saturation specific humidity. !
- ! Optionally return various temperature derivatives or enthalpy !
- ! at saturation. !
- !------------------------------------------------------------------!
-
- use wv_sat_methods, only: wv_sat_qsat_ice
-
- ! Inputs
- real(r8), intent(in) :: t ! Temperature
- real(r8), intent(in) :: p ! Pressure
- ! Outputs
- real(r8), intent(out) :: es ! Saturation vapor pressure
- real(r8), intent(out) :: qs ! Saturation specific humidity
-
- real(r8), intent(out), optional :: gam ! (l/cpair)*(d(qs)/dt)
- real(r8), intent(out), optional :: dqsdt ! (d(qs)/dt)
- real(r8), intent(out), optional :: enthalpy ! cpair*t + hltalt*q
-
- ! Local variables
- real(r8) :: hltalt ! Modified latent heat for T derivatives
-
- call wv_sat_qsat_ice(t, p, es, qs)
-
- if (present(gam) .or. present(dqsdt) .or. present(enthalpy)) then
-
- ! For pure ice, just add latent heats.
- hltalt = latvap + latice
-
- if (present(enthalpy)) enthalpy = tq_enthalpy(t, qs, hltalt)
-
- ! For pure water/ice transition term is 0.
- call deriv_outputs_line(t, p, es, qs, hltalt, 0._r8, &
- gam=gam, dqsdt=dqsdt)
-
- end if
-
-end subroutine qsat_ice_line
-
-subroutine qsat_ice_vect(t, p, es, qs, vlen, gam, dqsdt, enthalpy)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate SVP over ice at a given temperature, and then !
- ! calculate and return saturation specific humidity. !
- ! Optionally return various temperature derivatives or enthalpy !
- ! at saturation. !
- !------------------------------------------------------------------!
-
- use wv_sat_methods, only: wv_sat_qsat_ice_vect
-
- ! Inputs
- integer, intent(in) :: vlen
- real(r8), dimension(vlen), intent(in) :: t ! Temperature
- real(r8), dimension(vlen), intent(in) :: p ! Pressure
- ! Outputs
- real(r8), dimension(vlen), intent(out) :: es ! Saturation vapor pressure
- real(r8), dimension(vlen), intent(out) :: qs ! Saturation specific humidity
-
- real(r8), dimension(vlen), intent(out), optional :: gam ! (l/cpair)*(d(qs)/dt)
- real(r8), dimension(vlen), intent(out), optional :: dqsdt ! (d(qs)/dt)
- real(r8), dimension(vlen), intent(out), optional :: enthalpy ! cpair*t + hltalt*q
-
- ! Local variables
- real(r8), dimension(vlen) :: hltalt ! Modified latent heat for T derivatives
- real(r8), dimension(vlen) :: tterm
- integer :: i
- logical :: present_gam, present_dqsdt, present_enthalpy
-
- present_gam = present(gam)
- present_dqsdt = present(dqsdt)
- present_enthalpy = present(enthalpy)
-
- !$acc data copyin (t,p) &
- !$acc copyout (es,qs,gam,dqsdt,enthalpy) &
- !$acc create (hltalt,tterm)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i = 1, vlen
- tterm(i) = 0._r8
- end do
- !$acc end parallel
-
- call wv_sat_qsat_ice_vect(t, p, es, qs, vlen)
-
- if (present_gam .or. present_dqsdt .or. present_enthalpy) then
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector
- do i = 1, vlen
- ! For pure ice, just add latent heats.
- hltalt(i) = latvap + latice
- end do
- !$acc end parallel
-
- if (present_enthalpy) call tq_enthalpy_vect(t, qs, hltalt, enthalpy, vlen)
-
- ! For pure water/ice transition term is 0.
- call deriv_outputs_vect(t, p, es, qs, hltalt, tterm, vlen, &
- gam=gam, dqsdt=dqsdt)
-
- end if
-
- !$acc end data
-end subroutine qsat_ice_vect
-
-subroutine qsat_ice_2D(t, p, es, qs, dim1, dim2, gam, dqsdt, enthalpy)
- !------------------------------------------------------------------!
- ! Purpose: !
- ! Calculate SVP over ice at a given temperature, and then !
- ! calculate and return saturation specific humidity. !
- ! Optionally return various temperature derivatives or enthalpy !
- ! at saturation. !
- !------------------------------------------------------------------!
-
- use wv_sat_methods, only: wv_sat_qsat_ice_vect
-
- ! Inputs
- integer, intent(in) :: dim1, dim2
- real(r8), dimension(dim1,dim2), intent(in) :: t ! Temperature
- real(r8), dimension(dim1,dim2), intent(in) :: p ! Pressure
- ! Outputs
- real(r8), dimension(dim1,dim2), intent(out) :: es ! Saturation vapor pressure
- real(r8), dimension(dim1,dim2), intent(out) :: qs ! Saturation specific humidity
-
- real(r8), dimension(dim1,dim2), intent(out), optional :: gam ! (l/cpair)*(d(qs)/dt)
- real(r8), dimension(dim1,dim2), intent(out), optional :: dqsdt ! (d(qs)/dt)
- real(r8), dimension(dim1,dim2), intent(out), optional :: enthalpy ! cpair*t + hltalt*q
-
- ! Local variables
- real(r8), dimension(dim1,dim2) :: hltalt ! Modified latent heat for T derivatives
- real(r8), dimension(dim1,dim2) :: tterm
- integer :: i, k, vlen
- logical :: present_gam, present_dqsdt, present_enthalpy
-
- vlen = dim1 * dim2
- present_gam = present(gam)
- present_dqsdt = present(dqsdt)
- present_enthalpy = present(enthalpy)
-
- !$acc data copyin (t,p) &
- !$acc copyout (es,qs,gam,dqsdt,enthalpy) &
- !$acc create (hltalt,tterm)
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector collapse(2)
- do k = 1, dim2
- do i = 1, dim1
- tterm(i,k) = 0._r8
- end do
- end do
- !$acc end parallel
-
- call wv_sat_qsat_ice_vect(t, p, es, qs, vlen)
-
- if (present_gam .or. present_dqsdt .or. present_enthalpy) then
-
- !$acc parallel vector_length(VLENS) default(present)
- !$acc loop gang vector collapse(2)
- do k = 1, dim2
- do i = 1, dim1
- ! For pure ice, just add latent heats.
- hltalt(i,k) = latvap + latice
- end do
- end do
- !$acc end parallel
-
- if (present_enthalpy) call tq_enthalpy_vect(t, qs, hltalt, enthalpy, vlen)
-
- ! For pure water/ice transition term is 0.
- call deriv_outputs_vect(t, p, es, qs, hltalt, tterm, vlen, &
- gam=gam, dqsdt=dqsdt)
-
- end if
-
- !$acc end data
-end subroutine qsat_ice_2D
-
-!---------------------------------------------------------------------
-! FINDSP (WET BULB TEMPERATURE) PROCEDURES
-!---------------------------------------------------------------------
-
-subroutine findsp_vc(q, t, p, use_ice, tsp, qsp)
-
- use cam_logfile, only: iulog
- use cam_abortutils, only: endrun
-
- ! Wrapper for findsp which is 1D and handles the output status.
- ! Changing findsp to elemental restricted debugging output.
- ! If that output is needed again, it's preferable *not* to copy findsp,
- ! but to change the existing version.
-
- ! input arguments
- real(r8), intent(in) :: q(:) ! water vapor (kg/kg)
- real(r8), intent(in) :: t(:) ! temperature (K)
- real(r8), intent(in) :: p(:) ! pressure (Pa)
- logical, intent(in) :: use_ice ! flag to include ice phase in calculations
-
- ! output arguments
- real(r8), intent(out) :: tsp(:) ! saturation temp (K)
- real(r8), intent(out) :: qsp(:) ! saturation mixing ratio (kg/kg)
-
- integer :: status(size(q)) ! flag representing state of output
- ! 0 => Successful convergence
- ! 1 => No calculation done: pressure or specific
- ! humidity not within usable range
- ! 2 => Run failed to converge
- ! 4 => Temperature fell below minimum
- ! 8 => Enthalpy not conserved
-
- integer :: n, i
-
- n = size(q)
-
- ! Currently, only 2 and 8 seem to be treated as fatal errors.
- do i = 1,n
- call findsp(q(i), t(i), p(i), use_ice, tsp(i), qsp(i), status(i))
- if (status(i) == 2) then
- write(iulog,*) ' findsp not converging at i = ', i
- write(iulog,*) ' t, q, p ', t(i), q(i), p(i)
- write(iulog,*) ' tsp, qsp ', tsp(i), qsp(i)
- call endrun ('wv_saturation::FINDSP -- not converging')
- else if (status(i) == 8) then
- write(iulog,*) ' the enthalpy is not conserved at i = ', i
- write(iulog,*) ' t, q, p ', t(i), q(i), p(i)
- write(iulog,*) ' tsp, qsp ', tsp(i), qsp(i)
- call endrun ('wv_saturation::FINDSP -- enthalpy is not conserved')
- endif
- end do
-
-end subroutine findsp_vc
-
-subroutine findsp (q, t, p, use_ice, tsp, qsp, status)
-!-----------------------------------------------------------------------
-!
-! Purpose:
-! find the wet bulb temperature for a given t and q
-! in a longitude height section
-! wet bulb temp is the temperature and spec humidity that is
-! just saturated and has the same enthalpy
-! if q > qs(t) then tsp > t and qsp = qs(tsp) < q
-! if q < qs(t) then tsp < t and qsp = qs(tsp) > q
-!
-! Method:
-! a Newton method is used
-! first guess uses an algorithm provided by John Petch from the UKMO
-! we exclude points where the physical situation is unrealistic
-! e.g. where the temperature is outside the range of validity for the
-! saturation vapor pressure, or where the water vapor pressure
-! exceeds the ambient pressure, or the saturation specific humidity is
-! unrealistic
-!
-! Author: P. Rasch
-!
-!-----------------------------------------------------------------------
-!
-! input arguments
-!
-
- real(r8), intent(in) :: q ! water vapor (kg/kg)
- real(r8), intent(in) :: t ! temperature (K)
- real(r8), intent(in) :: p ! pressure (Pa)
- logical, intent(in) :: use_ice ! flag to include ice phase in calculations
-!
-! output arguments
-!
- real(r8), intent(out) :: tsp ! saturation temp (K)
- real(r8), intent(out) :: qsp ! saturation mixing ratio (kg/kg)
- integer, intent(out) :: status ! flag representing state of output
- ! 0 => Successful convergence
- ! 1 => No calculation done: pressure or specific
- ! humidity not within usable range
- ! 2 => Run failed to converge
- ! 4 => Temperature fell below minimum
- ! 8 => Enthalpy not conserved
-!
-! local variables
-!
- integer, parameter :: iter = 8 ! max number of times to iterate the calculation
- integer :: l ! iterator
-
- real(r8) es ! sat. vapor pressure
- real(r8) gam ! change in sat spec. hum. wrt temperature (times hltalt/cpair)
- real(r8) dgdt ! work variable
- real(r8) g ! work variable
- real(r8) hltalt ! lat. heat. of vap.
- real(r8) qs ! spec. hum. of water vapor
-
-! work variables
- real(r8) t1, q1, dt, dq
- real(r8) qvd
- real(r8) r1b, c1, c2
- real(r8), parameter :: dttol = 1.e-4_r8 ! the relative temp error tolerance required to quit the iteration
- real(r8), parameter :: dqtol = 1.e-4_r8 ! the relative moisture error tolerance required to quit the iteration
- real(r8) enin, enout
-
- ! Saturation specific humidity at this temperature
- if (use_ice) then
- call qsat(t, p, es, qs)
- else
- call qsat_water(t, p, es, qs)
- end if
-
- ! make sure a meaningful calculation is possible
- if (p <= 5._r8*es .or. qs <= 0._r8 .or. qs >= 0.5_r8 &
- .or. t < tmin .or. t > tmax) then
- status = 1
- ! Keep initial parameters when conditions aren't suitable
- tsp = t
- qsp = q
- enin = 1._r8
- enout = 1._r8
-
- return
- end if
-
- ! Prepare to iterate
- status = 2
-
- ! Get initial enthalpy
- if (use_ice) then
- call calc_hltalt(t,hltalt)
- else
- call no_ip_hltalt(t,hltalt)
- end if
- enin = tq_enthalpy(t, q, hltalt)
-
- ! make a guess at the wet bulb temp using a UKMO algorithm (from J. Petch)
- c1 = hltalt*c3
- c2 = (t + 36._r8)**2
- r1b = c2/(c2 + c1*qs)
- qvd = r1b * (q - qs)
- tsp = t + ((hltalt/cpair)*qvd)
-
- ! Generate qsp, gam, and enout from tsp.
- if (use_ice) then
- call qsat(tsp, p, es, qsp, gam=gam, enthalpy=enout)
- else
- call qsat_water(tsp, p, es, qsp, gam=gam, enthalpy=enout)
- end if
-
- ! iterate on first guess
- do l = 1, iter
-
- g = enin - enout
- dgdt = -cpair * (1 + gam)
-
- ! New tsp
- t1 = tsp - g/dgdt
- dt = abs(t1 - tsp)/t1
- tsp = t1
-
- ! bail out if past end of temperature range
- if ( tsp < tmin ) then
- tsp = tmin
- ! Get latent heat and set qsp to a value
- ! that preserves enthalpy.
- if (use_ice) then
- call calc_hltalt(tsp,hltalt)
- else
- call no_ip_hltalt(tsp,hltalt)
- end if
- qsp = (enin - cpair*tsp)/hltalt
- enout = tq_enthalpy(tsp, qsp, hltalt)
- status = 4
- exit
- end if
-
- ! Re-generate qsp, gam, and enout from new tsp.
- if (use_ice) then
- call qsat(tsp, p, es, q1, gam=gam, enthalpy=enout)
- else
- call qsat_water(tsp, p, es, q1, gam=gam, enthalpy=enout)
- end if
- dq = abs(q1 - qsp)/max(q1,1.e-12_r8)
- qsp = q1
-
- ! if converged at this point, exclude it from more iterations
- if (dt < dttol .and. dq < dqtol) then
- status = 0
- exit
- endif
- end do
-
- ! Test for enthalpy conservation
- if (abs((enin-enout)/(enin+enout)) > 1.e-4_r8) status = 8
-
-end subroutine findsp
-
-end module wv_saturation
diff --git a/src/physics/cam/zm_conv_intr.F90 b/src/physics/cam/zm_conv_intr.F90
index 5db6d1bc03..44488ac737 100644
--- a/src/physics/cam/zm_conv_intr.F90
+++ b/src/physics/cam/zm_conv_intr.F90
@@ -14,6 +14,7 @@ module zm_conv_intr
use zm_convr, only: zm_convr_init, zm_convr_run
use zm_conv_convtran, only: zm_conv_convtran_run
use zm_conv_momtran, only: zm_conv_momtran_run
+ use cloud_fraction_fice, only: cloud_fraction_fice_run
use rad_constituents, only: rad_cnst_get_info, rad_cnst_get_mode_num, rad_cnst_get_aer_mmr, &
rad_cnst_get_aer_props, rad_cnst_get_mode_props !, &
@@ -23,6 +24,8 @@ module zm_conv_intr
use perf_mod
use cam_logfile, only: iulog
use constituents, only: cnst_add
+ use ref_pres, only: trop_cloud_top_lev
+ use phys_control, only: phys_getopts
implicit none
private
@@ -37,7 +40,7 @@ module zm_conv_intr
zm_conv_tend, &! return tendencies
zm_conv_tend_2 ! return tendencies
- public zmconv_ke, zmconv_ke_lnd, zmconv_org ! needed by convect_shallow
+ public zmconv_ke, zmconv_ke_lnd ! needed by convect_shallow
integer ::& ! indices for fields in the physics buffer
zm_mu_idx, &
@@ -52,11 +55,9 @@ module zm_conv_intr
zm_ideep_idx, &
dp_flxprc_idx, &
dp_flxsnw_idx, &
- ixorg, &
+ dp_cldliq_idx, &
+ dp_cldice_idx, &
dlfzm_idx, & ! detrained convective cloud water mixing ratio.
- difzm_idx, & ! detrained convective cloud ice mixing ratio.
- dnlfzm_idx, & ! detrained convective cloud water num concen.
- dnifzm_idx, & ! detrained convective cloud ice num concen.
prec_dp_idx, &
snow_dp_idx, &
mconzm_idx ! convective mass flux
@@ -70,8 +71,6 @@ module zm_conv_intr
real(r8) :: zmconv_momcd = unset_r8
integer :: zmconv_num_cin ! Number of negative buoyancy regions that are allowed
! before the convection top and CAPE calculations are completed.
- logical :: zmconv_org ! Parameterization for sub-grid scale convective organization for the ZM deep
- ! convective scheme based on Mapes and Neale (2011)
real(r8) :: zmconv_dmpdz = unset_r8 ! Parcel fractional mass entrainment rate
real(r8) :: zmconv_tiedke_add = unset_r8 ! Convective parcel temperature perturbation
real(r8) :: zmconv_capelmt = unset_r8 ! Triggering thereshold for ZM convection
@@ -141,15 +140,9 @@ subroutine zm_conv_register
! detrained convective cloud water mixing ratio.
call pbuf_add_field('DLFZM', 'physpkg', dtype_r8, (/pcols,pver/), dlfzm_idx)
- ! detrained convective cloud ice mixing ratio.
- call pbuf_add_field('DIFZM', 'physpkg', dtype_r8, (/pcols,pver/), difzm_idx)
! convective mass fluxes
call pbuf_add_field('CMFMC_DP', 'physpkg', dtype_r8, (/pcols,pverp/), mconzm_idx)
-!CACNOTE - Is zm_org really a constituent or was it just a handy structure to use for an allocatable which persists in the run?
- if (zmconv_org) then
- call cnst_add('ZM_ORG',0._r8,0._r8,0._r8,ixorg,longname='organization parameter')
- endif
end subroutine zm_conv_register
@@ -159,7 +152,6 @@ subroutine zm_conv_readnl(nlfile)
use spmd_utils, only: mpicom, masterproc, masterprocid, mpi_real8, mpi_integer, mpi_logical
use namelist_utils, only: find_group_name
- use units, only: getunit, freeunit
character(len=*), intent(in) :: nlfile ! filepath for file containing namelist input
@@ -168,15 +160,14 @@ subroutine zm_conv_readnl(nlfile)
character(len=*), parameter :: subname = 'zm_conv_readnl'
namelist /zmconv_nl/ zmconv_c0_lnd, zmconv_c0_ocn, zmconv_num_cin, &
- zmconv_ke, zmconv_ke_lnd, zmconv_org, &
+ zmconv_ke, zmconv_ke_lnd, &
zmconv_momcu, zmconv_momcd, &
zmconv_dmpdz, zmconv_tiedke_add, zmconv_capelmt, &
zmconv_parcel_pbl, zmconv_tau
!-----------------------------------------------------------------------------
if (masterproc) then
- unitn = getunit()
- open( unitn, file=trim(nlfile), status='old' )
+ open( newunit=unitn, file=trim(nlfile), status='old' )
call find_group_name(unitn, 'zmconv_nl', status=ierr)
if (ierr == 0) then
read(unitn, zmconv_nl, iostat=ierr)
@@ -185,7 +176,6 @@ subroutine zm_conv_readnl(nlfile)
end if
end if
close(unitn)
- call freeunit(unitn)
end if
@@ -204,8 +194,6 @@ subroutine zm_conv_readnl(nlfile)
if (ierr /= 0) call endrun("zm_conv_readnl: FATAL: mpi_bcast: zmconv_momcu")
call mpi_bcast(zmconv_momcd, 1, mpi_real8, masterprocid, mpicom, ierr)
if (ierr /= 0) call endrun("zm_conv_readnl: FATAL: mpi_bcast: zmconv_momcd")
- call mpi_bcast(zmconv_org, 1, mpi_logical, masterprocid, mpicom, ierr)
- if (ierr /= 0) call endrun("zm_conv_readnl: FATAL: mpi_bcast: zmconv_org")
call mpi_bcast(zmconv_dmpdz, 1, mpi_real8, masterprocid, mpicom, ierr)
if (ierr /= 0) call endrun("zm_conv_readnl: FATAL: mpi_bcast: zmconv_dmpdz")
call mpi_bcast(zmconv_tiedke_add, 1, mpi_real8, masterprocid, mpicom, ierr)
@@ -260,10 +248,6 @@ subroutine zm_conv_init(pref_edge)
! Register fields with the output buffer
!
- if (zmconv_org) then
- call addfld ('ZM_ORG ', (/ 'lev' /), 'A', '- ','Organization parameter')
- call addfld ('ZM_ORG2D ', (/ 'lev' /), 'A', '- ','Organization parameter 2D')
- endif
call addfld ('PRECZ', horiz_only, 'A', 'm/s','total precipitation from ZM convection')
call addfld ('ZMDT', (/ 'lev' /), 'A', 'K/s','T tendency - Zhang-McFarlane moist convection')
call addfld ('ZMDQ', (/ 'lev' /), 'A', 'kg/kg/s','Q tendency - Zhang-McFarlane moist convection')
@@ -306,16 +290,11 @@ subroutine zm_conv_init(pref_edge)
call addfld ('ZMICVU', (/ 'lev' /), 'A', 'm/s', 'ZM in-cloud V updrafts')
call addfld ('ZMICVD', (/ 'lev' /), 'A', 'm/s', 'ZM in-cloud V downdrafts')
- call addfld ('DIFZM' ,(/ 'lev' /), 'A','kg/kg/s ','Detrained ice water from ZM convection')
call addfld ('DLFZM' ,(/ 'lev' /), 'A','kg/kg/s ','Detrained liquid water from ZM convection')
call phys_getopts( history_budget_out = history_budget, &
history_budget_histfile_num_out = history_budget_histfile_num)
- if (zmconv_org) then
- call add_default('ZM_ORG', 1, ' ')
- call add_default('ZM_ORG2D', 1, ' ')
- endif
if ( history_budget ) then
call add_default('EVAPTZM ', history_budget_histfile_num, ' ')
call add_default('EVAPQZM ', history_budget_histfile_num, ' ')
@@ -362,14 +341,17 @@ subroutine zm_conv_init(pref_edge)
end if
no_deep_pbl = phys_deepconv_pbl()
-!CACNOTE - Need to check errflg and report errors
- call zm_convr_init(cpair, epsilo, gravit, latvap, tmelt, rair, &
- limcnv,zmconv_c0_lnd, zmconv_c0_ocn, zmconv_ke, zmconv_ke_lnd, &
- zmconv_momcu, zmconv_momcd, zmconv_num_cin, zmconv_org, &
+ call zm_convr_init(plev, plevp, cpair, epsilo, gravit, latvap, tmelt, rair, &
+ pref_edge,zmconv_c0_lnd, zmconv_c0_ocn, zmconv_ke, zmconv_ke_lnd, &
+ zmconv_momcu, zmconv_momcd, zmconv_num_cin, &
no_deep_pbl, zmconv_tiedke_add, &
zmconv_capelmt, zmconv_dmpdz,zmconv_parcel_pbl, zmconv_tau, &
masterproc, iulog, errmsg, errflg)
+ if (errflg /= 0) then
+ call endrun('From zm_convr_init:' // errmsg)
+ end if
+
cld_idx = pbuf_get_index('CLD')
fracis_idx = pbuf_get_index('FRACIS')
@@ -394,8 +376,10 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
use physics_buffer, only : pbuf_get_field, physics_buffer_desc, pbuf_old_tim_idx
use constituents, only: pcnst, cnst_get_ind, cnst_is_convtran1
use physconst, only: gravit, latice, latvap, tmelt, cpwv, cpliq, rh2o
+ use phys_grid, only: get_rlat_all_p, get_rlon_all_p
use phys_control, only: cam_physpkg_is
+ use ccpp_constituent_prop_mod, only: ccpp_const_props
! Arguments
@@ -451,9 +435,6 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
real(r8), pointer, dimension(:,:) :: flxprec ! Convective-scale flux of precip at interfaces (kg/m2/s)
real(r8), pointer, dimension(:,:) :: flxsnow ! Convective-scale flux of snow at interfaces (kg/m2/s)
real(r8), pointer :: dlf(:,:) ! detrained convective cloud water mixing ratio.
- real(r8), pointer :: dif(:,:) ! detrained convective cloud ice mixing ratio.
- real(r8), pointer :: dnlf(:,:) ! detrained convective cloud water num concen.
- real(r8), pointer :: dnif(:,:) ! detrained convective cloud ice num concen.
real(r8), pointer :: lambdadpcu(:,:) ! slope of cloud liquid size distr
real(r8), pointer :: mudpcu(:,:) ! width parameter of droplet size distr
real(r8), pointer :: mconzm(:,:) !convective mass fluxes
@@ -475,10 +456,13 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
real(r8) :: pcont(pcols), pconb(pcols), freqzm(pcols)
+ real(r8) :: lat_all(pcols), long_all(pcols)
+
! history output fields
real(r8) :: cape(pcols) ! w convective available potential energy.
real(r8) :: mu_out(pcols,pver)
real(r8) :: md_out(pcols,pver)
+ real(r8) :: dif(pcols,pver)
! used in momentum transport calculation
real(r8) :: pguallu(pcols, pver)
@@ -490,17 +474,18 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
real(r8) :: icwdu(pcols,pver)
real(r8) :: icwdv(pcols,pver)
real(r8) :: seten(pcols, pver)
- logical :: l_windt(2)
+ logical :: l_windt
real(r8) :: tfinal1, tfinal2
integer :: ii
- real(r8),pointer :: zm_org2d(:,:)
- real(r8),allocatable :: orgt_alloc(:,:), org_alloc(:,:)
-
- real(r8) :: zm_org2d_ncol(state%ncol,pver)
- real(r8) :: orgt_ncol(state%ncol,pver), org_ncol(state%ncol,pver)
+ real(r8) :: fice(pcols,pver)
+ real(r8) :: fsnow_conv(pcols,pver)
logical :: lq(pcnst)
+ character(len=16) :: macrop_scheme
+ character(len=40) :: scheme_name
+ character(len=40) :: str
+ integer :: top_lev
!----------------------------------------------------------------------
@@ -517,9 +502,6 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
lq(:) = .FALSE.
lq(1) = .TRUE.
- if (zmconv_org) then
- lq(ixorg) = .TRUE.
- endif
call physics_ptend_init(ptend_loc, state%psetcols, 'zm_convr_run', ls=.true., lq=lq)! initialize local ptend type
!
@@ -547,35 +529,22 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
call pbuf_get_field(pbuf, zm_ideep_idx, ideep)
call pbuf_get_field(pbuf, dlfzm_idx, dlf)
- call pbuf_get_field(pbuf, difzm_idx, dif)
call pbuf_get_field(pbuf, mconzm_idx, mconzm)
- allocate(dnlf(pcols,pver), dnif(pcols,pver))
-
-!
! Begin with Zhang-McFarlane (1996) convection parameterization
!
call t_startf ('zm_convr_run')
- if (zmconv_org) then
- allocate(zm_org2d(pcols,pver))
- allocate(org_alloc(ncol,pver))
- allocate(orgt_alloc(ncol,pver))
- org_ncol(:ncol,:) = state%q(1:ncol,:,ixorg)
- endif
-
!REMOVECAM - no longer need these when CAM is retired and pcols no longer exists
ptend_loc%q(:,:,1) = 0._r8
ptend_loc%s(:,:) = 0._r8
+ dif(:,:) = 0._r8
mcon(:,:) = 0._r8
dlf(:,:) = 0._r8
cme(:,:) = 0._r8
cape(:) = 0._r8
zdu(:,:) = 0._r8
rprd(:,:) = 0._r8
- dif(:,:) = 0._r8
- dnlf(:,:) = 0._r8
- dnif(:,:) = 0._r8
mu(:,:) = 0._r8
eu(:,:) = 0._r8
du(:,:) = 0._r8
@@ -591,29 +560,27 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
ideep(:) = 0._r8
!REMOVECAM_END
-!CACNOTE - Need to check errflg and report errors
+
+ call get_rlat_all_p(lchnk, ncol, lat_all)
+ call get_rlon_all_p(lchnk, ncol, long_all)
+
call zm_convr_run(ncol, pver, &
pverp, gravit, latice, cpwv, cpliq, rh2o, &
+ lat_all, long_all, &
state%t(:ncol,:), state%q(:ncol,:,1), prec(:ncol), &
pblh(:ncol), state%zm(:ncol,:), state%phis(:ncol), state%zi(:ncol,:), ptend_loc%q(:ncol,:,1), &
ptend_loc%s(:ncol,:), state%pmid(:ncol,:), state%pint(:ncol,:), state%pdel(:ncol,:), &
ztodt, mcon(:ncol,:), cme(:ncol,:), cape(:ncol), &
- tpert(:ncol), dlf(:ncol,:), zdu(:ncol,:), rprd(:ncol,:), &
+ tpert(:ncol), dlf(:ncol,:), dif(:ncol,:), zdu(:ncol,:), rprd(:ncol,:), &
mu(:ncol,:), md(:ncol,:), du(:ncol,:), eu(:ncol,:), ed(:ncol,:), &
dp(:ncol,:), dsubcld(:ncol), jt(:ncol), maxg(:ncol), ideep(:ncol), &
ql(:ncol,:), rliq(:ncol), landfrac(:ncol), &
- org_ncol(:ncol,:), orgt_ncol(:ncol,:), zm_org2d_ncol(:ncol,:), &
- dif(:ncol,:), dnlf(:ncol,:), dnif(:ncol,:), &
- rice(:ncol), errmsg, errflg)
+ rice(:ncol), lengath, scheme_name, errmsg, errflg)
-
- if (zmconv_org) then
- ptend_loc%q(:,:,ixorg)=orgt_ncol(:ncol,:)
- zm_org2d(:ncol,:) = zm_org2d_ncol(:ncol,:)
- endif
-
- lengath = count(ideep > 0)
- if (lengath > ncol) lengath = ncol ! should not happen, but force it to not be larger than ncol for safety sake
+ if (errflg /= 0) then
+ write(str,*) 'From zm_convr_run: at chunk ',lchnk, ' : '
+ call endrun(str // errmsg)
+ end if
jctop(:) = real(pver,r8)
jcbot(:) = 1._r8
@@ -657,7 +624,6 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
call outfld('ZMDQ ',ptend_loc%q(1,1,1) ,pcols ,lchnk )
call t_stopf ('zm_convr_run')
- call outfld('DIFZM' ,dif ,pcols, lchnk)
call outfld('DLFZM' ,dlf ,pcols, lchnk)
pcont(:ncol) = state%ps(:ncol)
@@ -683,9 +649,6 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
! initialize ptend for next process
lq(:) = .FALSE.
lq(1) = .TRUE.
- if (zmconv_org) then
- lq(ixorg) = .TRUE.
- endif
call physics_ptend_init(ptend_loc, state1%psetcols, 'zm_conv_evap_run', ls=.true., lq=lq)
call t_startf ('zm_conv_evap_run')
@@ -702,25 +665,28 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
flxprec(:,:) = 0._r8
flxsnow(:,:) = 0._r8
snow(:) = 0._r8
+ fice(:,:) = 0._r8
+ fsnow_conv(:,:) = 0._r8
!REMOVECAM_END
+ top_lev = 1
+ call phys_getopts (macrop_scheme_out = macrop_scheme)
+ if ( .not. (macrop_scheme == "rk" .or. macrop_scheme == "SPCAM_sam1mom")) top_lev = trop_cloud_top_lev
+
+ call cloud_fraction_fice_run(ncol, state1%t(:ncol,:), tmelt, top_lev, pver, fice(:ncol,:), fsnow_conv(:ncol,:))
+
call zm_conv_evap_run(state1%ncol, pver, pverp, &
gravit, latice, latvap, tmelt, &
- cpair, zmconv_ke, zmconv_ke_lnd, zmconv_org, &
+ cpair, zmconv_ke, zmconv_ke_lnd, &
state1%t(:ncol,:),state1%pmid(:ncol,:),state1%pdel(:ncol,:),state1%q(:ncol,:pver,1), &
landfrac(:ncol), &
ptend_loc%s(:ncol,:), tend_s_snwprd(:ncol,:), tend_s_snwevmlt(:ncol,:), ptend_loc%q(:ncol,:pver,1), &
rprd(:ncol,:), cld(:ncol,:), ztodt, &
- prec(:ncol), snow(:ncol), ntprprd(:ncol,:), ntsnprd(:ncol,:), flxprec(:ncol,:), flxsnow(:ncol,:) )
+ prec(:ncol), snow(:ncol), ntprprd(:ncol,:), ntsnprd(:ncol,:), fsnow_conv(:ncol,:), flxprec(:ncol,:), flxsnow(:ncol,:),&
+ scheme_name, errmsg, errflg)
evapcdp(:ncol,:pver) = ptend_loc%q(:ncol,:pver,1)
- if (zmconv_org) then
- ptend_loc%q(:ncol,:pver,ixorg) = min(1._r8,max(0._r8,(50._r8*1000._r8*1000._r8*abs(evapcdp(:ncol,:pver))) &
- -(state%q(:ncol,:pver,ixorg)/10800._r8)))
- ptend_loc%q(:ncol,:pver,ixorg) = (ptend_loc%q(:ncol,:pver,ixorg) - state%q(:ncol,:pver,ixorg))/ztodt
- endif
-
!
! Write out variables from zm_conv_evap_run
!
@@ -755,8 +721,7 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
call physics_ptend_init(ptend_loc, state1%psetcols, 'zm_conv_momtran_run', ls=.true., lu=.true., lv=.true.)
- l_windt(1) = .true.
- l_windt(2) = .true.
+ l_windt = .true.
!REMOVECAM - no longer need these when CAM is retired and pcols no longer exists
ptend_loc%s(:,:) = 0._r8
ptend_loc%u(:,:) = 0._r8
@@ -765,15 +730,16 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
call t_startf ('zm_conv_momtran_run')
- call zm_conv_momtran_run (ncol, pver, pverp, &
- l_windt,state1%u(:ncol,:), state1%v(:ncol,:), 2, mu(:ncol,:), md(:ncol,:), &
+ call zm_conv_momtran_run (ncol, pver, pverp, &
+ l_windt,state1%u(:ncol,:), state1%v(:ncol,:), mu(:ncol,:), md(:ncol,:), &
zmconv_momcu, zmconv_momcd, &
du(:ncol,:), eu(:ncol,:), ed(:ncol,:), dp(:ncol,:), dsubcld(:ncol), &
jt(:ncol), maxg(:ncol), ideep(:ncol), 1, lengath, &
nstep, ptend_loc%u(:ncol,:), ptend_loc%v(:ncol,:),&
pguallu(:ncol,:), pguallv(:ncol,:), pgdallu(:ncol,:), pgdallv(:ncol,:), &
- icwuu(:ncol,:), icwuv(:ncol,:), icwdu(:ncol,:), icwdv(:ncol,:), ztodt, seten(:ncol,:) )
- call t_stopf ('zm_conv_momtran_run')
+ icwuu(:ncol,:), icwuv(:ncol,:), icwdu(:ncol,:), icwdv(:ncol,:), ztodt, seten(:ncol,:), &
+ scheme_name, errmsg, errflg)
+ call t_stopf ('zm_conv_momtran_run')
ptend_loc%s(:ncol,:pver) = seten(:ncol,:pver)
@@ -786,12 +752,8 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
! update physics state type state1 with ptend_loc
call physics_update(state1, ptend_loc, ztodt)
- ftem(:ncol,:pver) = seten(:ncol,:pver)/cpair
- if (zmconv_org) then
- call outfld('ZM_ORG', state%q(:,:,ixorg), pcols, lchnk)
- call outfld('ZM_ORG2D', zm_org2d, pcols, lchnk)
- endif
- call outfld('ZMMTT', ftem , pcols, lchnk)
+ ftem(:ncol,:pver) = seten(:ncol,:pver)/cpair
+ call outfld('ZMMTT', ftem , pcols, lchnk)
! Output apparent force from pressure gradient
call outfld('ZMUPGU', pguallu, pcols, lchnk)
@@ -828,7 +790,8 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
ptend_loc%lq,state1%q(:ncol,:,:), pcnst, mu(:ncol,:), md(:ncol,:), &
du(:ncol,:), eu(:ncol,:), ed(:ncol,:), dp(:ncol,:), dsubcld(:ncol), &
jt(:ncol), maxg(:ncol), ideep(:ncol), 1, lengath, &
- nstep, fracis(:ncol,:,:), ptend_loc%q(:ncol,:,:), fake_dpdry(:ncol,:))
+ nstep, fracis(:ncol,:,:), ptend_loc%q(:ncol,:,:), fake_dpdry(:ncol,:), ccpp_const_props, &
+ scheme_name, errmsg, errflg)
call t_stopf ('convtran1')
call outfld('ZMDICE ',ptend_loc%q(1,1,ixcldice) ,pcols ,lchnk )
@@ -840,11 +803,7 @@ subroutine zm_conv_tend(pblh ,mcon ,cme , &
call physics_state_dealloc(state1)
call physics_ptend_dealloc(ptend_loc)
- if (zmconv_org) then
- deallocate(zm_org2d)
- end if
- deallocate(dnlf, dnif)
end subroutine zm_conv_tend
!=========================================================================================
@@ -856,6 +815,8 @@ subroutine zm_conv_tend_2( state, ptend, ztodt, pbuf)
use time_manager, only: get_nstep
use physics_buffer, only: pbuf_get_index, pbuf_get_field, physics_buffer_desc
use constituents, only: pcnst, cnst_is_convtran2
+ use ccpp_constituent_prop_mod, only: ccpp_const_props
+
! Arguments
type(physics_state), intent(in ) :: state ! Physics state variables
@@ -885,6 +846,11 @@ subroutine zm_conv_tend_2( state, ptend, ztodt, pbuf)
integer, pointer :: jt(:) ! (pcols)
integer, pointer :: maxg(:) ! (pcols)
integer, pointer :: ideep(:) ! (pcols)
+
+ character(len=40) :: scheme_name
+ character(len=512) :: errmsg
+ integer :: errflg
+
!-----------------------------------------------------------------------------------
@@ -928,7 +894,13 @@ subroutine zm_conv_tend_2( state, ptend, ztodt, pbuf)
ptend%lq,state%q(:ncol,:,:), pcnst, mu(:ncol,:), md(:ncol,:), &
du(:ncol,:), eu(:ncol,:), ed(:ncol,:), dp(:ncol,:), dsubcld(:ncol), &
jt(:ncol), maxg(:ncol), ideep(:ncol), 1, lengath, &
- nstep, fracis(:ncol,:,:), ptend%q(:ncol,:,:), dpdry(:ncol,:))
+ nstep, fracis(:ncol,:,:), ptend%q(:ncol,:,:), dpdry(:ncol,:), ccpp_const_props, &
+ scheme_name, errmsg, errflg)
+
+ if (errflg /= 0) then
+ call endrun('From zm_conv_convtran_run:' // errmsg)
+ end if
+
call t_stopf ('convtran2')
end if
diff --git a/src/physics/cam7/physpkg.F90 b/src/physics/cam7/physpkg.F90
index abcd3ea7c2..d5ad2ee3e1 100644
--- a/src/physics/cam7/physpkg.F90
+++ b/src/physics/cam7/physpkg.F90
@@ -785,6 +785,7 @@ subroutine phys_init( phys_state, phys_tend, pbuf2d, cam_in, cam_out )
! local variables
integer :: lchnk
integer :: ierr
+ integer :: ixq
logical :: history_budget ! output tendencies and state variables for
! temperature, water vapor, cloud
@@ -955,7 +956,8 @@ subroutine phys_init( phys_state, phys_tend, pbuf2d, cam_in, cam_out )
! Initialize CAM CCPP constituent properties array
! for use in CCPP-ized physics schemes:
- call ccpp_const_props_init()
+ call cnst_get_ind('Q', ixq)
+ call ccpp_const_props_init(ixq)
! Initialize qneg3 and qneg4
call qneg_init()
diff --git a/src/physics/simple/physpkg.F90 b/src/physics/simple/physpkg.F90
index c236f7f021..48bc0200f7 100644
--- a/src/physics/simple/physpkg.F90
+++ b/src/physics/simple/physpkg.F90
@@ -208,6 +208,8 @@ subroutine phys_init( phys_state, phys_tend, pbuf2d, cam_in, cam_out )
use nudging, only: Nudge_Model, nudging_init
use cam_snapshot, only: cam_snapshot_init
use cam_budget, only: cam_budget_init
+ use constituents, only: cnst_get_ind
+
use ccpp_constituent_prop_mod, only: ccpp_const_props_init
@@ -220,7 +222,7 @@ subroutine phys_init( phys_state, phys_tend, pbuf2d, cam_in, cam_out )
type(cam_out_t),intent(inout) :: cam_out(begchunk:endchunk)
! local variables
- integer :: lchnk
+ integer :: lchnk, ixq
!-----------------------------------------------------------------------
call physics_type_alloc(phys_state, phys_tend, begchunk, endchunk, pcols)
@@ -281,7 +283,8 @@ subroutine phys_init( phys_state, phys_tend, pbuf2d, cam_in, cam_out )
! Initialize CAM CCPP constituent properties array
! for use in CCPP-ized physics schemes:
- call ccpp_const_props_init()
+ call cnst_get_ind('Q', ixq)
+ call ccpp_const_props_init(ixq)
! Initialize qneg3 and qneg4
call qneg_init()
diff --git a/src/utils/cam_ccpp/ccpp_constituent_prop_mod.F90 b/src/utils/cam_ccpp/ccpp_constituent_prop_mod.F90
index edbc4d59e9..29ae61fc53 100644
--- a/src/utils/cam_ccpp/ccpp_constituent_prop_mod.F90
+++ b/src/utils/cam_ccpp/ccpp_constituent_prop_mod.F90
@@ -8,12 +8,19 @@ module ccpp_constituent_prop_mod
type, public :: ccpp_constituent_prop_ptr_t
logical, private :: thermo_active = .false.
logical, private :: water_species = .false.
- contains
+ logical, private :: species_is_dry
+ character(len=256) :: std_name = ''
+
+ contains
procedure :: standard_name => ccp_get_standard_name
+ procedure :: set_standard_name => ccp_set_standard_name
procedure :: is_thermo_active => ccp_is_thermo_active
procedure :: is_water_species => ccp_is_water_species
procedure :: set_thermo_active => ccp_set_thermo_active
procedure :: set_water_species => ccp_set_water_species
+ procedure :: is_dry => ccp_is_dry
+ procedure :: set_dry => ccp_set_dry
+
end type ccpp_constituent_prop_ptr_t
! CCPP properties init routine
@@ -37,20 +44,43 @@ subroutine ccp_get_standard_name(this, std_name, errcode, errmsg)
integer, optional, intent(out) :: errcode
character(len=*), optional, intent(out) :: errmsg
- std_name = 'Not Used!'
+ std_name = this%std_name
! Provide err values if requested:
if(present(errcode)) then
errcode = 0
end if
if(present(errmsg)) then
- errmsg = 'Still Not Used!'
+ errmsg = ''
end if
end subroutine ccp_get_standard_name
!------
+ subroutine ccp_set_standard_name(this, std_name, errcode, errmsg)
+ ! Set this constituent's standard name
+
+ ! Dummy arguments
+ class(ccpp_constituent_prop_ptr_t), intent(inout) :: this
+ character(len=*), intent(in) :: std_name
+ integer, optional, intent(out) :: errcode
+ character(len=*), optional, intent(out) :: errmsg
+
+ this%std_name = std_name
+
+ ! Provide err values if requested:
+ if(present(errcode)) then
+ errcode = 0
+ end if
+ if(present(errmsg)) then
+ errmsg = ''
+ end if
+
+ end subroutine ccp_set_standard_name
+
+ !------
+
subroutine ccp_is_thermo_active(this, val_out, errcode, errmsg)
! Dummy arguments
@@ -97,6 +127,29 @@ end subroutine ccp_is_water_species
!------
+ subroutine ccp_is_dry(this, val_out, errcode, errmsg)
+
+ ! Dummy arguments
+ class(ccpp_constituent_prop_ptr_t), intent(in) :: this
+ logical, intent(out) :: val_out
+ integer, optional, intent(out) :: errcode
+ character(len=*), optional, intent(out) :: errmsg
+
+ ! Pass back water species property:
+ val_out = this%species_is_dry
+
+ ! Provide err values if requested:
+ if(present(errcode)) then
+ errcode = 0
+ end if
+ if(present(errmsg)) then
+ errmsg = ''
+ end if
+
+ end subroutine ccp_is_dry
+
+ !------
+
subroutine ccp_set_thermo_active(this, thermo_flag, errcode, errmsg)
! Set whether this constituent is thermodynamically active, which
! means that certain physics schemes will use this constitutent
@@ -147,18 +200,41 @@ subroutine ccp_set_water_species(this, water_flag, errcode, errmsg)
end subroutine ccp_set_water_species
+ subroutine ccp_set_dry(this, dry_flag, errcode, errmsg)
+ ! Set whether this constituent is a dry species or not using the dry_flag which is passed in
+
+ ! Dummy arguments
+ class(ccpp_constituent_prop_ptr_t), intent(inout) :: this
+ logical, intent(in) :: dry_flag
+ integer, optional, intent(out) :: errcode
+ character(len=*), optional, intent(out) :: errmsg
+
+ ! Set dry_flag for this constituent:
+ this%species_is_dry = dry_flag
+
+ ! Provide err values if requested:
+ if(present(errcode)) then
+ errcode = 0
+ end if
+ if(present(errmsg)) then
+ errmsg = ''
+ end if
+
+ end subroutine ccp_set_dry
+
!+++++++++++++++++++++++++++++++++++++++++++++++++++++++
!CAM-equivalent CCPP constituents initialization routine
!+++++++++++++++++++++++++++++++++++++++++++++++++++++++
-subroutine ccpp_const_props_init()
+subroutine ccpp_const_props_init(ix_qv)
! Use statements:
- use constituents, only: pcnst
+ use constituents, only: pcnst, cnst_get_type_byind
use cam_abortutils, only: handle_allocate_error
use air_composition, only: dry_air_species_num
use air_composition, only: thermodynamic_active_species_idx
+ integer, intent(in) :: ix_qv
! Local variables:
integer :: ierr
integer :: m
@@ -185,6 +261,18 @@ subroutine ccpp_const_props_init()
end if
end do
+ ! Set "set_dry" property:
+ do m=1,pcnst
+ if (cnst_get_type_byind(m).eq.'dry') then
+ call ccpp_const_props(m)%set_dry(.true.)
+ else
+ call ccpp_const_props(m)%set_dry(.false.)
+ end if
+ end do
+
+ ! Set "std_name" property:
+ call ccpp_const_props(ix_qv)%set_standard_name('water_vapor_mixing_ratio_wrt_moist_air_and_condensed_water')
+
end subroutine ccpp_const_props_init
end module ccpp_constituent_prop_mod
diff --git a/src/utils/error_messages.F90 b/src/utils/error_messages.F90
deleted file mode 100644
index a2a64bca91..0000000000
--- a/src/utils/error_messages.F90
+++ /dev/null
@@ -1,151 +0,0 @@
-module error_messages
-
- !-----------------------------------------------------------------------
- !
- ! Purpose:
- ! General purpose routines for issuing error messages.
- !
- ! Author: B. Eaton
- !
- !-----------------------------------------------------------------------
- use cam_abortutils, only: endrun
- use cam_logfile, only: iulog
-
- implicit none
- save
- private
- public :: &
- alloc_err, &! Issue error message after non-zero return from an allocate statement.
- handle_err, &! Issue error message after non-zero return from anything
- handle_ncerr ! Handle error returns from netCDF library procedures.
-
- ! If an error message string is not empty, abort with that string as the
- ! error message.
- public :: handle_errmsg
-
-!##############################################################################
-contains
-!##############################################################################
-
- subroutine alloc_err( istat, routine, name, nelem )
-
- !-----------------------------------------------------------------------
- ! Purpose:
- ! Issue error message after non-zero return from an allocate statement.
- !
- ! Author: B. Eaton
- !-----------------------------------------------------------------------
-
- integer, intent(in) ::&
- istat ! status from allocate statement
- character(len=*), intent(in) ::&
- routine, &! routine that called allocate
- name ! name of array
- integer, intent(in) ::&
- nelem ! number of elements attempted to allocate
- !-----------------------------------------------------------------------
-
- if ( istat .ne. 0 ) then
- write(iulog,*)'ERROR trying to allocate memory in routine: ' &
- //trim(routine)
- write(iulog,*)' Variable name: '//trim(name)
- write(iulog,*)' Number of elements: ',nelem
- call endrun ('ALLOC_ERR')
- end if
-
- return
-
- end subroutine alloc_err
-
-!##############################################################################
-
- subroutine handle_err( istat, msg )
-
- !-----------------------------------------------------------------------
- ! Purpose:
- ! Issue error message after non-zero return from anything.
- !
- ! Author: T. Henderson
- !-----------------------------------------------------------------------
-
- integer, intent(in) :: istat ! status, zero = "no error"
- character(len=*), intent(in) :: msg ! error message to print
- !-----------------------------------------------------------------------
-
- if ( istat .ne. 0 ) then
- call endrun (trim(msg))
- end if
-
- return
-
- end subroutine handle_err
-
-!##############################################################################
-
- subroutine handle_ncerr( ret, mes, line )
-
- !-----------------------------------------------------------------------
- ! Purpose:
- ! Check netCDF library function return code. If error detected
- ! issue error message then abort.
- !
- ! Author: B. Eaton
- !-----------------------------------------------------------------------
-
-!-----------------------------------------------------------------------
- use netcdf
-!-----------------------------------------------------------------------
-
- integer, intent(in) ::&
- ret ! return code from netCDF library routine
- character(len=*), intent(in) ::&
- mes ! message to be printed if error detected
- integer, intent(in), optional :: line
- !-----------------------------------------------------------------------
-
- if ( ret .ne. NF90_NOERR ) then
- if(present(line)) then
- write(iulog,*) mes, line
- else
- write(iulog,*) mes
- end if
- write(iulog,*) nf90_strerror( ret )
- call endrun ('HANDLE_NCERR')
- endif
-
- return
-
- end subroutine handle_ncerr
-
-!##############################################################################
-
- subroutine handle_errmsg(errmsg, subname, extra_msg)
-
- ! String that is asserted to be null.
- character(len=*), intent(in) :: errmsg
- ! Name of procedure generating the message.
- character(len=*), intent(in), optional :: subname
- ! Additional message from the procedure calling this one.
- character(len=*), intent(in), optional :: extra_msg
-
- if (trim(errmsg) /= "") then
-
- if (present(extra_msg)) &
- write(iulog,*) "handle_errmsg: &
- &Message from caller: ",trim(extra_msg)
-
- if (present(subname)) then
- call endrun("ERROR: handle_errmsg: "// &
- trim(subname)//": "//trim(errmsg))
- else
- call endrun("ERROR: handle_errmsg: "// &
- "Error message received from routine: "//trim(errmsg))
- end if
-
- end if
-
- end subroutine handle_errmsg
-
-!##############################################################################
-
-end module error_messages
diff --git a/src/utils/namelist_utils.F90 b/src/utils/namelist_utils.F90
deleted file mode 100644
index c12dfad2d6..0000000000
--- a/src/utils/namelist_utils.F90
+++ /dev/null
@@ -1,6 +0,0 @@
-module namelist_utils
-
-use shr_nl_mod, only: &
- find_group_name => shr_nl_find_group_name
-
-end module namelist_utils