diff --git a/cime_config/tests.py b/cime_config/tests.py index bb93575208dc..c6c97ac9e1f0 100644 --- a/cime_config/tests.py +++ b/cime_config/tests.py @@ -91,7 +91,8 @@ "SMS_R_Ld5.ne4_ne4.FSCM-ARM97.eam-scm", "SMS_D_Ln5.ne4_oQU240.F2010", "SMS_Ln5.ne4pg2_oQU480.F2010", - "ERS_D.ne4_oQU240.F2010.eam-hommexx" + "ERS_D.ne4_oQU240.F2010.eam-hommexx", + "SMS_Ln9_P24x1.ne4_ne4.FDPSCREAM-ARM97", ) }, @@ -210,6 +211,7 @@ "SMS_D_Ld1.T62_oEC60to30v3.DTESTM", "SMS_D_Ld3.T62_oQU120.CMPASO-IAF", "SMS_D_Ld1.ne30pg2_r05_EC30to60E2r2.WCYCL1850", + "SMS_Ln5.ne30pg2_ne30pg2.F2010-SCREAM-LR-DYAMOND2", ) }, diff --git a/components/eam/bld/configure b/components/eam/bld/configure index b95e213439c0..5e3aa616a3ca 100755 --- a/components/eam/bld/configure +++ b/components/eam/bld/configure @@ -1752,6 +1752,19 @@ if (defined $opts{'cppdefs'}) { } $cfg_ref->set('cppdefs', $usr_cppdefs); +# Define is_scream_config flag to select which P3 version codes to use +# If the config uses both p3 and shoc, it is recognized as a scream case. + +my $is_scream_config = 0; +if ($microphys_pkg =~ /^p3/) { + if ($shoc_sgs) { + $is_scream_config = 1; + print "The case uses SCREAM config and SCREAM's version of P3 codes in Fortran 90.\n"; + } else { + print "The case use EAM's version of P3 codes in Fortran 90.\n"; + } +} + if ($usr_cppdefs and $print>=2) { print "Default and user CPP definitions: \'$usr_cppdefs\'$eol";} # The following CPP macro definitions are used to implement the compile-time options. They are @@ -2463,7 +2476,7 @@ my $fp_filename = 'Filepath'; # name of output filepath file my $cpp_filename = 'CIME_cppdefs'; # name of output file for eam's cppdefs in cime # Write the filepath file. -write_filepath("$cam_bld/$fp_filename", $cfg_ref); +write_filepath("$cam_bld/$fp_filename", $cfg_ref, $is_scream_config); if ($print) { print "creating $cam_bld/$fp_filename\n"; } if (($ccsm_seq)) { @@ -2507,7 +2520,7 @@ exit; sub write_filepath { - my ($file, $cfg_ref) = @_; + my ($file, $cfg_ref, $is_scream_config) = @_; my $fh = new IO::File; $fh->open(">$file") or die "** can't open filepath file: $file\n"; @@ -2701,6 +2714,13 @@ sub write_filepath print $fh "$camsrcdir/eam/src/physics/silhs\n"; } + # set Fortran P3 path, needed even if microphys is not P3 + if ($is_scream_config eq '1') { + print $fh "$camsrcdir/eam/src/physics/p3/scream\n"; + } else { + print $fh "$camsrcdir/eam/src/physics/p3/eam\n"; + } + print $fh "$camsrcdir/eam/src/dynamics/$dyn\n"; if($dyn eq 'se') { print $fh "$camsrcdir/homme/src/share\n"; diff --git a/components/eam/src/physics/p3/eam/micro_p3.F90 b/components/eam/src/physics/p3/eam/micro_p3.F90 new file mode 100644 index 000000000000..2e8c073404cf --- /dev/null +++ b/components/eam/src/physics/p3/eam/micro_p3.F90 @@ -0,0 +1,280 @@ +!__________________________________________________________________________________________ +! This module contains the Predicted Particle Property (P3) bulk microphysics scheme. ! +! ! +! This code was originally written by H. Morrison, MMM Division, NCAR (Dec 2012). ! +! Modifications were made by J. Milbrandt, RPN, Environment Canada (July 2014). ! +! Peter Caldwell (caldwell19@llnl.gov) further modified this code to remove multiple ! +! ice categories and to clean up/simplify the code for conversion to C++ (9/11/18) ! +! Jacob Shpund (jacob.shpund@pnnl.gov) implemented and further modified/cleaned ! +! the scheme in E3SMv2 ! +! ! +! Three configurations of the P3 scheme are currently available: ! +! 1) specified droplet number (i.e. 1-moment cloud water), 1 ice category ! +! 2) predicted droplet number (i.e. 2-moment cloud water), 1 ice category ! +! ! +! The 2-moment cloud version is based on a specified aerosol distribution and ! +! does not include a subgrid-scale vertical velocity for droplet activation. Hence, ! +! this version should only be used for high-resolution simulations that resolve ! +! vertical motion driving droplet activation. ! +! ! +! For details see: Morrison and Milbrandt (2015) [J. Atmos. Sci., 72, 287-311] ! +! Milbrandt and Morrison (2016) [J. Atmos. Sci., 73, 975-995] ! +! ! +! For questions or bug reports, please contact: ! +! Hugh Morrison (morrison@ucar.edu), or ! +! Jason Milbrandt (jason.milbrandt@canada.ca) ! +!__________________________________________________________________________________________! +! ! +! Version: 2.8.2.4 + Peter/Aaron's fixes ! +! Last updated: 2018-02-04 ! +!__________________________________________________________________________________________! +! Comments from Aaron Donahue: ! +! 1) Need to change the dz coordinate system in sedimentation to be consistent ! +! with E3SM's pressure based coordinate system, i.e. dp. ! +! 2) Move all physical constants into a micro_p3_util module and match them to ! +! universal constants in E3SM for consistency. ! +! 3) Need to include extra in/out values which correspond with microphysics PBUF ! +! variables and outputs expected in E3SM. ! +!__________________________________________________________________________________________! + +! Include bit-for-bit math macros. +! #include "bfb_math.inc" + +module micro_p3 + + ! get real kind from utils + use physics_utils, only: rtype,rtype8,btype + + use phys_control, only: use_hetfrz_classnuc + + ! physical and mathematical constants + use micro_p3_utils, only: rho_1000mb,rho_600mb,ar,br,f1r,f2r,rho_h2o,kr,kc,aimm,mi0,nccnst, & + eci,eri,bcn,cpw,cons1,cons3,cons4,cons5,cons6,cons7, & + inv_rho_h2o,inv_dropmass,qsmall,nsmall,cp,g,rd,rv,ep_2,inv_cp, & + thrd,sxth,piov6,rho_rimeMin, & + rho_rimeMax,inv_rho_rimeMax,max_total_ni,dbrk,nmltratio,clbfact_sub, & + clbfact_dep,iparam, isize, densize, rimsize, rcollsize, ice_table_size, collect_table_size, & + get_latent_heat, T_zerodegc, pi=>pi_e3sm, dnu, & + T_rainfrz, T_icenuc, T_homogfrz, iulog=>iulog_e3sm, & + masterproc=>masterproc_e3sm, calculate_incloud_mixingratios, mu_r_constant, & + lookup_table_1a_dum1_c, rho_h2o, & + do_Cooper_inP3 + + use wv_sat_scream, only:qv_sat + + ! Bit-for-bit math functions. +#ifdef SCREAM_CONFIG_IS_CMAKE + use physics_share_f2c, only: cxx_pow, cxx_sqrt, cxx_cbrt, cxx_gamma, cxx_log, & + cxx_log10, cxx_exp, cxx_expm1, cxx_tanh +#endif + + implicit none + save + + public :: p3_init,p3_main + + ! protected items should be treated as private for everyone except tests + + real(rtype), protected, dimension(densize,rimsize,isize,ice_table_size) :: ice_table_vals !ice lookup table values + + !ice lookup table values for ice-rain collision/collection + real(rtype), protected, dimension(densize,rimsize,isize,rcollsize,collect_table_size) :: collect_table_vals + + ! lookup table values for rain shape parameter mu_r + real(rtype), protected, dimension(150) :: mu_r_table_vals + + ! lookup table values for rain number- and mass-weighted fallspeeds and ventilation parameters + real(rtype), protected, dimension(300,10) :: vn_table_vals,vm_table_vals,revap_table_vals + +! + + + type realptr + real(rtype), dimension(:), pointer :: p + end type realptr + +contains + + SUBROUTINE p3_init(lookup_file_dir,version_p3) + !------------------------------------------------------------------------------------------! + ! This subroutine initializes all physical constants and parameters needed by the P3 ! + ! scheme, including reading in two lookup table files and creating a third. ! + ! 'P3_INIT' be called at the first model time step, prior to first call to 'P3_MAIN'. ! + !------------------------------------------------------------------------------------------! + + implicit none + + ! Passed arguments: + character*(*), intent(in) :: lookup_file_dir !directory of the lookup tables + character(len=16), intent(in) :: version_p3 !version number of P3 package + + END SUBROUTINE p3_init + + SUBROUTINE p3_main(qc,nc,qr,nr,th_atm,qv,dt,qi,qm,ni,bm, & + pres,dz,nc_nuceat_tend,nccn_prescribed,ni_activated,frzimm,frzcnt,frzdep,inv_qc_relvar,it,precip_liq_surf,precip_ice_surf,its,ite,kts,kte,diag_eff_radius_qc, & + diag_eff_radius_qi,rho_qi,do_predict_nc, do_prescribed_CCN,p3_autocon_coeff,p3_accret_coeff,p3_qc_autocon_expon,p3_nc_autocon_expon,p3_qc_accret_expon, & + p3_wbf_coeff,p3_max_mean_rain_size,p3_embryonic_rain_size, & + dpres,exner,qv2qi_depos_tend,precip_total_tend,nevapr,qr_evap_tend,precip_liq_flux,precip_ice_flux,rflx,sflx,cflx,cld_frac_r,cld_frac_l,cld_frac_i, & + p3_tend_out,mu_c,lamc,liq_ice_exchange,vap_liq_exchange, & + vap_ice_exchange,qv_prev,t_prev,col_location,diag_equiv_reflectivity,diag_ze_rain,diag_ze_ice & +#ifdef SCREAM_CONFIG_IS_CMAKE + ,elapsed_s & +#endif + ) + + !----------------------------------------------------------------------------------------! + ! ! + ! This is the main subroutine for the P3 microphysics scheme. It is called from the ! + ! wrapper subroutine ('MP_P3_WRAPPER') and is passed i,k slabs of all prognostic ! + ! variables -- hydrometeor fields, potential temperature, and water vapor mixing ratio. ! + ! Microphysical process rates are computed first. These tendencies are then used to ! + ! computed updated values of the prognostic variables. The hydrometeor variables are ! + ! then updated further due to sedimentation. ! + ! ! + ! Several diagnostic values are also computed and returned to the wrapper subroutine, ! + ! including precipitation rates. ! + ! ! + !----------------------------------------------------------------------------------------! + + implicit none + + !----- Input/ouput arguments: ----------------------------------------------------------! + + real(rtype), intent(inout), dimension(its:ite,kts:kte) :: qc ! cloud, mass mixing ratio kg kg-1 + ! note: Nc may be specified or predicted (set by do_predict_nc) + real(rtype), intent(inout), dimension(its:ite,kts:kte) :: nc ! cloud, number mixing ratio # kg-1 + real(rtype), intent(inout), dimension(its:ite,kts:kte) :: qr ! rain, mass mixing ratio kg kg-1 + real(rtype), intent(inout), dimension(its:ite,kts:kte) :: nr ! rain, number mixing ratio # kg-1 + real(rtype), intent(inout), dimension(its:ite,kts:kte) :: qi ! ice, total mass mixing ratio kg kg-1 + real(rtype), intent(inout), dimension(its:ite,kts:kte) :: qm ! ice, rime mass mixing ratio kg kg-1 + real(rtype), intent(inout), dimension(its:ite,kts:kte) :: ni ! ice, total number mixing ratio # kg-1 + real(rtype), intent(inout), dimension(its:ite,kts:kte) :: bm ! ice, rime volume mixing ratio m3 kg-1 + + real(rtype), intent(inout), dimension(its:ite,kts:kte) :: qv ! water vapor mixing ratio kg kg-1 + real(rtype), intent(inout), dimension(its:ite,kts:kte) :: th_atm ! potential temperature K + real(rtype), intent(in), dimension(its:ite,kts:kte) :: pres ! pressure Pa + real(rtype), intent(in), dimension(its:ite,kts:kte) :: dz ! vertical grid spacing m + real(rtype), intent(in), dimension(its:ite,kts:kte) :: nc_nuceat_tend ! IN ccn activated number tendency kg-1 s-1 + real(rtype), intent(in), dimension(its:ite,kts:kte) :: nccn_prescribed + real(rtype), intent(in), dimension(its:ite,kts:kte) :: ni_activated ! IN actived ice nuclei concentration 1/kg + real(rtype), intent(in), dimension(its:ite,kts:kte) :: frzimm,frzcnt,frzdep ! From macrophysics aerop (CNT scheme) [#/cm3] + real(rtype), intent(in) :: dt ! model time step s + + real(rtype), intent(out), dimension(its:ite) :: precip_liq_surf ! precipitation rate, liquid m s-1 + real(rtype), intent(out), dimension(its:ite) :: precip_ice_surf ! precipitation rate, solid m s-1 + real(rtype), intent(out), dimension(its:ite,kts:kte) :: diag_eff_radius_qc ! effective radius, cloud m + real(rtype), intent(out), dimension(its:ite,kts:kte) :: diag_eff_radius_qi ! effective radius, ice m + real(rtype), intent(out), dimension(its:ite,kts:kte) :: rho_qi ! bulk density of ice kg m-3 + real(rtype), intent(out), dimension(its:ite,kts:kte) :: mu_c ! Size distribution shape parameter for radiation + real(rtype), intent(out), dimension(its:ite,kts:kte) :: lamc ! Size distribution slope parameter for radiation + + integer, intent(in) :: its,ite ! array bounds (horizontal) + integer, intent(in) :: kts,kte ! array bounds (vertical) + integer, intent(in) :: it ! time step counter NOTE: starts at 1 for first time step + + logical(btype), intent(in) :: do_predict_nc ! .T. (.F.) for prediction (specification) of Nc + + real(rtype), intent(in), dimension(its:ite,kts:kte) :: dpres ! pressure thickness Pa + real(rtype), intent(in), dimension(its:ite,kts:kte) :: exner ! Exner expression + + ! OUTPUT for PBUF variables used by other parameterizations + real(rtype), intent(out), dimension(its:ite,kts:kte) :: qv2qi_depos_tend ! qitend due to deposition/sublimation + real(rtype), intent(out), dimension(its:ite,kts:kte) :: precip_total_tend ! Total precipitation (rain + snow) + real(rtype), intent(out), dimension(its:ite,kts:kte) :: nevapr ! evaporation of total precipitation (rain + snow) + real(rtype), intent(out), dimension(its:ite,kts:kte) :: qr_evap_tend ! evaporation of rain + real(rtype), intent(out), dimension(its:ite,kts:kte+1) :: precip_liq_flux ! grid-box average rain flux (kg m^-2 s^-1) pverp + real(rtype), intent(out), dimension(its:ite,kts:kte+1) :: precip_ice_flux ! grid-box average ice/snow flux (kg m^-2 s^-1) pverp + real(rtype), intent(out), dimension(its:ite,kts:kte+1) :: rflx ! grid-box average rain flux (kg m^-2 s^-1) pverp + real(rtype), intent(out), dimension(its:ite,kts:kte+1) :: sflx ! grid-box average ice/snow flux (kg m^-2 s^-1) pverp + real(rtype), intent(out), dimension(its:ite,kts:kte+1) :: cflx ! grid-box average cloud droplets flux (kg m^-2 s^-1) pverp + real(rtype), intent(out), dimension(its:ite,kts:kte) :: liq_ice_exchange ! sum of liq-ice phase change tendenices + real(rtype), intent(out), dimension(its:ite,kts:kte) :: vap_liq_exchange ! sum of vap-liq phase change tendenices + real(rtype), intent(out), dimension(its:ite,kts:kte) :: vap_ice_exchange ! sum of vap-ice phase change tendenices + + ! INPUT for prescribed CCN option + logical(btype), intent(in) :: do_prescribed_CCN + + ! INPUT for p3 tuning parameters + real(rtype), intent(in) :: p3_autocon_coeff ! autconversion coefficient + real(rtype), intent(in) :: p3_accret_coeff ! accretion coefficient + real(rtype), intent(in) :: p3_qc_autocon_expon ! autconversion qc exponent + real(rtype), intent(in) :: p3_nc_autocon_expon ! autconversion nc exponent + real(rtype), intent(in) :: p3_qc_accret_expon ! accretion qc and qr exponent + real(rtype), intent(in) :: p3_wbf_coeff ! WBF coefficient + real(rtype), intent(in) :: p3_max_mean_rain_size ! max mean rain size allowed + real(rtype), intent(in) :: p3_embryonic_rain_size ! embryonic rain size from autoconversion + + ! INPUT needed for PBUF variables used by other parameterizations + + real(rtype), intent(in), dimension(its:ite,kts:kte) :: cld_frac_i, cld_frac_l, cld_frac_r ! Ice, Liquid and Rain cloud fraction + real(rtype), intent(in), dimension(its:ite,kts:kte) :: qv_prev, t_prev ! qv and t from previous p3_main call + ! AaronDonahue, the following variable (p3_tend_out) is a catch-all for passing P3-specific variables outside of p3_main + ! so that they can be written as ouput. NOTE TO C++ PORT: This variable is entirely optional and doesn't need to be + ! included in the port to C++, or can be changed if desired. + real(rtype), intent(out), dimension(its:ite,kts:kte,49) :: p3_tend_out ! micro physics tendencies + real(rtype), intent(in), dimension(its:ite,3) :: col_location + real(rtype), intent(in), dimension(its:ite,kts:kte) :: inv_qc_relvar + real(rtype), intent(out), dimension(its:ite,kts:kte) :: diag_equiv_reflectivity,diag_ze_rain,diag_ze_ice ! equivalent reflectivity [dBZ] + +#ifdef SCREAM_CONFIG_IS_CMAKE + real(rtype), optional, intent(out) :: elapsed_s ! duration of main loop in seconds +#endif + + ! + !----- Local variables and parameters: -------------------------------------------------! + ! + + ! These outputs are no longer provided by p3_main. + real(rtype), dimension(its:ite,kts:kte) :: diag_vm_qi ! mass-weighted fall speed of ice m s-1 + real(rtype), dimension(its:ite,kts:kte) :: diag_diam_qi ! mean diameter of ice m + real(rtype), dimension(its:ite,kts:kte) :: pratot ! accretion of cloud by rain + real(rtype), dimension(its:ite,kts:kte) :: prctot ! autoconversion of cloud to rain + + real(rtype), dimension(its:ite,kts:kte) :: mu_r ! shape parameter of rain + real(rtype), dimension(its:ite,kts:kte) :: t_atm ! temperature at the beginning of the microhpysics step [K] + + ! 2D size distribution and fallspeed parameters: + + real(rtype), dimension(its:ite,kts:kte) :: lamr + real(rtype), dimension(its:ite,kts:kte) :: logn0r + + real(rtype), dimension(its:ite,kts:kte) :: nu + real(rtype), dimension(its:ite,kts:kte) :: cdist + real(rtype), dimension(its:ite,kts:kte) :: cdist1 + real(rtype), dimension(its:ite,kts:kte) :: cdistr + + ! Variables needed for in-cloud calculations + real(rtype), dimension(its:ite,kts:kte) :: inv_cld_frac_i, inv_cld_frac_l, inv_cld_frac_r ! Inverse cloud fractions (1/cld) + real(rtype), dimension(its:ite,kts:kte) :: qc_incld, qr_incld, qi_incld, qm_incld ! In cloud mass-mixing ratios + real(rtype), dimension(its:ite,kts:kte) :: nc_incld, nr_incld, ni_incld, bm_incld ! In cloud number concentrations + + real(rtype), dimension(its:ite,kts:kte) :: inv_dz,inv_rho,ze_ice,ze_rain,prec,rho, & + rhofacr,rhofaci,acn,latent_heat_sublim,latent_heat_vapor,latent_heat_fusion,qv_sat_l,qv_sat_i,qv_supersat_i, & + tmparr1,inv_exner + + ! -- scalar locals -- ! + + real(rtype) :: inv_dt, timeScaleFactor + + integer :: ktop,kbot,kdir,i + + logical(btype) :: is_nucleat_possible, is_hydromet_present + + !--These will be added as namelist parameters in the future + logical(btype), parameter :: debug_ON = .true. !.true. to switch on debugging checks/traps throughout code TODO: Turn this back off as default once the tlay error is found. + logical(btype), parameter :: debug_ABORT = .false. !.true. will result in forced abort in s/r 'check_values' + + real(rtype),dimension(its:ite,kts:kte) :: qc_old, nc_old, qr_old, nr_old, qi_old, ni_old, qv_old, th_atm_old + +#ifdef SCREAM_CONFIG_IS_CMAKE + integer :: clock_count1, clock_count_rate, clock_count_max, clock_count2, clock_count_diff +#endif + + END SUBROUTINE p3_main + +end module micro_p3 diff --git a/components/eam/src/physics/p3/eam/micro_p3_interface.F90 b/components/eam/src/physics/p3/eam/micro_p3_interface.F90 new file mode 100644 index 000000000000..b3beaf01248b --- /dev/null +++ b/components/eam/src/physics/p3/eam/micro_p3_interface.F90 @@ -0,0 +1,426 @@ +module micro_p3_interface + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + !! + !! Interface between E3SM and P3 microphysics + !! + !! Author: Peter Caldwell + !! + !! Last updated: 2018-09-12 + !! +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + + use shr_kind_mod, only: rtype=>shr_kind_r8 + use ppgrid, only: pcols,pver,pverp + +!comment: I think Kai added handle_errmsg. It would be better to +!use standard E3SM libraries if possible. + use error_messages, only: handle_errmsg + + use physics_types, only: physics_state, & + physics_ptend, & + physics_ptend_init + use physconst, only: mwdry, cpair, mwh2o, gravit, rair, cpliq, pi, & + rh2o, latvap, latice, tmelt, rhoh2o, rairv + use constituents, only: cnst_add, pcnst, sflxnam, apcnst, bpcnst, pcnst,& + cnst_name, cnst_get_ind,cnst_longname + use physics_buffer, only: physics_buffer_desc, dtype_r8, & + pbuf_get_field, pbuf_add_field,dyn_time_lvls,dtype_i4, & + pbuf_set_field, pbuf_get_index, & + pbuf_old_tim_idx + use ref_pres, only: top_lev=>trop_cloud_top_lev + use phys_control, only: phys_getopts,use_hetfrz_classnuc + use cam_abortutils, only: endrun + use spmd_utils, only: masterproc + use cam_logfile, only: iulog + use time_manager, only: is_first_step, get_curr_date + use perf_mod, only: t_startf, t_stopf + use micro_p3_utils, only: do_Cooper_inP3 + use pio, only: file_desc_t, pio_nowrite + use cam_pio_utils, only: cam_pio_openfile,cam_pio_closefile + use cam_grid_support, only: cam_grid_check, cam_grid_id, cam_grid_get_dim_names + use ncdio_atm, only: infld + use ppgrid, only: begchunk, endchunk, pcols, pver, pverp,psubcols + use cam_history_support, only: add_hist_coord + + implicit none + save + + public :: micro_p3_init, micro_p3_register, micro_p3_tend, & + micro_p3_init_cnst, micro_p3_implements_cnst, & + micro_p3_readnl + + + character(len=16), parameter :: unset_str = 'UNSET' + + private + + !Define indices for state%q constituents at module level so + !defining them in micro_p3_register makes them permanently + !available. + CHARACTER(len=16) :: precip_frac_method = 'max_overlap' ! AaronDonahue, Hard-coded for now, should be fixed in the future + + integer, public :: & + ixcldliq = -1, & ! cloud liquid amount index + ixcldice = -1, & ! ice index + ixnumliq = -1, & ! cloud liquid number index + ixnumice = -1, & ! cloud ice number index + ixrain = -1, & ! rain index + ixnumrain= -1, & ! rain number index + ixcldrim = -1, & ! rime index ?? + ixrimvol = -1, & ! rime volume index ?? + ixqm = -1 ! ?? index ?? + +!! pbuf + integer :: & + cldo_idx, & + qme_idx, & + precip_total_tend_idx, & + nevapr_idx, & + dei_idx, & + rate1_cw2pr_st_idx, & + mu_idx, & + lambdac_idx, & + rei_idx, & + rel_idx, & + ls_flxprc_idx, & + ls_flxsnw_idx, & + ls_reffrain_idx, & + ls_reffsnow_idx, & + cv_reffliq_idx, & + cv_reffice_idx, & + qr_evap_tend_idx, & + cmeliq_idx, & + relvar_idx, & + qv_prev_idx, & + t_prev_idx, & + accre_enhan_idx, & + mon_ccn_1_idx, & + mon_ccn_2_idx, & + current_month !Needed for prescribed CCN option + +! Physics buffer indices for fields registered by other modules + integer :: & + ast_idx = -1 + + integer :: & + ni_activated_idx = -1, & + npccn_idx = -1, & + prec_str_idx = -1, & + prec_pcw_idx = -1, & + prec_sed_idx = -1, & + snow_str_idx = -1, & + snow_pcw_idx = -1, & + snow_sed_idx = -1 + + integer :: & + frzimm_idx = -1, & + frzcnt_idx = -1, & + frzdep_idx = -1 + + real(rtype) :: & + micro_mg_accre_enhan_fac = huge(1.0_rtype), & !Accretion enhancement factor from namelist + prc_coef1_in = huge(1.0_rtype), & + prc_exp_in = huge(1.0_rtype), & + prc_exp1_in = huge(1.0_rtype), & + p3_autocon_coeff = huge(1.0_rtype), & + p3_accret_coeff = huge(1.0_rtype), & + p3_qc_autocon_expon = huge(1.0_rtype), & + p3_nc_autocon_expon = huge(1.0_rtype), & + p3_qc_accret_expon = huge(1.0_rtype), & + p3_wbf_coeff = huge(1.0_rtype), & + p3_max_mean_rain_size = huge(1.0_rtype), & + p3_embryonic_rain_size = huge(1.0_rtype) + + + integer :: ncnst + + character(len=8), parameter :: & ! Constituent names + cnst_names(8) = (/'CLDLIQ', 'CLDICE','NUMLIQ','NUMICE', & + 'RAINQM', 'CLDRIM','NUMRAI','BVRIM '/) + + character(len=128) :: micro_p3_lookup_dir = unset_str ! location of p3 input files + character(len=16) :: micro_p3_tableversion = unset_str ! P3 table version + logical :: micro_aerosolactivation = .false. ! Use aerosol activation + logical :: micro_subgrid_cloud = .false. ! Use subgrid cloudiness + logical :: micro_tend_output = .false. ! Default microphysics tendencies to output file + logical :: do_prescribed_CCN = .false. ! Use prescribed CCN + + contains +!=============================================================================== +subroutine micro_p3_readnl(nlfile) + + use namelist_utils, only: find_group_name + use units, only: getunit, freeunit + use mpishorthand + + character(len=*), intent(in) :: nlfile ! filepath for file containing namelist input + + ! Local variables + integer :: unitn, ierr + character(len=*), parameter :: subname = 'micro_p3_cam_readnl' + + +end subroutine micro_p3_readnl + !================================================================================================ + + subroutine micro_p3_register() + + logical :: prog_modal_aero ! prognostic aerosols + + end subroutine micro_p3_register + + !================================================================================================ + function micro_p3_implements_cnst(name) + + ! Return true if specified constituent is implemented by the + ! microphysics package + + character(len=*), intent(in) :: name ! constituent name + logical :: micro_p3_implements_cnst ! return value + + + end function micro_p3_implements_cnst + + + !================================================================================================ + + subroutine micro_p3_init_cnst(name, q) + + ! Initialize the microphysics constituents, if they are + ! not read from the initial file. + + character(len=*), intent(in) :: name ! constituent name + real(rtype), intent(out) :: q(:,:) ! mass mixing ratio (gcol, plev) + + + end subroutine micro_p3_init_cnst + + !================================================================================================ + + subroutine get_prescribed_CCN_from_file(micro_p3_lookup_dir,month_int, ccn_values) + + character*(*), intent(in) :: micro_p3_lookup_dir !directory of the lookup tables + real(rtype), intent(inout) :: ccn_values(pcols,pver,begchunk:endchunk) + integer, intent(in) :: month_int + + !internal variables + character(len=100) :: base_file_name + character(len=500) :: filename + character(len=20) :: dim1name, dim2name + character(len=20) :: mon_str + type(file_desc_t) :: nccn_ncid + integer :: year, month, day, tod, next_month, grid_id + logical :: found = .false. + + + end subroutine get_prescribed_CCN_from_file + + !================================================================================================ + + subroutine micro_p3_init(pbuf2d) + use micro_p3, only: p3_init + use cam_history, only: addfld, add_default, horiz_only + use micro_p3_utils, only: micro_p3_utils_init + + type(physics_buffer_desc), pointer :: pbuf2d(:,:) + integer :: m, mm + integer :: ierr + logical :: history_amwg ! output the variables used by the AMWG diag package + logical :: history_verbose ! produce verbose history output + logical :: history_budget ! Output tendencies and state variables for CAM4 + integer :: budget_histfile ! output history file number for budget fields + ! temperature, water vapor, cloud ice and cloud + + !needed for prescribed CCN option: + character(len=20) :: base_file_name + character(len=500) :: filename, filename_next_month + character(len=20) :: dim1name, dim2name + type(file_desc_t) :: nccn_ncid + integer :: year, month, day, tod, next_month, grid_id + logical :: found = .false. + real(rtype), pointer :: ccn_values(:,:,:) + + + end subroutine micro_p3_init + + !================================================================================================ + subroutine get_cloud_fraction(its,ite,kts,kte,ast,qc,qr,qi,method, & + cld_frac_i,cld_frac_l,cld_frac_r) + + use micro_p3_utils, only: mincld, qsmall + + integer,intent(in) :: its,ite,kts,kte + real(rtype),dimension(its:ite,kts:kte),intent(in) :: ast, qc, qr, qi + character(len=16),intent(in) :: method + real(rtype),dimension(its:ite,kts:kte),intent(out) :: cld_frac_i, cld_frac_l, cld_frac_r + real(rtype),dimension(its:ite,kts:kte) :: cldm + + integer :: i,k + integer :: ktop, kbot, kdir + + end subroutine get_cloud_fraction + + !================================================================================================ + + subroutine get_prescribed_CCN(nccn_prescribed,micro_p3_lookup_dir,its,ite,kts,kte,pbuf,lchnk) + + !INOUT/OUTPUT VARIABLES + integer,intent(in) :: its,ite,kts,kte,lchnk + real(rtype),dimension(its:ite,kts:kte),intent(inout) :: nccn_prescribed + character*(*), intent(in) :: micro_p3_lookup_dir !directory of the lookup tables + type(physics_buffer_desc), pointer :: pbuf(:) + + + + end subroutine get_prescribed_CCN + + !================================================================================================ + subroutine micro_p3_tend(state, ptend, dtime, pbuf) + + use phys_grid, only: get_rlat_all_p, get_rlon_all_p, get_gcol_all_p + use time_manager, only: get_nstep + use cam_history, only: outfld + use time_manager, only: get_nstep + use micro_p3, only: p3_main + use micro_p3_utils, only: avg_diameter, & + rho_h2o, & + rho_h2os, & + qsmall, & + mincld, & + inv_cp + + !INPUT/OUTPUT VARIABLES + type(physics_state), intent(in) :: state + type(physics_ptend), intent(out) :: ptend + real(rtype), intent(in) :: dtime + type(physics_buffer_desc), pointer :: pbuf(:) + logical :: lq(pcnst) !list of what constituents to update + + !INTERNAL VARIABLES + real(rtype) :: dz(pcols,pver) !geometric layer thickness m + real(rtype) :: cldliq(pcols,pver) !cloud liquid water mixing ratio kg/kg + real(rtype) :: numliq(pcols,pver) !cloud liquid water drop concentraiton #/kg + real(rtype) :: rain(pcols,pver) !rain water mixing ratio kg/kg + real(rtype) :: numrain(pcols,pver) !rain water number concentration #/kg + real(rtype) :: qv(pcols,pver) !water vapor mixing ratio kg/kg + real(rtype) :: ice(pcols,pver) !total ice water mixing ratio kg/kg + real(rtype) :: qm(pcols,pver) !rime ice mixing ratio kg/kg + real(rtype) :: numice(pcols,pver) !total ice crystal number concentration #/kg + real(rtype) :: rimvol(pcols,pver) !rime volume mixing ratio m3/kg + real(rtype) :: temp(pcols,pver) !temperature copy needed for tendency K + real(rtype) :: th(pcols,pver) !potential temperature K + real(rtype) :: precip_liq_surf(pcols) !precipitation rate, liquid m s-1 + real(rtype) :: precip_ice_surf(pcols) !precipitation rate, solid m s-1 + + real(rtype) :: rho_qi(pcols,pver) !bulk density of ice kg m-1 + real(rtype) :: pres(pcols,pver) !pressure at midlevel hPa + real(rtype) :: qv2qi_depos_tend(pcols,pver) + real(rtype) :: precip_liq_flux(pcols,pver+1) !grid-box average rain flux (kg m^-2s^-1) pverp + real(rtype) :: precip_ice_flux(pcols,pver+1) !grid-box average ice/snow flux (kg m^-2s^-1) pverp + real(rtype) :: rflx(pcols,pver+1) !grid-box average rain flux (kg m^-2s^-1) pverp + real(rtype) :: sflx(pcols,pver+1) !grid-box average ice/snow flux (kg m^-2s^-1) pverp + real(rtype) :: cflx(pcols,pver+1) !grid-box average cloud flux (kg m^-2s^-1) pverp + real(rtype) :: exner(pcols,pver) !exner formula for converting between potential and normal temp + real(rtype) :: cld_frac_r(pcols,pver) !rain cloud fraction + real(rtype) :: cld_frac_l(pcols,pver) !liquid cloud fraction + real(rtype) :: cld_frac_i(pcols,pver) !ice cloud fraction + real(rtype) :: tend_out(pcols,pver,49) !microphysical tendencies + real(rtype), dimension(pcols,pver) :: liq_ice_exchange ! sum of liq-ice phase change tendenices + real(rtype), dimension(pcols,pver) :: vap_liq_exchange ! sum of vap-liq phase change tendenices + real(rtype), dimension(pcols,pver) :: vap_ice_exchange ! sum of vap-ice phase change tendenices + real(rtype), dimension(pcols,pver) :: diag_equiv_reflectivity,diag_ze_rain,diag_ze_ice ! equivalent reflectivity [dBz] + + !Prescribed CCN concentration + real(rtype), dimension(pcols,pver) :: nccn_prescribed + + ! PBUF Variables + real(rtype), pointer :: ast(:,:) ! Relative humidity cloud fraction + real(rtype), pointer :: ni_activated(:,:) ! ice nucleation number + real(rtype), pointer :: npccn(:,:) ! liquid activation number tendency + real(rtype), pointer :: cmeliq(:,:) + !! + real(rtype), pointer :: prec_str(:) ! [Total] Sfc flux of precip from stratiform [ m/s ] + real(rtype), pointer :: prec_sed(:) ! Surface flux of total cloud water from sedimentation + real(rtype), pointer :: prec_pcw(:) ! Sfc flux of precip from microphysics [ m/s ] + real(rtype), pointer :: snow_str(:) ! [Total] Sfc flux of snow from stratiform [ m/s ] + real(rtype), pointer :: snow_pcw(:) ! Sfc flux of snow from microphysics [ m/s ] + real(rtype), pointer :: snow_sed(:) ! Surface flux of cloud ice from sedimentation + real(rtype), pointer :: relvar(:,:) ! cloud liquid relative variance [-] + real(rtype), pointer :: cldo(:,:) ! Old cloud fraction + real(rtype), pointer :: qr_evap_tend(:,:) ! precipitation evaporation rate + real(rtype), pointer :: qv_prev(:,:) ! qv from previous p3_main call + real(rtype), pointer :: t_prev(:,:) ! t from previous p3_main call + !! wetdep + real(rtype), pointer :: qme(:,:) + real(rtype), pointer :: precip_total_tend(:,:) ! Total precipitation (rain + snow) + real(rtype), pointer :: nevapr(:,:) ! Evaporation of total precipitation (rain + snow) + !! COSP simulator + real(rtype), pointer :: rel(:,:) ! Liquid effective drop radius (microns) + real(rtype), pointer :: rei(:,:) ! Ice effective drop size (microns) + real(rtype), pointer :: flxprc(:,:) ! P3 grid-box mean flux_large_scale_cloud_rain+snow at interfaces (kg/m2/s) + real(rtype), pointer :: flxsnw(:,:) ! P3 grid-box mean flux_large_scale_cloud_snow at interfaces (kg/m2/s) + real(rtype), pointer :: reffrain(:,:) ! P3 diagnostic rain effective radius (um) + real(rtype), pointer :: reffsnow(:,:) ! P3 diagnostic snow effective radius (um) + real(rtype), pointer :: cvreffliq(:,:) ! convective cloud liquid effective radius (um) + real(rtype), pointer :: cvreffice(:,:) ! convective cloud ice effective radius (um) + !! radiation + real(rtype), pointer :: dei(:,:) ! Ice effective diameter (um) + real(rtype), pointer :: mu(:,:) ! Size distribution shape parameter for radiation + real(rtype), pointer :: lambdac(:,:) ! Size distribution slope parameter for radiation + ! DONE PBUF + ! For recording inputs/outputs to p3_main + real(rtype) :: p3_main_inputs(pcols,pver+1,17) ! Record of inputs for p3_main + real(rtype) :: p3_main_outputs(pcols,pver+1,31) ! Record of outputs for p3_main + + ! Derived Variables + real(rtype) :: icimrst(pcols,pver) ! stratus ice mixing ratio - on grid + real(rtype) :: icwmrst(pcols,pver) ! stratus water mixing ratio - on grid + real(rtype) :: rho(pcols,pver) + real(rtype) :: drout2(pcols,pver) + real(rtype) :: reff_rain(pcols,pver) + real(rtype) :: col_location(pcols,3),tmp_loc(pcols) ! Array of column lon (index 1) and lat (index 2) + integer :: tmpi_loc(pcols) ! Global column index temp array + + ! Variables used for microphysics output + real(rtype) :: aqrain(pcols,pver) + real(rtype) :: anrain(pcols,pver) + real(rtype) :: nfice(pcols,pver) + real(rtype) :: efcout(pcols,pver) + real(rtype) :: efiout(pcols,pver) + real(rtype) :: ncout(pcols,pver) + real(rtype) :: niout(pcols,pver) + real(rtype) :: freqr(pcols,pver) + real(rtype) :: freql(pcols,pver) + real(rtype) :: freqi(pcols,pver) + real(rtype) :: cdnumc(pcols) + real(rtype) :: icinc(pcols,pver) + real(rtype) :: icwnc(pcols,pver) + + ! variables for the CNT primary / heterogeneous freezing + real(rtype), pointer :: frzimm(:,:) + real(rtype), pointer :: frzcnt(:,:) + real(rtype), pointer :: frzdep(:,:) + real(rtype) :: frzimm_in(pcols,pver) + real(rtype) :: frzcnt_in(pcols,pver) + real(rtype) :: frzdep_in(pcols,pver) + + integer :: it !timestep counter - + integer :: its, ite !horizontal bounds (column start,finish) + integer :: kts !closest level to TOM - + integer :: kte !near surface level - + + logical :: do_predict_nc !prognostic droplet concentration or not? + logical :: do_subgrid_clouds !use subgrid cloudiness in tendency calculations? + integer :: icol, ncol, k + integer :: psetcols, lchnk + integer :: itim_old + + ! For rrtmg optics. specified distribution. + real(rtype), parameter :: dcon = 25.e-6_rtype ! Convective size distribution effective radius (um) + real(rtype), parameter :: mucon = 5.3_rtype ! Convective size distribution shape parameter + real(rtype), parameter :: deicon = 50._rtype ! Convective ice effective diameter (um) + + end subroutine micro_p3_tend + +end module micro_p3_interface diff --git a/components/eam/src/physics/p3/eam/micro_p3_utils.F90 b/components/eam/src/physics/p3/eam/micro_p3_utils.F90 new file mode 100644 index 000000000000..45a6920e196d --- /dev/null +++ b/components/eam/src/physics/p3/eam/micro_p3_utils.F90 @@ -0,0 +1,147 @@ +module micro_p3_utils + + use physics_utils, only: rtype, rtype8, itype, btype + + implicit none + private + save + + public :: get_latent_heat, micro_p3_utils_init, & + avg_diameter, calculate_incloud_mixingratios + + integer, public :: iulog_e3sm + logical(btype), public :: masterproc_e3sm + + ! Signaling NaN bit pattern that represents a limiter that's turned off. + integer(itype), parameter :: limiter_off = int(Z'7FF1111111111111', itype) + + real(rtype), public, parameter :: qsmall = 1.e-14_rtype + real(rtype), public, parameter :: nsmall = 1.e-16_rtype + + real(rtype) :: latent_heat_vapor, latent_heat_sublim, latent_heat_fusion + + real(rtype),public :: rho_1000mb,rho_600mb,ar,br,f1r,f2r,ecr,rho_h2o,kr,kc,aimm,bimm,rin,mi0,nccnst, & + eci,eri,bcn,cpw,cons1,cons2,cons3,cons4,cons5,cons6,cons7, & + inv_rho_h2o,inv_dropmass,cp,g,rd,rv,ep_2,inv_cp, & + thrd,sxth,piov3,piov6,rho_rimeMin, & + rho_rimeMax,inv_rho_rimeMax,max_total_ni,dbrk,nmltratio,clbfact_sub, & + clbfact_dep + + logical,public :: do_Cooper_inP3 ! Use prescribed CCN + + real(rtype),dimension(16), public :: dnu + + real(rtype), public, parameter :: mu_r_constant = 0.0_rtype !1.0_rtype + real(rtype), public, parameter :: lookup_table_1a_dum1_c = 4.135985029041767d+00 ! 1.0/(0.1*log10(261.7)) + + real(rtype),public :: T_zerodegc ! Temperature at zero degree celcius ~K + real(rtype),public :: T_rainfrz ! Contact and immersion freexing temp, -4C ~K + real(rtype),public :: T_homogfrz ! Homogeneous freezing temperature, -40C ~K + real(rtype),public :: T_icenuc ! Ice nucleation temperature, -5C ~K + + real(rtype),public :: pi_e3sm + ! ice microphysics lookup table array dimensions + integer, public,parameter :: isize = 50 + integer, public,parameter :: densize = 5 + integer, public,parameter :: rimsize = 4 + integer, public,parameter :: rcollsize = 30 + integer, public,parameter :: ice_table_size = 12 ! number of quantities used from lookup table + integer, public,parameter :: collect_table_size = 2 ! number of ice-rain collection quantities used from lookup table + ! switch for warm-rain parameterization + ! = 1 Seifert and Beheng 2001 + ! = 2 Beheng 1994 + ! = 3 Khairoutdinov and Kogan 2000 + integer, public,parameter :: iparam = 3 + + real(rtype), parameter, public :: mincld=0.0001_rtype + real(rtype), parameter, public :: rho_h2os = 917._rtype ! bulk density water solid + real(rtype), parameter, public :: dropmass = 5.2e-7_rtype + + ! particle mass-diameter relationship + ! currently we assume spherical particles for cloud ice/snow + ! m = cD^d + ! exponent + real(rtype), parameter :: dsph = 3._rtype + + ! Bounds for mean diameter for different constituents. + ! real(rtype), parameter :: lam_bnd_rain(2) = 1._rtype/[500.e-6_rtype, 20.e-6_rtype] + ! real(rtype), parameter :: lam_bnd_snow(2) = 1._rtype/[2000.e-6_rtype, 10.e-6_rtype] + + ! Minimum average mass of particles. + real(rtype), parameter :: min_mean_mass_liq = 1.e-20_rtype + real(rtype), parameter :: min_mean_mass_ice = 1.e-20_rtype + + ! in-cloud values + REAL(rtype), PARAMETER :: min_cld_frac = 1.e-20_rtype !! threshold min value for cloud fraction + real(rtype), parameter :: incloud_limit = 5.1E-3 + real(rtype), parameter :: precip_limit = 1.0E-2 + + contains +!__________________________________________________________________________________________! +! ! +!__________________________________________________________________________________________! + subroutine micro_p3_utils_init(cpair,rair,rh2o,rhoh2o,mwh2o,mwdry,gravit,latvap,latice, & + cpliq,tmelt,pi,iulog,masterproc) + + real(rtype), intent(in) :: cpair + real(rtype), intent(in) :: rair + real(rtype), intent(in) :: rh2o + real(rtype), intent(in) :: rhoh2o + real(rtype), intent(in) :: mwh2o + real(rtype), intent(in) :: mwdry + real(rtype), intent(in) :: gravit + real(rtype), intent(in) :: latvap + real(rtype), intent(in) :: latice + real(rtype), intent(in) :: cpliq + real(rtype), intent(in) :: tmelt + real(rtype), intent(in) :: pi + integer, intent(in) :: iulog + logical(btype), intent(in) :: masterproc + + + return + end subroutine micro_p3_utils_init +!__________________________________________________________________________________________! +! ! +!__________________________________________________________________________________________! + subroutine get_latent_heat(its,ite,kts,kte,v,s,f) + + integer,intent(in) :: its,ite,kts,kte + real(rtype),dimension(its:ite,kts:kte),intent(out) :: v,s,f + +! integer i,k + + end subroutine get_latent_heat + +!__________________________________________________________________________________________! +! ! +!__________________________________________________________________________________________! + subroutine calculate_incloud_mixingratios(qc,qr,qi,qm,nc,nr,ni,bm, & + inv_cld_frac_l,inv_cld_frac_i,inv_cld_frac_r, & + qc_incld,qr_incld,qi_incld,qm_incld,nc_incld,nr_incld,ni_incld,bm_incld) + + real(rtype),intent(in) :: qc, qr, qi, qm + real(rtype),intent(in) :: nc, nr, ni, bm + real(rtype),intent(in) :: inv_cld_frac_l, inv_cld_frac_i, inv_cld_frac_r + real(rtype),intent(out) :: qc_incld, qr_incld, qi_incld, qm_incld + real(rtype),intent(out) :: nc_incld, nr_incld, ni_incld, bm_incld + + end subroutine calculate_incloud_mixingratios +!__________________________________________________________________________________________! +! ! +!__________________________________________________________________________________________! + +real(rtype) elemental function avg_diameter(q, n, rho_air, rho_sub) + ! Finds the average diameter of particles given their density, and + ! mass/number concentrations in the air. + ! Assumes that diameter follows an exponential distribution. + real(rtype), intent(in) :: q ! mass mixing ratio + real(rtype), intent(in) :: n ! number concentration (per volume) + real(rtype), intent(in) :: rho_air ! local density of the air + real(rtype), intent(in) :: rho_sub ! density of the particle substance + + +end function avg_diameter + + +end module micro_p3_utils diff --git a/components/eam/src/physics/cam/micro_p3.F90 b/components/eam/src/physics/p3/scream/micro_p3.F90 similarity index 100% rename from components/eam/src/physics/cam/micro_p3.F90 rename to components/eam/src/physics/p3/scream/micro_p3.F90 diff --git a/components/eam/src/physics/cam/micro_p3_interface.F90 b/components/eam/src/physics/p3/scream/micro_p3_interface.F90 similarity index 100% rename from components/eam/src/physics/cam/micro_p3_interface.F90 rename to components/eam/src/physics/p3/scream/micro_p3_interface.F90 diff --git a/components/eam/src/physics/cam/micro_p3_utils.F90 b/components/eam/src/physics/p3/scream/micro_p3_utils.F90 similarity index 100% rename from components/eam/src/physics/cam/micro_p3_utils.F90 rename to components/eam/src/physics/p3/scream/micro_p3_utils.F90 diff --git a/components/eamxx/src/physics/p3/CMakeLists.txt b/components/eamxx/src/physics/p3/CMakeLists.txt index 23b3eb558e2d..7719cbfdf156 100644 --- a/components/eamxx/src/physics/p3/CMakeLists.txt +++ b/components/eamxx/src/physics/p3/CMakeLists.txt @@ -2,7 +2,7 @@ set(P3_SRCS p3_f90.cpp p3_ic_cases.cpp p3_iso_c.f90 - ${SCREAM_BASE_DIR}/../eam/src/physics/cam/micro_p3.F90 + ${SCREAM_BASE_DIR}/../eam/src/physics/p3/scream/micro_p3.F90 atmosphere_microphysics.cpp atmosphere_microphysics_run.cpp )