module rk_stratiform !------------------------------------------------------------------------------------------------------- ! ! Provides the CAM interface to the Rasch and Kristjansson (RK) ! prognostic cloud microphysics, and the cam3/4 macrophysics. ! !------------------------------------------------------------------------------------------------------- use shr_kind_mod, only: r8=>shr_kind_r8 use ppgrid, only: pcols, pver, pverp use physconst, only: gravit, latvap, latice use phys_control, only: phys_getopts use constituents, only: pcnst use spmd_utils, only: masterproc use cam_logfile, only: iulog use cam_abortutils, only: endrun use perf_mod implicit none private save public :: rk_stratiform_register, rk_stratiform_init_cnst, rk_stratiform_implements_cnst public :: rk_stratiform_init public :: rk_stratiform_tend public :: rk_stratiform_readnl ! Physics buffer indices integer :: landm_idx = 0 integer :: qcwat_idx = 0 integer :: lcwat_idx = 0 integer :: tcwat_idx = 0 integer :: cld_idx = 0 integer :: ast_idx = 0 integer :: concld_idx = 0 integer :: fice_idx = 0 integer :: qme_idx = 0 integer :: prain_idx = 0 integer :: nevapr_idx = 0 integer :: wsedl_idx = 0 integer :: rei_idx = 0 integer :: rel_idx = 0 integer :: shfrc_idx = 0 integer :: cmfmc_sh_idx = 0 integer :: prec_str_idx = 0 integer :: snow_str_idx = 0 integer :: prec_sed_idx = 0 integer :: snow_sed_idx = 0 integer :: prec_pcw_idx = 0 integer :: snow_pcw_idx = 0 integer :: ls_flxprc_idx = 0 integer :: ls_flxsnw_idx = 0 integer, parameter :: ncnst = 2 ! Number of constituents character(len=8), dimension(ncnst), parameter & ! Constituent names :: cnst_names = (/'CLDLIQ', 'CLDICE'/) logical :: use_shfrc ! Local copy of flag from convect_shallow_use_shfrc logical :: do_cnst = .false. ! True when this module has registered constituents. integer :: & ixcldliq, &! cloud liquid amount index ixcldice ! cloud ice amount index real(r8), parameter :: unset_r8 = huge(1.0_r8) logical :: do_psrhmin !=============================================================================== contains !=============================================================================== subroutine rk_stratiform_readnl(nlfile) use namelist_utils, only: find_group_name use units, only: getunit, freeunit use cldwat, only: cldwat_init use mpishorthand character(len=*), intent(in) :: nlfile ! filepath for file containing namelist input ! Local variables integer :: unitn, ierr character(len=*), parameter :: subname = 'rk_stratiform_readnl' ! Namelist variables real(r8) :: rk_strat_icritw = unset_r8 ! icritw = threshold for autoconversion of warm ice real(r8) :: rk_strat_icritc = unset_r8 ! icritc = threshold for autoconversion of cold ice real(r8) :: rk_strat_conke = unset_r8 ! conke = tunable constant for evaporation of precip real(r8) :: rk_strat_r3lcrit = unset_r8 ! r3lcrit = critical radius where liq conversion begins real(r8) :: rk_strat_polstrat_rhmin = unset_r8 ! condensation threadhold in polar stratosphere namelist /rk_stratiform_nl/ rk_strat_icritw, rk_strat_icritc, rk_strat_conke, rk_strat_r3lcrit namelist /rk_stratiform_nl/ rk_strat_polstrat_rhmin !----------------------------------------------------------------------------- if (masterproc) then unitn = getunit() open( unitn, file=trim(nlfile), status='old' ) call find_group_name(unitn, 'rk_stratiform_nl', status=ierr) if (ierr == 0) then read(unitn, rk_stratiform_nl, iostat=ierr) if (ierr /= 0) then call endrun(subname // ':: ERROR reading namelist') end if end if close(unitn) call freeunit(unitn) end if #ifdef SPMD ! Broadcast namelist variables call mpibcast(rk_strat_icritw, 1, mpir8, 0, mpicom) call mpibcast(rk_strat_icritc, 1, mpir8, 0, mpicom) call mpibcast(rk_strat_conke, 1, mpir8, 0, mpicom) call mpibcast(rk_strat_r3lcrit, 1, mpir8, 0, mpicom) call mpibcast(rk_strat_polstrat_rhmin, 1, mpir8, 0, mpicom) #endif do_psrhmin = rk_strat_polstrat_rhmin .ne. unset_r8 call cldwat_init(rk_strat_icritw,rk_strat_icritc,rk_strat_conke,& rk_strat_r3lcrit,rk_strat_polstrat_rhmin,do_psrhmin) end subroutine rk_stratiform_readnl subroutine rk_stratiform_register !---------------------------------------------------------------------- ! ! ! ! Register the constituents (cloud liquid and cloud ice) and the fields ! ! in the physics buffer. ! ! ! !---------------------------------------------------------------------- ! use constituents, only: cnst_add, pcnst use physconst, only: mwh2o, cpair use physics_buffer, only : pbuf_add_field, dtype_r8, dyn_time_lvls !----------------------------------------------------------------------- ! Take note of the fact that we are registering constituents. do_cnst = .true. ! Register cloud water and save indices. call cnst_add(cnst_names(1), mwh2o, cpair, 0._r8, ixcldliq, & longname='Grid box averaged cloud liquid amount', is_convtran1=.true.) call cnst_add(cnst_names(2), mwh2o, cpair, 0._r8, ixcldice, & longname='Grid box averaged cloud ice amount', is_convtran1=.true.) call pbuf_add_field('QCWAT', 'global', dtype_r8, (/pcols,pver,dyn_time_lvls/), qcwat_idx) call pbuf_add_field('LCWAT', 'global', dtype_r8, (/pcols,pver,dyn_time_lvls/), lcwat_idx) call pbuf_add_field('TCWAT', 'global', dtype_r8, (/pcols,pver,dyn_time_lvls/), tcwat_idx) call pbuf_add_field('CLD', 'global', dtype_r8, (/pcols,pver,dyn_time_lvls/), cld_idx) call pbuf_add_field('AST', 'global', dtype_r8, (/pcols,pver,dyn_time_lvls/), ast_idx) call pbuf_add_field('CONCLD', 'global', dtype_r8, (/pcols,pver,dyn_time_lvls/), concld_idx) call pbuf_add_field('FICE', 'physpkg', dtype_r8, (/pcols,pver/), fice_idx) call pbuf_add_field('QME', 'physpkg', dtype_r8, (/pcols,pver/), qme_idx) call pbuf_add_field('PRAIN', 'physpkg', dtype_r8, (/pcols,pver/), prain_idx) call pbuf_add_field('NEVAPR', 'physpkg', dtype_r8, (/pcols,pver/), nevapr_idx) call pbuf_add_field('WSEDL', 'physpkg', dtype_r8, (/pcols,pver/), wsedl_idx) call pbuf_add_field('REI', 'physpkg', dtype_r8, (/pcols,pver/), rei_idx) call pbuf_add_field('REL', 'physpkg', dtype_r8, (/pcols,pver/), rel_idx) call pbuf_add_field('LS_FLXPRC', 'physpkg', dtype_r8, (/pcols,pverp/), ls_flxprc_idx) call pbuf_add_field('LS_FLXSNW', 'physpkg', dtype_r8, (/pcols,pverp/), ls_flxsnw_idx) end subroutine rk_stratiform_register !=============================================================================== function rk_stratiform_implements_cnst(name) !----------------------------------------------------------------------------- ! ! ! ! Return true if specified constituent is implemented by this package ! ! ! !----------------------------------------------------------------------------- ! character(len=*), intent(in) :: name ! constituent name logical :: rk_stratiform_implements_cnst ! return value !----------------------------------------------------------------------- rk_stratiform_implements_cnst = (do_cnst .and. any(name == cnst_names)) end function rk_stratiform_implements_cnst !=============================================================================== subroutine rk_stratiform_init_cnst(name, latvals, lonvals, mask, q) !----------------------------------------------------------------------- ! ! ! ! Initialize the cloud water mixing ratios (liquid and ice), if they are ! ! not read from the initial file ! ! ! !----------------------------------------------------------------------- ! character(len=*), intent(in) :: name ! constituent name real(r8), intent(in) :: latvals(:) ! lat in degrees (ncol) real(r8), intent(in) :: lonvals(:) ! lon in degrees (ncol) logical, intent(in) :: mask(:) ! Only initialize where .true. real(r8), intent(out) :: q(:,:) ! kg tracer/kg dry air (gcol, plev !----------------------------------------------------------------------- if (any(name == cnst_names)) q = 0.0_r8 end subroutine rk_stratiform_init_cnst !=============================================================================== subroutine rk_stratiform_init() !-------------------------------------------- ! ! ! ! Initialize the cloud water parameterization ! ! ! !-------------------------------------------- ! use physics_buffer, only: physics_buffer_desc, pbuf_get_index use constituents, only: cnst_get_ind, cnst_name, cnst_longname, sflxnam, apcnst, bpcnst use cam_history, only: addfld, add_default, horiz_only use convect_shallow, only: convect_shallow_use_shfrc use phys_control, only: cam_physpkg_is use physconst, only: tmelt, rhodair, rh2o use cldwat, only: inimc integer :: m, mm logical :: history_amwg ! output the variables used by the AMWG diag package logical :: history_aerosol ! Output the MAM aerosol tendencies logical :: history_budget ! Output tendencies and state variables for CAM4 ! temperature, water vapor, cloud ice and cloud ! liquid budgets. integer :: history_budget_histfile_num ! output history file number for budget fields !----------------------------------------------------------------------- call phys_getopts( history_aerosol_out = history_aerosol , & history_amwg_out = history_amwg , & history_budget_out = history_budget , & history_budget_histfile_num_out = history_budget_histfile_num) landm_idx = pbuf_get_index("LANDM") ! Find out whether shfrc from convect_shallow will be used in cldfrc if( convect_shallow_use_shfrc() ) then use_shfrc = .true. shfrc_idx = pbuf_get_index('shfrc') else use_shfrc = .false. endif ! Register history variables do m = 1, ncnst call cnst_get_ind( cnst_names(m), mm ) call addfld( cnst_name(mm), (/ 'lev' /), 'A', 'kg/kg', cnst_longname(mm) ) call addfld( sflxnam (mm), horiz_only, 'A', 'kg/m2/s', trim(cnst_name(mm))//' surface flux' ) if (history_amwg) then call add_default( cnst_name(mm), 1, ' ' ) call add_default( sflxnam (mm), 1, ' ' ) endif enddo call addfld (apcnst(ixcldliq), (/ 'lev' /), 'A', 'kg/kg', trim(cnst_name(ixcldliq))//' after physics' ) call addfld (apcnst(ixcldice), (/ 'lev' /), 'A', 'kg/kg', trim(cnst_name(ixcldice))//' after physics' ) call addfld (bpcnst(ixcldliq), (/ 'lev' /), 'A', 'kg/kg', trim(cnst_name(ixcldliq))//' before physics' ) call addfld (bpcnst(ixcldice), (/ 'lev' /), 'A', 'kg/kg', trim(cnst_name(ixcldice))//' before physics' ) if( history_budget) then call add_default (cnst_name(ixcldliq), history_budget_histfile_num, ' ') call add_default (cnst_name(ixcldice), history_budget_histfile_num, ' ') call add_default (apcnst (ixcldliq), history_budget_histfile_num, ' ') call add_default (apcnst (ixcldice), history_budget_histfile_num, ' ') call add_default (bpcnst (ixcldliq), history_budget_histfile_num, ' ') call add_default (bpcnst (ixcldice), history_budget_histfile_num, ' ') end if call addfld ('FWAUT', (/ 'lev' /), 'A', 'fraction', 'Relative importance of liquid autoconversion' ) call addfld ('FSAUT', (/ 'lev' /), 'A', 'fraction', 'Relative importance of ice autoconversion' ) call addfld ('FRACW', (/ 'lev' /), 'A', 'fraction', 'Relative importance of rain accreting liquid' ) call addfld ('FSACW', (/ 'lev' /), 'A', 'fraction', 'Relative importance of snow accreting liquid' ) call addfld ('FSACI', (/ 'lev' /), 'A', 'fraction', 'Relative importance of snow accreting ice' ) call addfld ('CME', (/ 'lev' /), 'A', 'kg/kg/s' , 'Rate of cond-evap within the cloud' ) call addfld ('CMEICE', (/ 'lev' /), 'A', 'kg/kg/s' , 'Rate of cond-evap of ice within the cloud' ) call addfld ('CMELIQ', (/ 'lev' /), 'A', 'kg/kg/s' , 'Rate of cond-evap of liq within the cloud' ) call addfld ('ICE2PR', (/ 'lev' /), 'A', 'kg/kg/s' , 'Rate of conversion of ice to precip' ) call addfld ('LIQ2PR', (/ 'lev' /), 'A', 'kg/kg/s' , 'Rate of conversion of liq to precip' ) call addfld ('ZMDLF', (/ 'lev' /), 'A', 'kg/kg/s' , 'Detrained liquid water from ZM convection' ) call addfld ('SHDLF', (/ 'lev' /), 'A', 'kg/kg/s' , 'Detrained liquid water from shallow convection' ) call addfld ('PRODPREC', (/ 'lev' /), 'A', 'kg/kg/s' , 'Rate of conversion of condensate to precip' ) call addfld ('EVAPPREC', (/ 'lev' /), 'A', 'kg/kg/s' , 'Rate of evaporation of falling precip' ) call addfld ('EVAPSNOW', (/ 'lev' /), 'A', 'kg/kg/s' , 'Rate of evaporation of falling snow' ) call addfld ('HPROGCLD', (/ 'lev' /), 'A', 'W/kg' , 'Heating from prognostic clouds' ) call addfld ('HCME', (/ 'lev' /), 'A', 'W/kg' , 'Heating from cond-evap within the cloud' ) call addfld ('HEVAP', (/ 'lev' /), 'A', 'W/kg' , 'Heating from evaporation of falling precip' ) call addfld ('HFREEZ', (/ 'lev' /), 'A', 'W/kg' , 'Heating rate due to freezing of precip' ) call addfld ('HMELT', (/ 'lev' /), 'A', 'W/kg' , 'Heating from snow melt' ) call addfld ('HREPART', (/ 'lev' /), 'A', 'W/kg' , 'Heating from cloud ice/liquid repartitioning' ) call addfld ('REPARTICE', (/ 'lev' /), 'A', 'kg/kg/s' , 'Cloud ice tendency from cloud ice/liquid repartitioning' ) call addfld ('REPARTLIQ', (/ 'lev' /), 'A', 'kg/kg/s' , 'Cloud liq tendency from cloud ice/liquid repartitioning' ) call addfld ('FICE', (/ 'lev' /), 'A', 'fraction', 'Fractional ice content within cloud' ) call addfld ('ICWMR', (/ 'lev' /), 'A', 'kg/kg' , 'Prognostic in-cloud water mixing ratio' ) call addfld ('ICIMR', (/ 'lev' /), 'A', 'kg/kg' , 'Prognostic in-cloud ice mixing ratio' ) call addfld ('PCSNOW', horiz_only , 'A', 'm/s' , 'Snow fall from prognostic clouds' ) call addfld ('DQSED', (/ 'lev' /), 'A', 'kg/kg/s' , 'Water vapor tendency from cloud sedimentation' ) call addfld ('DLSED', (/ 'lev' /), 'A', 'kg/kg/s' , 'Cloud liquid tendency from sedimentation' ) call addfld ('DISED', (/ 'lev' /), 'A', 'kg/kg/s' , 'Cloud ice tendency from sedimentation' ) call addfld ('HSED', (/ 'lev' /), 'A', 'W/kg' , 'Heating from cloud sediment evaporation' ) call addfld ('SNOWSED', horiz_only, 'A', 'm/s' , 'Snow from cloud ice sedimentation' ) call addfld ('RAINSED', horiz_only, 'A', 'm/s' , 'Rain from cloud liquid sedimentation' ) call addfld ('PRECSED', horiz_only, 'A', 'm/s' , 'Precipitation from cloud sedimentation' ) call addfld ('CNVCLD', horiz_only, 'A', 'fraction', 'Vertically integrated convective cloud amount' ) call addfld ('CLDST', (/ 'lev' /), 'A', 'fraction', 'Stratus cloud fraction' ) call addfld ('CONCLD', (/ 'lev' /), 'A', 'fraction', 'Convective cloud cover' ) call addfld ('AST', (/ 'lev' /), 'A','fraction' , 'Stratus cloud fraction' ) call addfld ('LIQCLDF', (/ 'lev' /), 'A', 'fraction', 'Stratus Liquid cloud fraction' ) call addfld ('ICECLDF', (/ 'lev' /), 'A', 'fraction', 'Stratus ICE cloud fraction' ) call addfld ('IWC', (/ 'lev' /), 'A', 'kg/m3' , 'Grid box average ice water content' ) call addfld ('LWC', (/ 'lev' /), 'A', 'kg/m3' , 'Grid box average liquid water content' ) call addfld ('ICWNC', (/ 'lev' /), 'A', 'm-3' , 'Prognostic in-cloud water number conc' ) call addfld ('ICINC', (/ 'lev' /), 'A', 'm-3' , 'Prognostic in-cloud ice number conc' ) call addfld ('REL', (/ 'lev' /), 'A', 'micron' , 'effective liquid drop radius' ) call addfld ('REI', (/ 'lev' /), 'A', 'micron' , 'effective ice particle radius' ) if ( history_budget ) then call add_default ('EVAPSNOW ', history_budget_histfile_num, ' ') call add_default ('EVAPPREC ', history_budget_histfile_num, ' ') call add_default ('CMELIQ ', history_budget_histfile_num, ' ') if( cam_physpkg_is('cam3') .or. cam_physpkg_is('cam4') ) then call add_default ('ZMDLF ', history_budget_histfile_num, ' ') call add_default ('CME ', history_budget_histfile_num, ' ') call add_default ('DQSED ', history_budget_histfile_num, ' ') call add_default ('DISED ', history_budget_histfile_num, ' ') call add_default ('DLSED ', history_budget_histfile_num, ' ') call add_default ('HSED ', history_budget_histfile_num, ' ') call add_default ('CMEICE ', history_budget_histfile_num, ' ') call add_default ('LIQ2PR ', history_budget_histfile_num, ' ') call add_default ('ICE2PR ', history_budget_histfile_num, ' ') call add_default ('HCME ', history_budget_histfile_num, ' ') call add_default ('HEVAP ', history_budget_histfile_num, ' ') call add_default ('HFREEZ ', history_budget_histfile_num, ' ') call add_default ('HMELT ', history_budget_histfile_num, ' ') call add_default ('HREPART ', history_budget_histfile_num, ' ') call add_default ('HPROGCLD ', history_budget_histfile_num, ' ') call add_default ('REPARTLIQ', history_budget_histfile_num, ' ') call add_default ('REPARTICE', history_budget_histfile_num, ' ') end if end if if (history_amwg) then call add_default ('ICWMR', 1, ' ') call add_default ('ICIMR', 1, ' ') call add_default ('CONCLD ', 1, ' ') call add_default ('FICE ', 1, ' ') endif ! History Variables for COSP/CFMIP call addfld ('LS_FLXPRC', (/ 'ilev' /), 'A', 'kg/m2/s', 'ls stratiform gbm interface rain+snow flux') call addfld ('LS_FLXSNW', (/ 'ilev' /), 'A', 'kg/m2/s', 'ls stratiform gbm interface snow flux') call addfld ('PRACWO', (/ 'lev' /), 'A', '1/s', 'Accretion of cloud water by rain') call addfld ('PSACWO', (/ 'lev' /), 'A', '1/s', 'Accretion of cloud water by snow') call addfld ('PSACIO', (/ 'lev' /), 'A', '1/s', 'Accretion of cloud ice by snow') call addfld ('CLDLIQSTR', (/ 'lev' /), 'A', 'kg/kg', 'Stratiform CLDLIQ') call addfld ('CLDICESTR', (/ 'lev' /), 'A', 'kg/kg', 'Stratiform CLDICE') call addfld ('CLDLIQCON', (/ 'lev' /), 'A', 'kg/kg', 'Convective CLDLIQ') call addfld ('CLDICECON', (/ 'lev' /), 'A', 'kg/kg', 'Convective CLDICE') cmfmc_sh_idx = pbuf_get_index('CMFMC_SH') prec_str_idx = pbuf_get_index('PREC_STR') snow_str_idx = pbuf_get_index('SNOW_STR') prec_pcw_idx = pbuf_get_index('PREC_PCW') snow_pcw_idx = pbuf_get_index('SNOW_PCW') prec_sed_idx = pbuf_get_index('PREC_SED') snow_sed_idx = pbuf_get_index('SNOW_SED') ! Initialize cldwat with constants. call inimc(tmelt, rhodair/1000.0_r8, gravit, rh2o) end subroutine rk_stratiform_init !=============================================================================== subroutine rk_stratiform_tend( & state, ptend_all, pbuf, dtime, icefrac, & landfrac, ocnfrac, snowh, dlf, & dlf2, rliq, cmfmc, ts, & sst, zdu) !-------------------------------------------------------- ! ! ! ! Interface to sedimentation, detrain, cloud fraction and ! ! cloud macro - microphysics subroutines ! ! ! !-------------------------------------------------------- ! use cloud_fraction, only: cldfrc, cldfrc_fice 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 use physics_types, only: physics_state_dealloc use cam_history, only: outfld use physics_buffer, only: physics_buffer_desc, pbuf_get_field, pbuf_old_tim_idx use pkg_cld_sediment, only: cld_sediment_vel, cld_sediment_tend use cldwat, only: pcond use pkg_cldoptics, only: cldefr use phys_control, only: cam_physpkg_is use tropopause, only: tropopause_find, TROP_ALG_TWMO, TROP_ALG_CLIMATE use phys_grid, only: get_rlat_all_p use physconst, only: pi ! Arguments type(physics_state), intent(in) :: state ! State variables type(physics_ptend), intent(out) :: ptend_all ! Package tendencies type(physics_buffer_desc), pointer :: pbuf(:) real(r8), intent(in) :: dtime ! Timestep real(r8), intent(in) :: icefrac (pcols) ! Sea ice fraction (fraction) real(r8), intent(in) :: landfrac(pcols) ! Land fraction (fraction) real(r8), intent(in) :: ocnfrac (pcols) ! Ocean fraction (fraction) real(r8), intent(in) :: snowh(pcols) ! Snow depth over land, water equivalent (m) real(r8), intent(in) :: dlf(pcols,pver) ! Detrained water from convection schemes real(r8), intent(in) :: dlf2(pcols,pver) ! Detrained water from shallow convection scheme real(r8), intent(in) :: rliq(pcols) ! Vertical integral of liquid not yet in q(ixcldliq) real(r8), intent(in) :: cmfmc(pcols,pverp) ! Deep + Shallow Convective mass flux [ kg /s/m^2 ] real(r8), intent(in) :: ts(pcols) ! Surface temperature real(r8), intent(in) :: sst(pcols) ! Sea surface temperature real(r8), intent(in) :: zdu(pcols,pver) ! Detrainment rate from deep convection ! Local variables type(physics_state) :: state1 ! Local copy of the state variable type(physics_ptend) :: ptend_loc ! Package tendencies integer :: i, k, m integer :: lchnk ! Chunk identifier integer :: ncol ! Number of atmospheric columns integer :: itim_old ! Physics buffer fields real(r8), pointer :: landm(:) ! Land fraction ramped over water real(r8), pointer :: prec_str(:) ! [Total] Sfc flux of precip from stratiform [ m/s ] real(r8), pointer :: snow_str(:) ! [Total] Sfc flux of snow from stratiform [ m/s ] real(r8), pointer :: prec_sed(:) ! Surface flux of total cloud water from sedimentation real(r8), pointer :: snow_sed(:) ! Surface flux of cloud ice from sedimentation real(r8), pointer :: prec_pcw(:) ! Sfc flux of precip from microphysics [ m/s ] real(r8), pointer :: snow_pcw(:) ! Sfc flux of snow from microphysics [ m/s ] real(r8), pointer, dimension(:,:) :: qcwat ! Cloud water old q real(r8), pointer, dimension(:,:) :: tcwat ! Cloud water old temperature real(r8), pointer, dimension(:,:) :: lcwat ! Cloud liquid water old q real(r8), pointer, dimension(:,:) :: cld ! Total cloud fraction real(r8), pointer, dimension(:,:) :: fice ! Cloud ice/water partitioning ratio. real(r8), pointer, dimension(:,:) :: ast ! Relative humidity cloud fraction real(r8), pointer, dimension(:,:) :: concld ! Convective cloud fraction real(r8), pointer, dimension(:,:) :: qme ! rate of cond-evap of condensate (1/s) real(r8), pointer, dimension(:,:) :: prain ! Total precipitation (rain + snow) real(r8), pointer, dimension(:,:) :: nevapr ! Evaporation of total precipitation (rain + snow) real(r8), pointer, dimension(:,:) :: rel ! Liquid effective drop radius (microns) real(r8), pointer, dimension(:,:) :: rei ! Ice effective drop size (microns) real(r8), pointer, dimension(:,:) :: wsedl ! Sedimentation velocity of liquid stratus cloud droplet [ m/s ] real(r8), pointer, dimension(:,:) :: shfrc ! Cloud fraction from shallow convection scheme real(r8), pointer, dimension(:,:) :: cmfmc_sh ! Shallow convective mass flux (pcols,pverp) [ kg/s/m^2 ] real(r8), target :: shfrc_local(pcols,pver) ! physics buffer fields for COSP simulator (RK only) real(r8), pointer, dimension(:,:) :: rkflxprc ! RK grid-box mean flux_large_scale_cloud_rain+snow at interfaces (kg/m2/s) real(r8), pointer, dimension(:,:) :: rkflxsnw ! RK grid-box mean flux_large_scale_cloud_snow at interfaces (kg/m2/s) ! Local variables for stratiform_sediment real(r8) :: rain(pcols) ! Surface flux of cloud liquid real(r8) :: pvliq(pcols,pverp) ! Vertical velocity of cloud liquid drops (Pa/s) real(r8) :: pvice(pcols,pverp) ! Vertical velocity of cloud ice particles (Pa/s) ! Local variables for cldfrc real(r8) :: cldst(pcols,pver) ! Stratus cloud fraction real(r8) :: rhcloud(pcols,pver) ! Relative humidity cloud (last timestep) real(r8) :: rhcloud2(pcols,pver) ! Relative humidity cloud (perturbation) real(r8) :: clc(pcols) ! Column convective cloud amount real(r8) :: relhum(pcols,pver) ! RH, output to determine drh/da real(r8) :: rhu00(pcols,pver) real(r8) :: rhu002(pcols,pver) ! Same as rhu00 but for perturbed rh real(r8) :: rhdfda(pcols,pver) real(r8) :: cld2(pcols,pver) ! Same as cld but for perturbed rh real(r8) :: concld2(pcols,pver) ! Same as concld but for perturbed rh real(r8) :: cldst2(pcols,pver) ! Same as cldst but for perturbed rh real(r8) :: relhum2(pcols,pver) ! RH after perturbation real(r8) :: icecldf(pcols,pver) ! Ice cloud fraction real(r8) :: liqcldf(pcols,pver) ! Liquid cloud fraction (combined into cloud) real(r8) :: icecldf_out(pcols,pver) ! Ice cloud fraction real(r8) :: liqcldf_out(pcols,pver) ! Liquid cloud fraction (combined into cloud) real(r8) :: icecldf2(pcols,pver) ! Ice cloud fraction real(r8) :: liqcldf2(pcols,pver) ! Liquid cloud fraction (combined into cloud) ! Local variables for microphysics real(r8) :: rdtime ! 1./dtime real(r8) :: qtend(pcols,pver) ! Moisture tendencies real(r8) :: ttend(pcols,pver) ! Temperature tendencies real(r8) :: ltend(pcols,pver) ! Cloud liquid water tendencies real(r8) :: evapheat(pcols,pver) ! Heating rate due to evaporation of precip real(r8) :: evapsnow(pcols,pver) ! Local evaporation of snow real(r8) :: prfzheat(pcols,pver) ! Heating rate due to freezing of precip (W/kg) real(r8) :: meltheat(pcols,pver) ! Heating rate due to phase change of precip real(r8) :: cmeheat (pcols,pver) ! Heating rate due to phase change of precip real(r8) :: prodsnow(pcols,pver) ! Local production of snow real(r8) :: totcw(pcols,pver) ! Total cloud water mixing ratio real(r8) :: fsnow(pcols,pver) ! Fractional snow production real(r8) :: repartht(pcols,pver) ! Heating rate due to phase repartition of input precip real(r8) :: icimr(pcols,pver) ! In cloud ice mixing ratio real(r8) :: icwmr(pcols,pver) ! In cloud water mixing ratio real(r8) :: fwaut(pcols,pver) real(r8) :: fsaut(pcols,pver) real(r8) :: fracw(pcols,pver) real(r8) :: fsacw(pcols,pver) real(r8) :: fsaci(pcols,pver) real(r8) :: cmeice(pcols,pver) ! Rate of cond-evap of ice within the cloud real(r8) :: cmeliq(pcols,pver) ! Rate of cond-evap of liq within the cloud real(r8) :: ice2pr(pcols,pver) ! Rate of conversion of ice to precip real(r8) :: liq2pr(pcols,pver) ! Rate of conversion of liquid to precip real(r8) :: liq2snow(pcols,pver) ! Rate of conversion of liquid to snow ! Local variables for CFMIP calculations real(r8) :: mr_lsliq(pcols,pver) ! mixing_ratio_large_scale_cloud_liquid (kg/kg) real(r8) :: mr_lsice(pcols,pver) ! mixing_ratio_large_scale_cloud_ice (kg/kg) real(r8) :: mr_ccliq(pcols,pver) ! mixing_ratio_convective_cloud_liquid (kg/kg) real(r8) :: mr_ccice(pcols,pver) ! mixing_ratio_convective_cloud_ice (kg/kg) real(r8) :: pracwo(pcols,pver) ! RK accretion of cloud water by rain (1/s) real(r8) :: psacwo(pcols,pver) ! RK accretion of cloud water by snow (1/s) real(r8) :: psacio(pcols,pver) ! RK accretion of cloud ice by snow (1/s) real(r8) :: iwc(pcols,pver) ! Grid box average ice water content real(r8) :: lwc(pcols,pver) ! Grid box average liquid water content logical :: lq(pcnst) integer :: troplev(pcols) real(r8) :: rlat(pcols) real(r8) :: dlat(pcols) real(r8), parameter :: rad2deg = 180._r8/pi ! ====================================================================== lchnk = state%lchnk ncol = state%ncol call physics_state_copy(state,state1) ! Copy state to local state1. ! Associate pointers with physics buffer fields call pbuf_get_field(pbuf, landm_idx, landm) itim_old = pbuf_old_tim_idx() call pbuf_get_field(pbuf, qcwat_idx, qcwat, start=(/1,1,itim_old/), kount=(/pcols,pver,1/) ) call pbuf_get_field(pbuf, tcwat_idx, tcwat, start=(/1,1,itim_old/), kount=(/pcols,pver,1/) ) call pbuf_get_field(pbuf, lcwat_idx, lcwat, start=(/1,1,itim_old/), kount=(/pcols,pver,1/) ) call pbuf_get_field(pbuf, cld_idx, cld, start=(/1,1,itim_old/), kount=(/pcols,pver,1/) ) call pbuf_get_field(pbuf, concld_idx, concld, start=(/1,1,itim_old/), kount=(/pcols,pver,1/) ) call pbuf_get_field(pbuf, ast_idx, ast, start=(/1,1,itim_old/), kount=(/pcols,pver,1/) ) call pbuf_get_field(pbuf, fice_idx, fice) call pbuf_get_field(pbuf, cmfmc_sh_idx, cmfmc_sh) call pbuf_get_field(pbuf, prec_str_idx, prec_str) call pbuf_get_field(pbuf, snow_str_idx, snow_str) call pbuf_get_field(pbuf, prec_sed_idx, prec_sed) call pbuf_get_field(pbuf, snow_sed_idx, snow_sed) call pbuf_get_field(pbuf, prec_pcw_idx, prec_pcw) call pbuf_get_field(pbuf, snow_pcw_idx, snow_pcw) call pbuf_get_field(pbuf, qme_idx, qme ) call pbuf_get_field(pbuf, prain_idx, prain) call pbuf_get_field(pbuf, nevapr_idx, nevapr) call pbuf_get_field(pbuf, rel_idx, rel) call pbuf_get_field(pbuf, rei_idx, rei) call pbuf_get_field(pbuf, wsedl_idx, wsedl) if ( do_psrhmin ) then call tropopause_find(state, troplev, primary=TROP_ALG_TWMO, backup=TROP_ALG_CLIMATE) call get_rlat_all_p(lchnk,ncol,rlat) dlat = rlat*rad2deg endif ! ------------- ! ! Sedimentation ! ! ------------- ! ! Allow the cloud liquid drops and ice particles to sediment. ! This is done before adding convectively detrained cloud water, ! because the phase of the detrained water is unknown. call t_startf('stratiform_sediment') call cld_sediment_vel( ncol, & icefrac, landfrac, ocnfrac, state1%pmid, state1%pdel, state1%t, & cld, state1%q(:,:,ixcldliq), state1%q(:,:,ixcldice), & pvliq, pvice, landm, snowh ) wsedl(:ncol,:pver) = pvliq(:ncol,:pver)/gravit/(state1%pmid(:ncol,:pver)/(287.15_r8*state1%t(:ncol,:pver))) lq(:) = .FALSE. lq(1) = .TRUE. lq(ixcldice) = .TRUE. lq(ixcldliq) = .TRUE. call physics_ptend_init(ptend_loc, state%psetcols, 'pcwsediment', ls=.true., lq=lq)! Initialize local ptend type call cld_sediment_tend( ncol, dtime , & state1%pint, state1%pmid, state1%pdel, state1%t, & cld, state1%q(:,:,ixcldliq), state1%q(:,:,ixcldice), pvliq, pvice, & ptend_loc%q(:,:,ixcldliq), ptend_loc%q(:,:,ixcldice), ptend_loc%q(:,:,1), & ptend_loc%s, rain, snow_sed ) ! Convert rain and snow fluxes at the surface from [kg/m2/s] to [m/s] ! Compute total precipitation flux at the surface in [m/s] snow_sed(:ncol) = snow_sed(:ncol)/1000._r8 rain(:ncol) = rain(:ncol)/1000._r8 prec_sed(:ncol) = rain(:ncol) + snow_sed(:ncol) ! Record history variables call outfld( 'DQSED' ,ptend_loc%q(:,:,1) , pcols,lchnk ) call outfld( 'DISED' ,ptend_loc%q(:,:,ixcldice), pcols,lchnk ) call outfld( 'DLSED' ,ptend_loc%q(:,:,ixcldliq), pcols,lchnk ) call outfld( 'HSED' ,ptend_loc%s , pcols,lchnk ) call outfld( 'PRECSED' ,prec_sed , pcols,lchnk ) call outfld( 'SNOWSED' ,snow_sed , pcols,lchnk ) call outfld( 'RAINSED' ,rain , pcols,lchnk ) ! Add tendency from this process to tend from other processes here call physics_ptend_init(ptend_all, state%psetcols, 'stratiform') call physics_ptend_sum( ptend_loc, ptend_all, ncol ) ! Update physics state type state1 with ptend_loc call physics_update( state1, ptend_loc, dtime ) call t_stopf('stratiform_sediment') ! Accumulate prec and snow flux at the surface [ m/s ] prec_str(:ncol) = prec_sed(:ncol) snow_str(:ncol) = snow_sed(:ncol) ! ----------------------------------------------------------------------------- ! ! Detrainment of convective condensate into the environment or stratiform cloud ! ! ----------------------------------------------------------------------------- ! ! Put all of the detraining cloud water from convection into the large scale cloud. ! It all goes in liquid for the moment. ! Strictly speaking, this approach is detraining all the cconvective water into ! the environment, not the large-scale cloud. lq(:) = .FALSE. lq(ixcldliq) = .TRUE. call physics_ptend_init( ptend_loc, state1%psetcols, 'pcwdetrain', lq=lq) do k = 1, pver do i = 1, state1%ncol ptend_loc%q(i,k,ixcldliq) = dlf(i,k) end do end do call outfld( 'ZMDLF', dlf, pcols, lchnk ) call outfld( 'SHDLF', dlf2, pcols, lchnk ) ! Add hie detrainment tendency to tend from the other prior processes call physics_ptend_sum( ptend_loc, ptend_all, ncol ) call physics_update( state1, ptend_loc, dtime ) ! Accumulate prec and snow, reserved liquid has now been used. prec_str(:ncol) = prec_str(:ncol) - rliq(:ncol) ! ( snow contribution is zero ) ! -------------------------------------- ! ! Computation of Various Cloud Fractions ! ! -------------------------------------- ! ! ----------------------------------------------------------------------------- ! ! Treatment of cloud fraction in CAM4 and CAM5 differs ! ! (1) CAM4 ! ! . Cumulus AMT = Deep Cumulus AMT ( empirical fcn of mass flux ) + ! ! Shallow Cumulus AMT ( empirical fcn of mass flux ) ! ! . Stratus AMT = max( RH stratus AMT, Stability Stratus AMT ) ! ! . Cumulus and Stratus are 'minimally' overlapped without hierarchy. ! ! . Cumulus LWC,IWC is assumed to be the same as Stratus LWC,IWC ! ! (2) CAM5 ! ! . Cumulus AMT = Deep Cumulus AMT ( empirical fcn of mass flux ) + ! ! Shallow Cumulus AMT ( internally fcn of mass flux and w ) ! ! . Stratus AMT = fcn of environmental-mean RH ( no Stability Stratus ) ! ! . Cumulus and Stratus are non-overlapped with higher priority on Cumulus ! ! . Cumulus ( both Deep and Shallow ) has its own LWC and IWC. ! ! ----------------------------------------------------------------------------- ! if( use_shfrc ) then call pbuf_get_field(pbuf, shfrc_idx, shfrc ) else shfrc=>shfrc_local shfrc(:,:) = 0._r8 endif ! Stratus ('ast' = max(alst,aist)) and total cloud fraction ('cld = ast + concld') ! will be computed using this updated 'concld' in the stratiform macrophysics ! scheme (mmacro_pcond) later below. call t_startf("cldfrc") call cldfrc( lchnk, ncol, pbuf, & state1%pmid, state1%t, state1%q(:,:,1), state1%omega, state1%phis, & shfrc, use_shfrc, & cld, rhcloud, clc, state1%pdel, & cmfmc, cmfmc_sh, landfrac,snowh, concld, cldst, & ts, sst, state1%pint(:,pverp), zdu, ocnfrac, rhu00, & state1%q(:,:,ixcldice), icecldf, liqcldf, & relhum, 0 ) ! Re-calculate cloud with perturbed rh add call cldfrc to estimate rhdfda. call cldfrc( lchnk, ncol, pbuf, & state1%pmid, state1%t, state1%q(:,:,1), state1%omega, state1%phis, & shfrc, use_shfrc, & cld2, rhcloud2, clc, state1%pdel, & cmfmc, cmfmc_sh, landfrac, snowh, concld2, cldst2, & ts, sst, state1%pint(:,pverp), zdu, ocnfrac, rhu002, & state1%q(:,:,ixcldice), icecldf2, liqcldf2, & relhum2, 1 ) call t_stopf("cldfrc") rhu00(:ncol,1) = 2.0_r8 do k = 1, pver do i = 1, ncol if( relhum(i,k) < rhu00(i,k) ) then rhdfda(i,k) = 0.0_r8 elseif( relhum(i,k) >= 1.0_r8 ) then rhdfda(i,k) = 0.0_r8 else ! Under certain circumstances, rh+ cause cld not to changed ! when at an upper limit, or w/ strong subsidence if( ( cld2(i,k) - cld(i,k) ) < 1.e-4_r8 ) then rhdfda(i,k) = 0.01_r8*relhum(i,k)*1.e+4_r8 else rhdfda(i,k) = 0.01_r8*relhum(i,k)/(cld2(i,k)-cld(i,k)) endif endif enddo enddo ! ---------------------------------------------- ! ! Stratiform Cloud Macrophysics and Microphysics ! ! ---------------------------------------------- ! call t_startf('stratiform_microphys') rdtime = 1._r8/dtime ! Define fractional amount of stratus condensate and precipitation in ice phase. ! This uses a ramp ( -30 ~ -10 for fice, -5 ~ 0 for fsnow ). ! The ramp within convective cloud may be different call cldfrc_fice(ncol, state1%t, fice, fsnow) ! Perform repartitioning of stratiform condensate. ! Corresponding heating tendency will be added later. lq(:) = .FALSE. lq(ixcldice) = .true. lq(ixcldliq) = .true. call physics_ptend_init( ptend_loc, state1%psetcols, 'cldwat-repartition', lq=lq ) totcw(:ncol,:pver) = state1%q(:ncol,:pver,ixcldice) + state1%q(:ncol,:pver,ixcldliq) repartht(:ncol,:pver) = state1%q(:ncol,:pver,ixcldice) ptend_loc%q(:ncol,:pver,ixcldice) = rdtime * ( totcw(:ncol,:pver)*fice(:ncol,:pver) - state1%q(:ncol,:pver,ixcldice) ) ptend_loc%q(:ncol,:pver,ixcldliq) = rdtime * ( totcw(:ncol,:pver)*(1.0_r8-fice(:ncol,:pver)) - state1%q(:ncol,:pver,ixcldliq) ) call outfld( 'REPARTICE', ptend_loc%q(:,:,ixcldice), pcols, lchnk ) call outfld( 'REPARTLIQ', ptend_loc%q(:,:,ixcldliq), pcols, lchnk ) call physics_ptend_sum( ptend_loc, ptend_all, ncol ) call physics_update( state1, ptend_loc, dtime ) ! Determine repartition heating from change in cloud ice. repartht(:ncol,:pver) = (latice/dtime) * ( state1%q(:ncol,:pver,ixcldice) - repartht(:ncol,:pver) ) ! Non-micro and non-macrophysical external advective forcings to compute net condensation rate. ! Note that advective forcing of condensate is aggregated into liquid phase. qtend(:ncol,:pver) = ( state1%q(:ncol,:pver,1) - qcwat(:ncol,:pver) ) * rdtime ttend(:ncol,:pver) = ( state1%t(:ncol,:pver) - tcwat(:ncol,:pver) ) * rdtime ltend(:ncol,:pver) = ( totcw (:ncol,:pver) - lcwat(:ncol,:pver) ) * rdtime ! Compute Stratiform Macro-Microphysical Tendencies ! Add rain and snow fluxes as output variables from pcond, and into physics buffer call pbuf_get_field(pbuf, ls_flxprc_idx, rkflxprc) call pbuf_get_field(pbuf, ls_flxsnw_idx, rkflxsnw) call t_startf('pcond') call pcond( lchnk, ncol, troplev, dlat, & state1%t, ttend, state1%q(1,1,1), qtend, state1%omega, & totcw, state1%pmid , state1%pdel, cld, fice, fsnow, & qme, prain, prodsnow, nevapr, evapsnow, evapheat, prfzheat, & meltheat, prec_pcw, snow_pcw, dtime, fwaut, & fsaut, fracw, fsacw, fsaci, ltend, & rhdfda, rhu00, landm, icefrac, state1%zi, ice2pr, liq2pr, & liq2snow, snowh, rkflxprc, rkflxsnw, pracwo, psacwo, psacio ) call t_stopf('pcond') lq(:) = .FALSE. lq(1) = .true. lq(ixcldice) = .true. lq(ixcldliq) = .true. call physics_ptend_init( ptend_loc, state1%psetcols, 'cldwat', ls=.true., lq=lq) do k = 1, pver do i = 1, ncol ptend_loc%s(i,k) = qme(i,k)*( latvap + latice*fice(i,k) ) + & evapheat(i,k) + prfzheat(i,k) + meltheat(i,k) + repartht(i,k) ptend_loc%q(i,k,1) = - qme(i,k) + nevapr(i,k) ptend_loc%q(i,k,ixcldice) = qme(i,k)*fice(i,k) - ice2pr(i,k) ptend_loc%q(i,k,ixcldliq) = qme(i,k)*(1._r8-fice(i,k)) - liq2pr(i,k) end do end do do k = 1, pver do i = 1, ncol ast(i,k) = cld(i,k) icimr(i,k) = (state1%q(i,k,ixcldice) + dtime*ptend_loc%q(i,k,ixcldice)) / max(0.01_r8,ast(i,k)) icwmr(i,k) = (state1%q(i,k,ixcldliq) + dtime*ptend_loc%q(i,k,ixcldliq)) / max(0.01_r8,ast(i,k)) end do end do ! Convert precipitation from [ kg/m2 ] to [ m/s ] snow_pcw(:ncol) = snow_pcw(:ncol)/1000._r8 prec_pcw(:ncol) = prec_pcw(:ncol)/1000._r8 do k = 1, pver do i = 1, ncol cmeheat(i,k) = qme(i,k) * ( latvap + latice*fice(i,k) ) cmeice (i,k) = qme(i,k) * fice(i,k) cmeliq (i,k) = qme(i,k) * ( 1._r8 - fice(i,k) ) end do end do ! Record history variables call outfld( 'FWAUT' , fwaut, pcols, lchnk ) call outfld( 'FSAUT' , fsaut, pcols, lchnk ) call outfld( 'FRACW' , fracw, pcols, lchnk ) call outfld( 'FSACW' , fsacw, pcols, lchnk ) call outfld( 'FSACI' , fsaci, pcols, lchnk ) call outfld( 'PCSNOW' , snow_pcw, pcols, lchnk ) call outfld( 'FICE' , fice, pcols, lchnk ) call outfld( 'CMEICE' , cmeice, pcols, lchnk ) call outfld( 'CMELIQ' , cmeliq, pcols, lchnk ) call outfld( 'ICE2PR' , ice2pr, pcols, lchnk ) call outfld( 'LIQ2PR' , liq2pr, pcols, lchnk ) call outfld( 'HPROGCLD', ptend_loc%s, pcols, lchnk ) call outfld( 'HEVAP ', evapheat, pcols, lchnk ) call outfld( 'HMELT' , meltheat, pcols, lchnk ) call outfld( 'HCME' , cmeheat , pcols, lchnk ) call outfld( 'HFREEZ' , prfzheat, pcols, lchnk ) call outfld( 'HREPART' , repartht, pcols, lchnk ) call outfld('LS_FLXPRC', rkflxprc, pcols, lchnk ) call outfld('LS_FLXSNW', rkflxsnw, pcols, lchnk ) call outfld('PRACWO' , pracwo, pcols, lchnk ) call outfld('PSACWO' , psacwo, pcols, lchnk ) call outfld('PSACIO' , psacio, pcols, lchnk ) ! initialize local variables mr_ccliq(1:ncol,1:pver) = 0._r8 mr_ccice(1:ncol,1:pver) = 0._r8 mr_lsliq(1:ncol,1:pver) = 0._r8 mr_lsice(1:ncol,1:pver) = 0._r8 do k=1,pver do i=1,ncol if (cld(i,k) .gt. 0._r8) then mr_ccliq(i,k) = (state%q(i,k,ixcldliq)/cld(i,k))*concld(i,k) mr_ccice(i,k) = (state%q(i,k,ixcldice)/cld(i,k))*concld(i,k) mr_lsliq(i,k) = (state%q(i,k,ixcldliq)/cld(i,k))*(cld(i,k)-concld(i,k)) mr_lsice(i,k) = (state%q(i,k,ixcldice)/cld(i,k))*(cld(i,k)-concld(i,k)) else mr_ccliq(i,k) = 0._r8 mr_ccice(i,k) = 0._r8 mr_lsliq(i,k) = 0._r8 mr_lsice(i,k) = 0._r8 end if end do end do call outfld( 'CLDLIQSTR ', mr_lsliq, pcols, lchnk ) call outfld( 'CLDICESTR ', mr_lsice, pcols, lchnk ) call outfld( 'CLDLIQCON ', mr_ccliq, pcols, lchnk ) call outfld( 'CLDICECON ', mr_ccice, pcols, lchnk ) ! ------------------------------- ! ! Update microphysical tendencies ! ! ------------------------------- ! call physics_ptend_sum( ptend_loc, ptend_all, ncol ) call physics_update( state1, ptend_loc, dtime ) if (.not. cam_physpkg_is('cam3')) then call t_startf("cldfrc") call cldfrc( lchnk, ncol, pbuf, & state1%pmid, state1%t, state1%q(:,:,1), state1%omega, state1%phis, & shfrc, use_shfrc, & cld, rhcloud, clc, state1%pdel, & cmfmc, cmfmc_sh, landfrac, snowh, concld, cldst, & ts, sst, state1%pint(:,pverp), zdu, ocnfrac, rhu00, & state1%q(:,:,ixcldice), icecldf, liqcldf, & relhum, 0 ) call t_stopf("cldfrc") endif call outfld( 'CONCLD ', concld, pcols, lchnk ) call outfld( 'CLDST ', cldst, pcols, lchnk ) call outfld( 'CNVCLD ', clc, pcols, lchnk ) call outfld( 'AST', ast, pcols, lchnk ) do k = 1, pver do i = 1, ncol iwc(i,k) = state1%q(i,k,ixcldice)*state1%pmid(i,k)/(287.15_r8*state1%t(i,k)) lwc(i,k) = state1%q(i,k,ixcldliq)*state1%pmid(i,k)/(287.15_r8*state1%t(i,k)) icimr(i,k) = state1%q(i,k,ixcldice) / max(0.01_r8,rhcloud(i,k)) icwmr(i,k) = state1%q(i,k,ixcldliq) / max(0.01_r8,rhcloud(i,k)) end do end do call outfld( 'IWC' , iwc, pcols, lchnk ) call outfld( 'LWC' , lwc, pcols, lchnk ) call outfld( 'ICIMR' , icimr, pcols, lchnk ) call outfld( 'ICWMR' , icwmr, pcols, lchnk ) call outfld( 'CME' , qme, pcols, lchnk ) call outfld( 'PRODPREC' , prain, pcols, lchnk ) call outfld( 'EVAPPREC' , nevapr, pcols, lchnk ) call outfld( 'EVAPSNOW' , evapsnow, pcols, lchnk ) call t_stopf('stratiform_microphys') prec_str(:ncol) = prec_str(:ncol) + prec_pcw(:ncol) snow_str(:ncol) = snow_str(:ncol) + snow_pcw(:ncol) ! Save variables for use in the macrophysics at the next time step do k = 1, pver qcwat(:ncol,k) = state1%q(:ncol,k,1) tcwat(:ncol,k) = state1%t(:ncol,k) lcwat(:ncol,k) = state1%q(:ncol,k,ixcldice) + state1%q(:ncol,k,ixcldliq) end do ! Cloud water and ice particle sizes, saved in physics buffer for radiation call cldefr( lchnk, ncol, landfrac, state1%t, rel, rei, state1%ps, state1%pmid, landm, icefrac, snowh ) call outfld('REL', rel, pcols, lchnk) call outfld('REI', rei, pcols, lchnk) call physics_state_dealloc(state1) end subroutine rk_stratiform_tend !============================================================================ ! ! ! !============================================================================ ! subroutine debug_microphys_1(state1,ptend,i,k, & dtime,qme,fice,snow_pcw,prec_pcw, & prain,nevapr,prodsnow, evapsnow, & ice2pr,liq2pr,liq2snow) use physics_types, only: physics_state, physics_ptend use physconst, only: tmelt implicit none integer, intent(in) :: i,k type(physics_state), intent(in) :: state1 ! local copy of the state variable type(physics_ptend), intent(in) :: ptend ! local copy of the ptend variable real(r8), intent(in) :: dtime ! timestep real(r8), intent(in) :: qme(pcols,pver) ! local condensation - evaporation of cloud water real(r8), intent(in) :: prain(pcols,pver) ! local production of precipitation real(r8), intent(in) :: nevapr(pcols,pver) ! local evaporation of precipitation real(r8), intent(in) :: prodsnow(pcols,pver) ! local production of snow real(r8), intent(in) :: evapsnow(pcols,pver) ! local evaporation of snow real(r8), intent(in) :: ice2pr(pcols,pver) ! rate of conversion of ice to precip real(r8), intent(in) :: liq2pr(pcols,pver) ! rate of conversion of liquid to precip real(r8), intent(in) :: liq2snow(pcols,pver) ! rate of conversion of liquid to snow real(r8), intent(in) :: fice (pcols,pver) ! Fractional ice content within cloud real(r8), intent(in) :: snow_pcw(pcols) real(r8), intent(in) :: prec_pcw(pcols) real(r8) hs1, qv1, ql1, qi1, qs1, qr1, fice2, pr1, w1, w2, w3, fliq, res real(r8) w4, wl, wv, wi, wlf, wvf, wif, qif, qlf, qvf pr1 = 0 hs1 = 0 qv1 = 0 ql1 = 0 qi1 = 0 qs1 = 0 qr1 = 0 w1 = 0 wl = 0 wv = 0 wi = 0 wlf = 0 wvf = 0 wif = 0 write(iulog,*) write(iulog,*) ' input state, t, q, l, i ', k, state1%t(i,k), state1%q(i,k,1), state1%q(i,k,ixcldliq), state1%q(i,k,ixcldice) write(iulog,*) ' rain, snow, total from components before accumulation ', qr1, qs1, qr1+qs1 write(iulog,*) ' total precip before accumulation ', k, pr1 wv = wv + state1%q(i,k,1 )*state1%pdel(i,k)/gravit wl = wl + state1%q(i,k,ixcldliq)*state1%pdel(i,k)/gravit wi = wi + state1%q(i,k,ixcldice)*state1%pdel(i,k)/gravit qvf = state1%q(i,k,1) + ptend%q(i,k,1)*dtime qlf = state1%q(i,k,ixcldliq) + ptend%q(i,k,ixcldliq)*dtime qif = state1%q(i,k,ixcldice) + ptend%q(i,k,ixcldice)*dtime if (qvf.lt.0._r8) then write(iulog,*) ' qvf is negative *******', qvf endif if (qlf.lt.0._r8) then write(iulog,*) ' qlf is negative *******', qlf endif if (qif.lt.0._r8) then write(iulog,*) ' qif is negative *******', qif endif write(iulog,*) ' qvf, qlf, qif ', qvf, qlf, qif wvf = wvf + qvf*state1%pdel(i,k)/gravit wlf = wlf + qlf*state1%pdel(i,k)/gravit wif = wif + qif*state1%pdel(i,k)/gravit hs1 = hs1 + ptend%s(i,k)*state1%pdel(i,k)/gravit pr1 = pr1 + state1%pdel(i,k)/gravit*(prain(i,k)-nevapr(i,k)) qv1 = qv1 - (qme(i,k)-nevapr(i,k))*state1%pdel(i,k)/gravit ! vdot w1 = w1 + (qme(i,k)-prain(i,k))*state1%pdel(i,k)/gravit ! cdot qi1 = qi1 + ((qme(i,k))*fice(i,k) -ice2pr(i,k) )*state1%pdel(i,k)/gravit ! idot ql1 = ql1 + ((qme(i,k))*(1._r8-fice(i,k))-liq2pr(i,k) )*state1%pdel(i,k)/gravit ! ldot qr1 = qr1 & + ( liq2pr(i,k)-liq2snow(i,k) & ! production of rain -(nevapr(i,k)-evapsnow(i,k)) & ! rain evaporation )*state1%pdel(i,k)/gravit qs1 = qs1 & + ( ice2pr(i,k) + liq2snow(i,k) & ! production of snow.Note last term has phase change -evapsnow(i,k) & ! snow evaporation )*state1%pdel(i,k)/gravit if (state1%t(i,k).gt.tmelt) then qr1 = qr1 + qs1 qs1 = 0._r8 endif write(iulog,*) ' rain, snow, total after accumulation ', qr1, qs1, qr1+qs1 write(iulog,*) ' total precip after accumulation ', k, pr1 write(iulog,*) write(iulog,*) ' layer prain, nevapr, pdel ', prain(i,k), nevapr(i,k), state1%pdel(i,k) write(iulog,*) ' layer prodsnow, ice2pr+liq2snow ', prodsnow(i,k), ice2pr(i,k)+liq2snow(i,k) write(iulog,*) ' layer prain-prodsnow, liq2pr-liq2snow ', prain(i,k)-prodsnow(i,k), liq2pr(i,k)-liq2snow(i,k) write(iulog,*) ' layer evapsnow, evaprain ', k, evapsnow(i,k), nevapr(i,k)-evapsnow(i,k) write(iulog,*) ' layer ice2pr, liq2pr, liq2snow ', ice2pr(i,k), liq2pr(i,k), liq2snow(i,k) write(iulog,*) ' layer ice2pr+liq2pr, prain ', ice2pr(i,k)+liq2pr(i,k), prain(i,k) write(iulog,*) write(iulog,*) ' qv1 vapor removed from col after accum (vdot) ', k, qv1 write(iulog,*) ' - (precip produced - vapor removed) after accum ', k, -pr1-qv1 write(iulog,*) ' condensate produce after accum ', k, w1 write(iulog,*) ' liq+ice tends accum ', k, ql1+qi1 write(iulog,*) ' change in total water after accum ', k, qv1+ql1+qi1 write(iulog,*) ' imbalance in colum after accum ', k, qs1+qr1+qv1+ql1+qi1 write(iulog,*) ' fice at this lev ', fice(i,k) write(iulog,*) res = abs((qs1+qr1+qv1+ql1+qi1)/max(abs(qv1),abs(ql1),abs(qi1),abs(qs1),abs(qr1),1.e-36_r8)) write(iulog,*) ' relative residual in column method 1 ', k, res write(iulog,*) ' relative residual in column method 2 ',& k, abs((qs1+qr1+qv1+ql1+qi1)/max(abs(qv1+ql1+qi1),1.e-36_r8)) ! if (abs((qs1+qr1+qv1+ql1+qi1)/(qs1+qr1+1.e-36)).gt.1.e-14) then if (res.gt.1.e-14_r8) then call endrun ('STRATIFORM_TEND') endif ! w3 = qme(i,k) * (latvap + latice*fice(i,k)) & ! + evapheat(i,k) + prfzheat(i,k) + meltheat(i,k) res = qs1+qr1-pr1 w4 = max(abs(qs1),abs(qr1),abs(pr1)) if (w4.gt.0._r8) then if (res/w4.gt.1.e-14_r8) then write(iulog,*) ' imbalance in precips calculated two ways ' write(iulog,*) ' res/w4, pr1, qr1, qs1, qr1+qs1 ', & res/w4, pr1, qr1, qs1, qr1+qs1 ! call endrun() endif endif if (k.eq.pver) then write(iulog,*) ' pcond returned precip, rain and snow rates ', prec_pcw(i), prec_pcw(i)-snow_pcw(i), snow_pcw(i) write(iulog,*) ' I calculate ', pr1, qr1, qs1 ! call endrun write(iulog,*) ' byrons water check ', wv+wl+wi-pr1*dtime, wvf+wlf+wif endif write(iulog,*) end subroutine debug_microphys_1 !============================================================================ ! ! ! !============================================================================ ! subroutine debug_microphys_2(state1,& snow_pcw,fsaut,fsacw ,fsaci, meltheat) use ppgrid, only: pver use physconst, only: tmelt use physics_types, only: physics_state implicit none type(physics_state), intent(in) :: state1 ! local copy of the state variable real(r8), intent(in) :: snow_pcw(pcols) real(r8), intent(in) :: fsaut(pcols,pver) real(r8), intent(in) :: fsacw(pcols,pver) real(r8), intent(in) :: fsaci(pcols,pver) real(r8), intent(in) :: meltheat(pcols,pver) ! heating rate due to phase change of precip integer i,ncol,lchnk ncol = state1%ncol lchnk = state1%lchnk do i = 1,ncol if (snow_pcw(i) .gt. 0.01_r8/8.64e4_r8 .and. state1%t(i,pver) .gt. tmelt) then write(iulog,*) ' stratiform: snow, temp, ', i, lchnk, & snow_pcw(i), state1%t(i,pver) write(iulog,*) ' t ', state1%t(i,:) write(iulog,*) ' fsaut ', fsaut(i,:) write(iulog,*) ' fsacw ', fsacw(i,:) write(iulog,*) ' fsaci ', fsaci(i,:) write(iulog,*) ' meltheat ', meltheat(i,:) call endrun ('STRATIFORM_TEND') endif if (snow_pcw(i)*8.64e4_r8 .lt. -1.e-5_r8) then write(iulog,*) ' neg snow ', snow_pcw(i)*8.64e4_r8 write(iulog,*) ' stratiform: snow_pcw, temp, ', i, lchnk, & snow_pcw(i), state1%t(i,pver) write(iulog,*) ' t ', state1%t(i,:) write(iulog,*) ' fsaut ', fsaut(i,:) write(iulog,*) ' fsacw ', fsacw(i,:) write(iulog,*) ' fsaci ', fsaci(i,:) write(iulog,*) ' meltheat ', meltheat(i,:) call endrun ('STRATIFORM_TEND') endif end do end subroutine debug_microphys_2 end module rk_stratiform