CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
m_micro.F90
1
6
8module m_micro
9
10 use machine, only: kind_phys
11
12 implicit none
13 public :: m_micro_init, m_micro_run
14 private
15 logical :: is_initialized = .false.
16
17 real, parameter :: one = 1.0_kind_phys, &
18 oneb3 = one/3.0_kind_phys, &
19 zero = 0.0_kind_phys, &
20 half = 0.5_kind_phys, &
21 qsmall = 1.0e-14_kind_phys, &
22 fourb3 = 4.0_kind_phys/3.0_kind_phys, &
23 rl_cub = 1.0e-15_kind_phys, &
24 nmin = 1.0_kind_phys
25
26 real(kind=kind_phys) :: grav, pi, rgas, cp, hvap, hfus, ttp, &
27 tice, eps, epsm1, vireps, onebcp, onebg, &
28 kapa, cpbg, lvbcp, lsbcp
29
30contains
31
36subroutine m_micro_init(imp_physics, imp_physics_mg, fprcp, gravit, rair, rh2o, cpair,&
37 eps_in, epsm1_in, tmelt, latvap, latice, pi_in, tice_in, &
38 VIREPS_in, mg_dcs, mg_qcvar, mg_ts_auto_ice, &
39 mg_rhmini, microp_uniform, do_cldice, hetfrz_classnuc, &
40 mg_precip_frac_method, mg_berg_eff_factor, sed_supersat, &
41 do_sb_physics, mg_do_hail, mg_do_graupel, mg_nccons, &
42 mg_nicons, mg_ngcons, mg_ncnst, mg_ninst, mg_ngnst, &
43 mg_do_ice_gmao, mg_do_liq_liu, errmsg, errflg)
44
45 use machine, only: kind_phys
46 use cldwat2m_micro, only: ini_micro
47 use micro_mg2_0, only: micro_mg_init2_0 => micro_mg_init
48 use micro_mg3_0, only: micro_mg_init3_0 => micro_mg_init
49 use aer_cloud, only: aer_cloud_init
50
51 integer, intent(in) :: imp_physics, imp_physics_mg, fprcp
52 logical, intent(in) :: microp_uniform, do_cldice, hetfrz_classnuc, &
53 sed_supersat, do_sb_physics, mg_do_hail, &
54 mg_do_graupel, mg_nccons, mg_nicons, mg_ngcons, &
55 mg_do_ice_gmao, mg_do_liq_liu
56 real(kind=kind_phys), intent(in) :: gravit, rair, rh2o, cpair, eps_in, epsm1_in, &
57 tmelt, latvap, latice, pi_in, tice_in, vireps_in
58 real(kind=kind_phys), intent(in) :: mg_dcs, mg_qcvar, mg_ts_auto_ice(:), mg_rhmini, &
59 mg_berg_eff_factor, mg_ncnst, mg_ninst, mg_ngnst
60 character(len=16), intent(in) :: mg_precip_frac_method
61 character(len=*), intent(out) :: errmsg
62 integer, intent(out) :: errflg
63
64 errmsg = ''
65 errflg = 0
66
67 if (is_initialized) return
68
69 if (imp_physics /= imp_physics_mg) then
70 write(errmsg,'(*(a))') "Logic error: namelist choice of microphysics is different from Morrison-Gettelman MP"
71 errflg = 1
72 return
73 end if
74
75 ! Assign constants
76 grav = gravit
77 pi = pi_in
78 rgas = rair
79 cp = cpair
80 hvap = latvap
81 hfus = latice
82 ttp = tmelt
83 tice = tice_in
84 eps = eps_in
85 epsm1 = epsm1_in
86 vireps = vireps_in
87 onebcp = one/cp
88 onebg = one/grav
89 kapa = rgas*onebcp
90 cpbg = cp/grav
91 lvbcp = hvap*onebcp
92 lsbcp = (hvap+hfus)*onebcp
93
94 if (fprcp <= 0) then
95 call ini_micro (mg_dcs, mg_qcvar, mg_ts_auto_ice(1))
96 elseif (fprcp == 1) then
97 call micro_mg_init2_0(kind_phys, gravit, rair, rh2o, cpair, &
98 eps, tmelt, latvap, latice, mg_rhmini,&
99 mg_dcs, mg_ts_auto_ice, &
100 mg_qcvar, &
101 microp_uniform, do_cldice, &
102 hetfrz_classnuc, &
103 mg_precip_frac_method, &
104 mg_berg_eff_factor, &
105 sed_supersat, do_sb_physics, &
106 mg_do_ice_gmao, mg_do_liq_liu, &
107 mg_nccons, mg_nicons, &
108 mg_ncnst, mg_ninst)
109 elseif (fprcp == 2) then
110 call micro_mg_init3_0(kind_phys, gravit, rair, rh2o, cpair, &
111 eps, tmelt, latvap, latice, mg_rhmini,&
112 mg_dcs, mg_ts_auto_ice, &
113 mg_qcvar, &
114 mg_do_hail, mg_do_graupel, &
115 microp_uniform, do_cldice, &
116 hetfrz_classnuc, &
117 mg_precip_frac_method, &
118 mg_berg_eff_factor, &
119 sed_supersat, do_sb_physics, &
120 mg_do_ice_gmao, mg_do_liq_liu, &
121 mg_nccons, mg_nicons, &
122 mg_ncnst, mg_ninst, &
123 mg_ngcons, mg_ngnst)
124 else
125 write(0,*)' fprcp = ',fprcp,' is not a valid option - aborting'
126 errflg = 1
127 errmsg = 'ERROR(m_micro_init): fprcp is not a valid option'
128 return
129 endif
130 call aer_cloud_init ()
131
132 is_initialized = .true.
133
134end subroutine m_micro_init
135
145 subroutine m_micro_run( im, lm, flipv, dt_i &
146 &, prsl_i, prsi_i, phil, phii &
147 &, omega_i, QLLS_i, QLCN_i, QILS_i, QICN_i&
148 &, lwheat_i, swheat_i, w_upi, cf_upi &
149 &, FRLAND, ZPBL, CNV_MFD_i &
150 &, CNV_DQLDT_i, CLCN_i, u_i, v_i &
151 &, TAUGWX, TAUGWY &
152 &, TAUOROX, TAUOROY, CNV_FICE_i &
153 &, CNV_NDROP_i,CNV_NICE_i, q_io, lwm_o &
154 &, qi_o, t_io, rn_o, sr_o &
155 &, ncpl_io, ncpi_io, fprcp, rnw_io, snw_io&
156 &, qgl_io, ncpr_io, ncps_io, ncgl_io &
157 &, CLLS_io, KCBL, rainmin &
158 &, CLDREFFL, CLDREFFI, CLDREFFR, CLDREFFS &
159 &, CLDREFFG, ntrcaer, aerfld_i &
160 &, naai_i, npccn_i, iccn &
161 &, skip_macro &
162 &, alf_fac, qc_min, pdfflag &
163 &, kdt, xlat, xlon, rhc_i, &
164 & errmsg, errflg)
165
166! use funcphys, only: fpvs !< saturation vapor pressure for water-ice mixed
167! use funcphys, only: fpvsl, fpvsi, fpvs !< saturation vapor pressure for water,ice & mixed
173 use micro_mg2_0, only: micro_mg_tend2_0 => micro_mg_tend, qcvar2 => qcvar
174 use micro_mg3_0, only: micro_mg_tend3_0 => micro_mg_tend, qcvar3 => qcvar
175
176! use wv_saturation, only: aqsat
177
178 implicit none
179! Anning Cheng July 2015 writing the interface for GSM. Based on GMAO version of M-2M,
180! and Donifan's nuclei activation, notice the vertical coordinate is top-down
181! opposite to the GSM dynamic core, much work is still needed to consistently
182! treat other parts of the model
183! Anning Cheng 9/29/2017 implemented the MG2 from NCAR
184! alphar8 for qc_var scaled from climatology value
185!
186! Feb 2018 : S. Moorthi Updated for MG3 with graupel as prognostic variable
187!------------------------------------
188! input
189! real, parameter :: r_air = 3.47d-3
190 integer, parameter :: kp = kind_phys
191 real(kind=kind_phys), intent(in ) :: rainmin
192
193 integer, parameter :: ncolmicro = 1
194 integer,intent(in) :: im, lm, kdt, fprcp, pdfflag, iccn, ntrcaer
195 logical,intent(in) :: flipv, skip_macro
196 real (kind=kind_phys), intent(in):: dt_i, alf_fac, qc_min(:)
197
198 real (kind=kind_phys), dimension(:,:),intent(in) :: &
199 & prsl_i,u_i,v_i,phil, omega_i, qlls_i,qils_i, &
200 & lwheat_i,swheat_i
201 real (kind=kind_phys), dimension(:,0:),intent(in):: prsi_i, phii
202 real (kind=kind_phys), dimension(:,:),intent(in), optional :: &
203 & cnv_dqldt_i, clcn_i, qlcn_i, qicn_i, &
204 & cnv_mfd_i, cf_upi, cnv_fice_i, cnv_ndrop_i, &
205 & cnv_nice_i, w_upi
206! *GJF
207 real (kind=kind_phys), dimension(:,:),intent(in) :: &
208 & rhc_i, naai_i, npccn_i
209 real (kind=kind_phys), dimension(:,:,:),intent(in) :: &
210 & aerfld_i
211 real (kind=kind_phys),dimension(:),intent(in):: taugwx, &
212 & taugwy, tauorox, tauoroy, frland,zpbl,xlat,xlon
213! & TAUGWY, TAUX, TAUY, TAUOROX, TAUOROY,ps_i,FRLAND,ZPBL
214! & CNVPRCP
215
216! output
217 real (kind=kind_phys),dimension(:,:), intent(out) :: lwm_o, qi_o
218 real (kind=kind_phys),dimension(:,:), intent(out), optional :: &
219 cldreffl, cldreffi, cldreffr, cldreffs, cldreffg
220 real (kind=kind_phys),dimension(:), intent(out) :: rn_o, sr_o
221 character(len=*), intent(out) :: errmsg
222 integer, intent(out) :: errflg
223
224! input and output
225! Anning Cheng 10/24/2016 twat for total water, diagnostic purpose
226 integer, dimension(:), intent(inout):: kcbl
227 real (kind=kind_phys),dimension(:,:),intent(inout):: q_io, t_io, &
228 & ncpi_io
229 real (kind=kind_phys),dimension(:,:),intent(inout), optional :: &
230 rnw_io, snw_io, ncpr_io, ncps_io, qgl_io, ncgl_io, ncpl_io, &
231 clls_io
232! *GJF
233!Moo real (kind=kind_phys),dimension(im,lm),intent(inout):: CLLS_io
234
235
236! Local variables
237 integer kcldtopcvn,i,k,ll, kbmin, naux, nbincontactdust,l
238 integer, dimension(im) :: kct
239 real (kind=kind_phys) t_ice_all, use_av_v,bkgtau,lccirrus, &
240 & npre_frac, nct, wct, fcn, ksa1, tauxr8, dt_moist, dt_kp, tem, &
241 & tmaxll, usurf,lts_up, lts_low, min_exp, fracover, c2_gw, est3
242
243 real(kind=kind_phys), allocatable, dimension(:,:) :: &
244 & cnv_mfd, cnv_fice,cnv_ndrop,cnv_nice
245! & CNV_MFD,CNV_PRC3,CNV_FICE,CNV_NDROP,CNV_NICE
246
247 real(kind=kind_phys), dimension(IM,LM)::ncpl,ncpi,omega,sc_ice, &
248 & rad_cf, radheat,q1,u1,v1, plo, zlo, temp, &
249 & qlls, qlcn, qils,qicn, cnv_cvw,cnv_updf, &
250! & QLLS, QLCN, QILS,QICN, CNV_CVW,CNV_UPDF,SMAXL,SMAXI, &
251! & NHET_NUC, NLIM_NUC, CDNC_NUC,INC_NUC,CNN01,CNN04,CNN1,DNHET_IMM, &
252 & nhet_nuc, nlim_nuc, cdnc_nuc,inc_nuc, dnhet_imm, &
253 & nhet_imm,nhet_dep,nhet_dhf,dust_imm,dust_dep, dust_dhf,wsub, &
254 & sigw_gw,sigw_cnv,sigw_turb, &
255! & SIGW_GW,SIGW_CNV,SIGW_TURB,SIGW_RC,REV_CN_X,REV_LS_X, &
256 & rnw,snw,ncpr,ncps,qgl,ncgl, &
257! & RSU_LS_X, ALPHT_X, DLPDF_X, DIPDF_X,rnw,snw,ncpr,ncps,qgl,ncgl, &
258! & ACLL_CN_X,ACIL_CN_X, PFRZ, FQA,QCNTOT,QTOT,QL_TOT,qi_tot,blk_l,rhc
259 & fqa,ql_tot,qi_tot,blk_l,rhc
260
261 real(kind=kind_phys) :: qcntot, qtot
262
263 real(kind=kind_phys), dimension(IM,LM):: cnv_dqldt, clcn, clls
264
265! real(kind=kind_phys), dimension(IM,LM):: DQRL_X, &
266! real(kind=kind_phys), dimension(IM,LM):: CNV_DQLDT, CLCN,CLLS, &
267! & CCN01,CCN04,CCN1
268
269! real(kind=kind_phys), allocatable, dimension(:,:) :: RHX_X &
270! &, CFPDF_X, VFALLSN_CN_X, QSNOW_CN, VFALLRN_CN_X, QRAIN_CN &
271
272 real(kind=kind_phys), allocatable, dimension(:,:) :: &
273 & alpht_x, pfrz
274! & QSNOW_CN, QRAIN_CN, ALPHT_X, PFRZ
275
276! real(kind=kind_phys), allocatable, dimension(:,:) :: &
277! & QSNOW_CN, QRAIN_CN &
278!! &, CFPDF_X, QSNOW_CN, QRAIN_CN &
279! &, ALPHT_X, PFRZ
280!! &, REV_CN_X, RSU_CN_X, DLPDF_X, DIPDF_X, ALPHT_X, PFRZ &
281!! &, ACLL_CN_X, ACIL_CN_X, DQRL_X &
282!! &, PFI_CN_X, PFL_CN_X, QST3, DZET, QDDF3
283!! real(kind=kind_phys), allocatable, dimension(:) :: vmip
284
285! real(kind=kind_phys), dimension(IM,LM) :: QDDF3
286! real(kind=kind_phys), dimension(IM,LM):: QST3, DZET, QDDF3
287! & MASS, RHX_X, CFPDF_X, &
288! & VFALLSN_CN_X, QSNOW_CN, &
289! & VFALLRN_CN_X, QRAIN_CN
290! & VFALLRN_CN_X, QRAIN_CN, dum
291
292 real(kind=kind_phys), dimension(IM,LM+1) :: zet
293 real(kind=kind_phys), dimension(IM,0:LM) :: ple, kh
294
295! real(kind=kind_phys), dimension(IM,0:LM) :: PLE, PKE, kh
296! &, PFI_CN_X, PFL_CN_X
297
298 real(kind=kind_phys),dimension(LM) :: rhdfdar8, rhu00r8, &
299 & ttendr8,qtendr8, cwtendr8,npre8, npccninr8,ter8, &
300 & plevr8,ndropr8,qir8,qcr8,wparc_turb,qvr8, nir8,ncr8, &
301 & nimmr8,nsootr8,rnsootr8,omegr8,qrr8,qsr8,nrr8,nsr8, &
302 & qgr8, ngr8
303
304 real(kind=kind_phys), dimension(1:LM,10) :: rndstr8,naconr8
305
306! real(kind=kind_phys), dimension(IM) :: CN_PRC2,CN_SNR,CN_ARFX,&
307! & LS_SNR, LS_PRC2, TPREC
308 real(kind=kind_phys), dimension(IM) :: ls_snr, ls_prc2
309! & VMIP, twat
310
311 real(kind=kind_phys), dimension (LM) :: uwind_gw,vwind_gw, &
312 & tm_gw, pm_gw, nm_gw, h_gw, rho_gw, khaux, qcaux, &
313 & dummyw , wparc_cgw, cfaux, dpre8, &
314 & wparc_ls,wparc_gw, swparc,smaxliq,smaxicer8,nheticer8, &
315 & nhet_immr8,dnhet_immr8,nhet_depr8,nhet_dhfr8,sc_icer8, &
316 & dust_immr8,dust_depr8,dust_dhfr8,nlimicer8,cldfr8,liqcldfr8, &
317 & icecldfr8,cldor8, pdelr8, &
318 & rpdelr8,lc_turb,zmr8,ficer8,rate1ord_cw2pr, tlatr8, qvlatr8, &
319 & qctendr8, qitendr8, nctendr8, nitendr8, effcr8, effc_fnr8, &
320 & effir8, nevaprr8, evapsnowr8, prainr8, &
321 & prodsnowr8, cmeoutr8, deffir8, pgamradr8, lamcradr8,qsoutr8, &
322 & qroutr8,droutr8, qcsevapr8,qisevapr8, qvresr8, &
323 & cmeioutr8, dsoutr8, qcsinksum_rate1ord,qrtend,nrtend, &
324 & qstend, nstend, alphar8, rhr8, &
325
326 & qgtend, ngtend, qgoutr8, ngoutr8, dgoutr8
327
328 real(kind=kind_phys), dimension(1) :: prectr8, precir8
329
330 real(kind=kind_phys), dimension (LM) :: vtrmcr8,vtrmir8, &
331 & qcsedtenr8,qisedtenr8, praor8,prcor8,mnucccor8, mnucctor8, &
332 & msacwior8,psacwsor8, bergsor8,bergor8,meltor8, homoor8, &
333 & qcresor8, &
334 & prcior8, praior8,qiresor8, mnuccror8,pracsor8, meltsdtr8, &
335 & frzrdtr8, &
336 & ncalr8, ncair8, mnuccdor8, nnucctor8, nsoutr8, nroutr8, &
337 & nnuccdor8, nnucccor8,naair8, &
338 & nsacwior8, nsubior8, nprcior8, npraior8, npccnor8, npsacwsor8, &
339 & nsubcor8, npraor8, nprc1or8, tlatauxr8,pfrz_inc_kp,sadice, &
340 & sadsnow, am_evp_st, reff_rain, reff_snow, &
341 & umr,ums,qrsedten,qssedten,refl,arefl,areflz,frefl,csrfl, &
342 & acsrfl,fcsrfl,rercld,qrout2,qsout2,nrout2,nsout2,drout2, &
343 & dsout2,freqs,freqr,nfice,qcrat,prer_evap, &
344! graupel related
345 & reff_grau, umg, qgsedtenr8, mnuccrior8, &
346 & pracgr8, psacwgr8, pgsacwr8, pgracsr8, prdgr8, qmultgr8,&
347 & qmultrgr8, psacrr8, npracgr8, nscngr8, ngracsr8, nmultgr8,&
348 & nmultrgr8, npsacwgr8, qgout2, ngout2, dgout2, freqg
349
350 real(kind=kind_phys), dimension (0:LM) :: pi_gw, rhoi_gw, &
351 & ni_gw, ti_gw
352
353 real(kind=kind_phys), dimension(LM+1) :: pintr8, kkvhr8
354 real(kind=kind_phys), dimension(2:LM+1) :: lflx, iflx, rflx, &
355 sflx, gflx
356
357! real (kind=kind_phys), parameter :: disp_liu=2., ui_scale=1.0 &
358! &, dcrit=20.0d-6 &
359 real (kind=kind_phys), parameter :: disp_liu=1.0_kp &
360 &, ui_scale=1.0_kp &
361 &, dcrit=1.0e-6_kp &
362! &, ts_autice=1800.0 &
363! &, ts_autice=3600.0 & !time scale
364 &, ninstr8 = 0.1e6_kp &
365 &, ncnstr8 = 100.0e6_kp
366
367 real(kind=kind_phys):: k_gw, maxkh, tausurf_gw, overscale, tx1, rh1_kp
368 real(kind=kind_phys):: t_ice_denom
369
370 integer, dimension(1) :: lev_sed_strt ! sedimentation start level
371 real(kind=kind_phys), parameter :: sig_sed_strt=0.05_kp ! normalized pressure at sedimentation start
372
373 real(kind=kind_phys),dimension(3) :: ccn_diag
374 real(kind=kind_phys),dimension(58) :: cloudparams
375
376 integer, parameter :: ccn_param=2, in_param=5
377
378 real(kind=kind_phys), parameter ::fdust_drop=1.0_kp, fsoot_drop=0.1_kp &
379 &, sigma_nuc_kp=0.28_kp,sclmfdfr=0.03_kp
380! &, sigma_nuc_kp=0.28,SCLMFDFR=0.1
381
382 type (aerprops), dimension (IM,LM) :: aeroprops
383 type (aerprops) :: aeroaux, aeroaux_b
384 real, allocatable, dimension(:,:,:) :: aermassmix
385
386 logical :: use_average_v, ltrue, lprint, lprnt
387 integer :: ipr
388
389!==================================
390!====2-moment Microhysics=
391!================== Start Stratiform cloud processes==========================================
392!set up initial values
393
394 data use_av_v/1.0_kp/, bkgtau/0.015_kp/, lccirrus/500.0_kp/, npre_frac/1.0_kp/, &
395 & tmaxll/296.0_kp/, fracover/1.0_kp/, lts_low/12.0_kp/, lts_up/24.0_kp/, &
396 & min_exp/0.5_kp/
397
398 data cloudparams/ &
399 & 10.0_kp, 4.0_kp , 4.0_kp , 1.0_kp , 2.e-3_kp, 8.e-4_kp, 2.0_kp , 1.0_kp , -1.0_kp &
400 &, 0.0_kp , 1.3_kp , 1.0e-9_kp, 3.3e-4_kp, 20.0_kp , 4.8_kp , 4.8_kp , 230.0_kp , 1.0_kp &
401 &, 1.0_kp , 230.0_kp, 14400._kp, 50.0_kp , 0.01_kp , 0.1_kp , 200.0_kp, 0.0_kp , 0.0_kp &
402 &, 0.5_kp , 0.5_kp , 2000.0_kp, 0.8_kp , 0.5_kp , -40.0_kp, 1.0_kp , 4.0_kp , 0.0_kp &
403 &, 0.0_kp , 0.0_kp , 1.0e-3_kp, 8.0e-4_kp, 1.0_kp , 0.95_kp , 1.0_kp , 0.0_kp , 900.0_kp&
404! &, 0.0_kp , 0.0_kp , 1.0e-3_kp, 8.0e-4_kp, 1.0_kp , 0.95_kp , 1.0_kp , 0.0_kp , 880.0_kp&
405! &, 0.0_kp , 0.0_kp , 1.0e-3_kp, 8.0e-4_kp, 1.0_kp , 0.95_kp , 1.0_kp , 0.0_kp , 980.0_kp&
406 &, 1.0_kp , 1.0_kp , 1.0_kp , 0.0_kp , 0.0_kp , 1.e-5_kp, 2.e-5_kp, 2.1e-5_kp, 4.e-5_kp&
407! &, 3e-5_kp, 0.1_kp , 4.0_kp , 250.0_kp/ ! Annings version
408 &, 3e-5_kp, 0.1_kp , 4.0_kp , 150.0_kp/ ! Annings version
409! &, 3e-5_kp, 0.1_kp , 1.0_kp , 150.0_kp/
410
411! Initialize CCPP error handling variables
412 errmsg = ''
413 errflg = 0
414
415 lprnt = .false.
416 ipr = 1
417
418! rhr8 = 1.0
419
420 if(flipv) then
421 DO k=1, lm
422 ll = lm-k+1
423 DO i = 1,im
424 q1(i,k) = q_io(i,ll)
425 u1(i,k) = u_i(i,ll)
426 v1(i,k) = v_i(i,ll)
427 omega(i,k) = omega_i(i,ll)
428 ncpl(i,k) = ncpl_io(i,ll)
429 ncpi(i,k) = ncpi_io(i,ll)
430 rnw(i,k) = rnw_io(i,ll)
431 snw(i,k) = snw_io(i,ll)
432 qgl(i,k) = qgl_io(i,ll)
433 ncpr(i,k) = ncpr_io(i,ll)
434 ncps(i,k) = ncps_io(i,ll)
435 ncgl(i,k) = ncgl_io(i,ll)
436! QLLS is the total cloud water
437 qlls(i,k) = qlls_i(i,ll)-qlcn_i(i,ll)
438 qlcn(i,k) = qlcn_i(i,ll)
439 qils(i,k) = qils_i(i,ll)-qicn_i(i,ll)
440 qicn(i,k) = qicn_i(i,ll)
441 cnv_cvw(i,k) = w_upi(i,ll)
442 cnv_updf(i,k) = cf_upi(i,ll)
443 cnv_dqldt(i,k) = cnv_dqldt_i(i,ll)
444 clcn(i,k) = clcn_i(i,ll)
445 clls(i,k) = max(clls_io(i,ll)-clcn_i(i,ll),zero)
446 plo(i,k) = prsl_i(i,ll)*0.01_kp
447 zlo(i,k) = phil(i,ll) * onebg
448 temp(i,k) = t_io(i,ll)
449 radheat(i,k) = lwheat_i(i,ll) + swheat_i(i,ll)
450 rhc(i,k) = rhc_i(i,ll)
451 if (iccn == 1) then
452 cdnc_nuc(i,k) = npccn_i(i,ll)
453 inc_nuc(i,k) = naai_i(i,ll)
454 endif
455
456 END DO
457 END DO
458 DO k=0, lm
459 ll = lm-k
460 DO i = 1,im
461 ple(i,k) = prsi_i(i,ll) * 0.01_kp ! interface pressure in hPa
462 zet(i,k+1) = phii(i,ll) * onebg
463 END DO
464 END DO
465 if (.not. skip_macro) then
466! allocate(CNV_MFD(im,lm), CNV_PRC3(im,lm), CNV_FICE(im,lm) &
467 allocate(cnv_mfd(im,lm), cnv_fice(im,lm) &
468 &, cnv_ndrop(im,lm), cnv_nice(im,lm))
469 DO k=1, lm
470 ll = lm-k+1
471 DO i = 1,im
472 cnv_mfd(i,k) = cnv_mfd_i(i,ll)
473! CNV_PRC3(i,k) = CNV_PRC3_i(i,ll)
474 cnv_fice(i,k) = cnv_fice_i(i,ll)
475 cnv_ndrop(i,k) = cnv_ndrop_i(i,ll)
476 cnv_nice(i,k) = cnv_nice_i(i,ll)
477 enddo
478 enddo
479 endif
480
481 else
482 DO k=1, lm
483 DO i = 1,im
484 q1(i,k) = q_io(i,k)
485 u1(i,k) = u_i(i,k)
486 v1(i,k) = v_i(i,k)
487 omega(i,k) = omega_i(i,k)
488 ncpl(i,k) = ncpl_io(i,k)
489 ncpi(i,k) = ncpi_io(i,k)
490 rnw(i,k) = rnw_io(i,k)
491 snw(i,k) = snw_io(i,k)
492 qgl(i,k) = qgl_io(i,k)
493 ncpr(i,k) = ncpr_io(i,k)
494 ncps(i,k) = ncps_io(i,k)
495 ncgl(i,k) = ncgl_io(i,k)
496! QLLS is the total cloud water
497 qlls(i,k) = qlls_i(i,k)-qlcn_i(i,k)
498 qlcn(i,k) = qlcn_i(i,k)
499 qils(i,k) = qils_i(i,k)-qicn_i(i,k)
500 qicn(i,k) = qicn_i(i,k)
501 cnv_cvw(i,k) = w_upi(i,k)
502 cnv_updf(i,k) = cf_upi(i,k)
503 cnv_dqldt(i,k) = cnv_dqldt_i(i,k)
504 clcn(i,k) = clcn_i(i,k)
505 clls(i,k) = max(clls_io(i,k)-clcn_i(i,k),zero)
506 plo(i,k) = prsl_i(i,k)*0.01_kp
507 zlo(i,k) = phil(i,k) * onebg
508 temp(i,k) = t_io(i,k)
509 radheat(i,k) = lwheat_i(i,k) + swheat_i(i,k)
510 rhc(i,k) = rhc_i(i,k)
511 if (iccn == 1) then
512 cdnc_nuc(i,k) = npccn_i(i,k)
513 inc_nuc(i,k) = naai_i(i,k)
514 endif
515
516 END DO
517 END DO
518 DO k=0, lm
519 DO i = 1,im
520 ple(i,k) = prsi_i(i,k) * 0.01_kp ! interface pressure in hPa
521 zet(i,k+1) = phii(i,k) * onebg
522 END DO
523 END DO
524 if (.not. skip_macro) then
525! allocate(CNV_MFD(im,lm), CNV_PRC3(im,lm), CNV_FICE(im,lm) &
526 allocate(cnv_mfd(im,lm), cnv_fice(im,lm) &
527 &, cnv_ndrop(im,lm), cnv_nice(im,lm))
528 DO k=1, lm
529 DO i = 1,im
530 cnv_mfd(i,k) = cnv_mfd_i(i,k)
531! CNV_PRC3(i,k) = CNV_PRC3_i(i,k)
532 cnv_fice(i,k) = cnv_fice_i(i,k)
533 cnv_ndrop(i,k) = cnv_ndrop_i(i,k)
534 cnv_nice(i,k) = cnv_nice_i(i,k)
535 enddo
536 enddo
537 endif
538 endif
539! if (lprnt) then
540! write(0,*)' inmic qlcn=',qlcn(ipr,:)
541! write(0,*)' inmic qlls=',qlls(ipr,:)
542! write(0,*)' inmic qicn=',qicn(ipr,:)
543! write(0,*)' inmic qils=',qils(ipr,:)
544! endif
545!
546 dt_moist = dt_i
547 dt_kp = dt_i
548
549 if (kdt == 1) then
550 DO k=1, lm
551 DO i = 1,im
552 CALL fix_up_clouds_2m(q1(i,k), temp(i,k), qlls(i,k), &
553 & qils(i,k), clls(i,k), qlcn(i,k), &
554 & qicn(i,k), clcn(i,k), ncpl(i,k), &
555 & ncpi(i,k), qc_min)
556 if (rnw(i,k) <= qc_min(1)) then
557 ncpr(i,k) = zero
558 rnw(i,k) = zero
559 elseif (ncpr(i,k) <= nmin) then ! make sure NL > 0 if Q >0
560 ncpr(i,k) = max(rnw(i,k) / (fourb3 * pi *rl_cub*997.0_kp), nmin)
561 endif
562 if (snw(i,k) <= qc_min(2)) then
563 ncps(i,k) = zero
564 snw(i,k) = zero
565 elseif (ncps(i,k) <= nmin) then
566 ncps(i,k) = max(snw(i,k) / (fourb3 * pi *rl_cub*500.0_kp), nmin)
567 endif
568 if (qgl(i,k) <= qc_min(2)) then
569 ncgl(i,k) = zero
570 qgl(i,k) = zero
571 elseif (ncgl(i,k) <= nmin) then
572 ncgl(i,k) = max(qgl(i,k) / (fourb3 * pi *rl_cub*500.0_kp), nmin)
573 endif
574
575 enddo
576 enddo
577 endif
578
579 do i=1,im
580 kcbl(i) = max(lm-kcbl(i),10)
581 kct(i) = 10
582 enddo
583
584 DO i=1, im
585 DO k = lm-2, 10, -1
586 If ((cnv_dqldt(i,k) <= 1.0e-9_kp) .and. &
587 & (cnv_dqldt(i,k+1) > 1.0e-9_kp)) then
588 kct(i) = k+1
589 exit
590 end if
591 END DO
592 END DO
593
594! do L=LM,1,-1
595! do i=1,im
596! DZET(i,L) = ZET(i,L) - ZET(i,L+1)
597! tx1 = plo(i,l)*100.0
598! est3 = min(tx1, fpvs(temp(i,l)))
599! qst3(i,l) = min(eps*est3/max(tx1+epsm1*est3,1.0e-10),1.0)
600! MASS(i,l) = (ple(i,l) - ple(i,l-1)) * (100.0/grav)
601! enddo
602! enddo
603!------------------------------------------------------------------------------
604! call aqsat(temp,plo*100.,est3,qst3,im,im,lm,1,lm)
605! do k=1,lm
606! do i=1,im
607! DZET(i,k) = TH1(i,k) * (pke(i,k)-pke(i,k-1)) &
608! & * cpbg * (1.0 + vireps*q1(i,k))
609! MASS(i,k) = (ple(i,k) - ple(i,k-1)) * (100.0/grav)
610! end do
611! end do
612
613! do k=1,lm
614! do i=1,im
615! temp(i,k) = th1(i,k) * PK(i,k)
616! est3 = fpvs(temp(i,k))
617! qst3(i,k) = min(eps*est3/max(plo(i,k)*100.0+epsm1*est3,1.0e-10),1.0)
618! enddo
619! enddo
620! call aqsat(temp,plo*100.,est3,qst3,im,im,lm,1,lm)
621! do k=1,lm
622! do i=1,im
623! DZET(i,k) = TH1(i,k) * (pke(i,k)-pke(i,k-1)) &
624! & * cpbg * (1.0 + vireps*q1(i,k))
625! MASS(i,k) = (ple(i,k) - ple(i,k-1)) * (100.0/grav)
626! enddo
627! enddo
628
629! do i=1,im
630! ZET(i,LM+1) = 0.0
631! vmip(i) = 0.0
632! enddo
633!------------------------------------------------------------------------------
634
635! if (.not. skip_macro) then
636! allocate(qddf3(im,lm))
637! allocate(vmip(im))
638! do i=1,im
639! vmip(i) = 0.0
640! enddo
641! DO K = LM, 1, -1
642! do i=1,im
643! if (zet(i,k) < 3000.0) then
644!! qddf3(i,k) = - (zet(i,k) - 3000.0) * zet(i,k) * mass(i,k)
645! qddf3(i,k) = - (zet(i,k) - 3000.0) * zet(i,k) &
646! & * (ple(i,k) - ple(i,k-1)) * (100.0/grav)
647! else
648! qddf3(i,k) = 0.0
649! endif
650! vmip(i) = vmip(i) + qddf3(i,k)
651! enddo
652! END DO
653! do i=1,im
654! if (vmip(i) /= 0.0) vmip(i) = 1.0 / vmip(i)
655! enddo
656! DO K = 1,LM
657! do i=1,im
658! QDDF3(i,K) = QDDF3(i,K) * VMIP(i)
659! enddo
660! END DO
661! deallocate (vmip)
662! endif
663
664 do l=lm-1,1,-1
665 do i=1,im
666 tx1 = half * (temp(i,l+1) + temp(i,l))
667 kh(i,l) = 3.55e-7_kp*tx1**2.5_kp*(rgas*0.01_kp) / ple(i,l) !kh molecule diff only needing refinement
668 enddo
669 end do
670 do i=1,im
671 kh(i,0) = kh(i,1)
672 kh(i,lm) = kh(i,lm-1)
673 enddo
674 do l=lm,1,-1
675 do i=1,im
676 blk_l(i,l) = one / ( one/max(0.15_kp*zpbl(i),0.4_kp*zlo(i,lm-1))&
677 & + one/(zlo(i,l)*0.4_kp) )
678
679 sc_ice(i,l) = one
680 ncpl(i,l) = max( ncpl(i,l), zero)
681 ncpi(i,l) = max( ncpi(i,l), zero)
682 rad_cf(i,l) = max(zero, min(clls(i,l)+clcn(i,l), one))
683 if (iccn /= 1) then
684 cdnc_nuc(i,l) = zero
685 inc_nuc(i,l) = zero
686 endif
687
688 enddo
689 end do
690! T_ICE_ALL = TICE - 40.0
691 t_ice_all = cloudparams(33) + tice
692 t_ice_denom = one / (tice - t_ice_all)
693
694
695 do l=1,lm
696 rhdfdar8(l) = 1.e-8_kp
697 rhu00r8(l) = 0.95_kp
698
699 ttendr8(l) = zero
700 qtendr8(l) = zero
701 cwtendr8(l) = zero
702
703 npccninr8(l) = zero
704 enddo
705 do k=1,10
706 do l=1,lm
707 rndstr8(l,k) = 2.0e-7_kp
708 enddo
709 enddo
710
711!need an estimate of convective area
712!=======================================================================================================================
713!=======================================================================================================================
719!=======================================================================================================================
720!=======================================================================================================================
721!=======================================================================================================================
722! if(aero_in) then
723! allocate(AERMASSMIX (IM,LM, 15))
724! AERMASSMIX = 1.e-15
725! call AerConversion1 (AERMASSMIX, AeroProps)
726! deallocate(AERMASSMIX)
727! end if
728
729!
731 do k=1,lm
732 do i=1,im
733 call init_aer(aeroprops(i, k))
734 enddo
735 enddo
736!
737
738 allocate(aermassmix(im,lm,15))
739 if ( iccn == 2) then
740 aermassmix(:,:,1:ntrcaer) = aerfld_i(:,:,1:ntrcaer)
741 else
742 aermassmix(:,:,1:5) = 1.0e-6_kp
743 aermassmix(:,:,6:15) = 2.0e-14_kp
744 end if
746 call aerconversion1 (aermassmix, aeroprops)
747 deallocate(aermassmix)
748
749 use_average_v = .false.
750 if (use_av_v > zero) then
751 use_average_v = .true.
752 end if
753
754 k_gw = (pi+pi) / lccirrus
755
756!-------------------------------------------------------------------------------
757 do i=1,im ! beginning of first big I loop
758
759 kcldtopcvn = kct(i)
760
761 tausurf_gw = min(half*sqrt(tauorox(i)*tauorox(i) &
762 & + tauoroy(i)*tauoroy(i)), 10.0_kp)
763 do k=1,lm
764
765 uwind_gw(k) = min(half*sqrt( u1(i,k)*u1(i,k) &
766 & + v1(i,k)*v1(i,k)), 50.0_kp)
767
768! tausurf_gw =tausurf_gw + max (tausurf_gw, min(0.5*SQRT(TAUX(I , J)**2+TAUY(I , J)**2), 10.0)*BKGTAU) !adds a minimum value from unresolved sources
769
770
771 pm_gw(k) = 100.0_kp*plo(i,k)
772 tm_gw(k) = temp(i,k)
773
774 nm_gw(k) = zero
775 rho_gw(k) = pm_gw(k) /(rgas*tm_gw(k))
776
777 ter8(k) = temp(i,k)
778 plevr8(k) = 100.0_kp*plo(i,k)
779 ndropr8(k) = ncpl(i,k)
780 qir8(k) = qils(i,k) + qicn(i,k)
781 qcr8(k) = qlls(i,k) + qlcn(i,k)
782 qcaux(k) = qcr8(k)
783
784 npccninr8(k) = zero
785 naair8(k) = zero
786
787 npre8(k) = zero
788
789 if (rad_cf(i,k) > 0.01_kp .and. qir8(k) > zero) then
790 npre8(k) = npre_frac*ncpi(i,k)
791 else
792 npre8(k) = zero
793 endif
794
795 omegr8(k) = omega(i,k)
796 lc_turb(k) = max(blk_l(i,k), 50.0_kp)
797! rad_cooling(k) = RADheat(I,k)
798
799 if (npre8(k) > zero .and. qir8(k) > zero) then
800 dpre8(k) = ( qir8(k)/(6.0_kp*npre8(k)*900.0_kp*pi))**(one/3.0_kp)
801 else
802 dpre8(k) = 1.0e-9_kp
803 endif
804 wparc_ls(k) = -omegr8(k) / (rho_gw(k)*grav) &
805 & + cpbg * radheat(i,k)
806! & + cpbg * rad_cooling(k)
807 enddo
808 do k=0,lm
809 pi_gw(k) = 100.0_kp*ple(i,k)
810 rhoi_gw(k) = zero
811 ni_gw(k) = zero
812 ti_gw(k) = zero
813 enddo
814
815
816! ====================================================================
818! ====================================================================
819
820
821 call gw_prof (1, lm, 1, tm_gw, pm_gw, pi_gw, rhoi_gw, ni_gw, &
822 & ti_gw, nm_gw, q1(i,:))
823
824 do k=1,lm
825 nm_gw(k) = max(nm_gw(k), 0.005_kp)
826 h_gw(k) = k_gw*rho_gw(k)*uwind_gw(k)*nm_gw(k)
827 if (h_gw(k) > zero) then
828 h_gw(k) = sqrt(2.0_kp*tausurf_gw/h_gw(k))
829 end if
830
831 wparc_gw(k) = k_gw*uwind_gw(k)*h_gw(k)*0.133_kp
832
833 wparc_cgw(k) = zero
834 end do
835
837
838 if (kcldtopcvn > 20) then
839
840 ksa1 = one
841 nct = nm_gw(kcldtopcvn)
842 wct = max(cnv_cvw(i,kcldtopcvn), zero)
843
844 fcn = maxval(cnv_updf(i,kcldtopcvn:lm))
845
846 do k=1,kcldtopcvn
847 c2_gw = (nm_gw(k) + nct) / nct
848 wparc_cgw(k) = sqrt(ksa1*fcn*fcn*12.56_kp* &
849 & 1.806_kp*c2_gw*c2_gw)*wct*0.133_kp
850 enddo
851
852 end if
853
854 do k=1,lm
855 dummyw(k) = 0.133_kp*k_gw*uwind_gw(k)/nm_gw(k)
856 enddo
857
858 do k=1, lm-5, 1
859 if (wparc_cgw(k)+wparc_gw(k) > dummyw(k)) then
860 exit
861 end if
862 end do
863
864 do l=1,min(k,lm-5)
865 wparc_cgw(l) = zero
866 wparc_gw(l) = zero
867 enddo
868
869
870
871 kbmin = kcbl(i)
872 kbmin = min(kbmin, lm-1) - 4
873 do k = 1, lm
874 wparc_turb(k) = kh(i,k) / lc_turb(k)
875 dummyw(k) = 10.0_kp
876 enddo
877
878 if (frland(i) < 0.1_kp .and. zpbl(i) < 800.0_kp .and. &
879 & temp(i,lm) < 298.0_kp .and. temp(i,lm) > 274.0_kp) then
880 do k = 1, lm
881 dummyw(k) = max(min((zet(i,k+1)-zpbl(i))*0.01_kp, 10.0_kp),-10.0_kp)
882 dummyw(k) = one / (one+exp(dummyw(k)))
883 enddo
884 maxkh = max(maxval(kh(i,kbmin:lm-1)*nm_gw(kbmin:lm-1)/ &
885 & 0.17_kp), 0.3_kp)
886 do k = 1, lm
887 wparc_turb(k) = (one-dummyw(k))*wparc_turb(k) &
888 & + dummyw(k)*maxkh
889 enddo
890
891 end if
892
893 wparc_turb(kbmin:lm) = max(wparc_turb(kbmin:lm), 0.2_kp)
894
895
896
898
899 do k = 1, lm
900 swparc(k) = sqrt(wparc_gw(k) * wparc_gw(k) &
901 & + wparc_turb(k) * wparc_turb(k) &
902 & + wparc_cgw(k) * wparc_cgw(k))
903 enddo
904
905
906! ==========================================================================================
907! ========================Activate the aerosols ============================================
908
909 do k = 1, lm
910
911 if (plevr8(k) > 70.0_kp) then
912
913 ccn_diag(1) = 0.001_kp
914 ccn_diag(2) = 0.004_kp
915 ccn_diag(3) = 0.01_kp
916
917 if (k > 2 .and. k <= lm-2) then
918 tauxr8 = (ter8(k-1) + ter8(k+1) + ter8(k)) * oneb3
919 else
920 tauxr8 = ter8(k)
921 endif
922
923! if(aero_in) then
924 aeroaux = aeroprops(i, k)
925! else
926! call init_Aer(AeroAux)
927! call init_Aer(AeroAux_b)
928! endif
929
930 pfrz_inc_kp(k) = zero
931 rh1_kp = zero !related to cnv_dql_dt, needed to changed soon
932
933! if (lprnt) write(0,*)' bef aero npccninr8=',npccninr8(k),' k=',k &
934! &,' ccn_param=',ccn_param,' in_param=',in_param &
935! &,' AeroAux%kap=',AeroAux%kap
936
938 call aerosol_activate(tauxr8, plevr8(k), swparc(k), &
939 & wparc_ls(k), aeroaux, npre8(k), dpre8(k), ccn_diag, &
940 & ndropr8(k), npccninr8(k), smaxliq(k), &
941! & ndropr8(k), qcr8(K), npccninr8(K), smaxliq(K), &
942 & naair8(k), smaxicer8(k), nheticer8(k), nhet_immr8(k), &
943 & dnhet_immr8(k), nhet_depr8(k), nhet_dhfr8(k), &
944 & sc_icer8(k), dust_immr8(k), dust_depr8(k), &
945 & dust_dhfr8(k), nlimicer8(k), use_average_v, &
946 & ccn_param, in_param, fdust_drop, &
947 & fsoot_drop,pfrz_inc_kp(k),sigma_nuc_kp, rh1_kp, &
948 & size(ccn_diag))
949! & size(ccn_diag), lprnt)
950! if (lprnt) write(0,*)' aft aero npccninr8=',npccninr8(k),' k=',k
951
952 if (npccninr8(k) < 1.0e-12_kp) npccninr8(k) = zero
953
954! CCN01(I,K) = max(ccn_diag(1)*1e-6, 0.0)
955! CCN04(I,K) = max(ccn_diag(2)*1e-6, 0.0)
956! CCN1 (I,K) = max(ccn_diag(3)*1e-6, 0.0)
957
958
959
960 else
961 ccn_diag(:) = zero
962 smaxliq(k) = zero
963 swparc(k) = zero
964 smaxicer8(k) = zero
965 nheticer8(k) = zero
966 sc_icer8(k) = 2.0_kp
967! sc_icer8(K) = 1.0d0
968 naair8(k) = zero
969 npccninr8(k) = zero
970 nlimicer8(k) = zero
971 nhet_immr8(k) = zero
972 dnhet_immr8(k) = zero
973 nhet_depr8(k) = zero
974 nhet_dhfr8(k) = zero
975 dust_immr8(k) = zero
976 dust_depr8(k) = zero
977 dust_dhfr8(k) = zero
978
979 end if
980
981! SMAXL(I,k) = smaxliq(k) * 100.0
982! SMAXI(I,k) = smaxicer8(k) * 100.0
983 nhet_nuc(i,k) = nheticer8(k) * 1.0e-6_kp
984 nlim_nuc(i,k) = nlimicer8(k) * 1.0e-6_kp
985 sc_ice(i,k) = min(max(sc_icer8(k),one),2.0_kp)
986! SC_ICE(I,k) = min(max(sc_icer8(k),1.0),1.2)
987! if(temp(i,k) < T_ICE_ALL) SC_ICE(i,k) = max(SC_ICE(I,k), 1.2)
988! if(temp(i,k) < T_ICE_ALL) SC_ICE(i,k) = max(SC_ICE(I,k), 1.5)
989! if(temp(i,k) > T_ICE_ALL) SC_ICE(i,k) = 1.0
990! if(temp(i,k) > TICE) SC_ICE(i,k) = rhc(i,k)
991!
992 if (iccn == 0) then
993 if(temp(i,k) < t_ice_all) then
994! SC_ICE(i,k) = max(SC_ICE(I,k), 1.2)
995 sc_ice(i,k) = max(sc_ice(i,k), 1.5_kp)
996 elseif(temp(i,k) > tice) then
997 sc_ice(i,k) = rhc(i,k)
998 else
999! SC_ICE(i,k) = 1.0
1000! tx1 = max(SC_ICE(I,k), 1.2)
1001 tx1 = max(sc_ice(i,k), 1.5_kp)
1002 sc_ice(i,k) = ((tice-temp(i,k))*tx1 + (temp(i,k)-t_ice_all)*rhc(i,k)) &
1003 * t_ice_denom
1004 endif
1005 endif
1006 if (iccn /= 1) then
1007 cdnc_nuc(i,k) = npccninr8(k)
1008 inc_nuc(i,k) = naair8(k)
1009 endif
1010 nhet_imm(i,k) = max(nhet_immr8(k), zero)
1011 dnhet_imm(i,k) = max(dnhet_immr8(k), zero)
1012 nhet_dep(i,k) = nhet_depr8(k) * 1.0e-6_kp
1013 nhet_dhf(i,k) = nhet_dhfr8(k) * 1.0e-6_kp
1014 dust_imm(i,k) = max(dust_immr8(k), zero)*1.0e-6_kp
1015 dust_dep(i,k) = max(dust_depr8(k), zero)*1.0e-6_kp
1016 dust_dhf(i,k) = max(dust_dhfr8(k), zero)*1.0e-6_kp
1017 wsub(i,k) = wparc_ls(k) + swparc(k)*0.8_kp
1018 sigw_gw(i,k) = wparc_gw(k) * wparc_gw(k)
1019 sigw_cnv(i,k) = wparc_cgw(k) * wparc_cgw(k)
1020 sigw_turb(i,k) = wparc_turb(k) * wparc_turb(k)
1021
1022 enddo ! end of K loop
1023 enddo ! end of first big I loop
1024!-------------------------------------------------------------------------------
1025
1026! SC_ICE=MIN(MAX(SC_ICE, 1.0), 2.0)
1027! WHERE (TEMP .gt. T_ICE_ALL)
1028! SC_ICE=1.0
1029! END WHERE
1030
1031!===========================End cloud particle nucleation=======================
1032! -----------------------------
1033!
1035
1036! do k=1,lm
1037! do i=1,im
1038! REV_CN_X(i,k) = 0.0
1039! REV_LS_X(i,k) = 0.0
1040! RSU_CN_X(i,k) = 0.0
1041! RSU_LS_X(i,k) = 0.0
1042! CFX(i,k) = INC_NUC(i,k) + NHET_IMM(i,k)
1043! enddo
1044! enddo
1045! do k=0,lm
1046! do i=1,im
1047! PFI_CN_X(i,k) = 0.0
1048! PFL_CN_X(i,k) = 0.0
1049! enddo
1050! enddo
1051
1052! if(lprnt) write(0,*)' skip_macro=',skip_macro
1053
1054 if (.not. skip_macro) then
1055
1056! if (lprnt) write(0,*) ' in micro qicn2=',qicn(ipr,25),' kdt=',kdt&
1057! &,' qils=',qils(ipr,25)
1058! if(lprnt) write(0,*)' bef macro_cloud clcn=',clcn(ipr,:)
1059! if(lprnt) write(0,*)' bef macro_cloud clls=',clls(ipr,:)
1060
1061! allocate(RHX_X(im,lm), CFPDF_X(im,lm), VFALLSN_CN_X(im,lm), &
1062 allocate( &
1063! & QSNOW_CN(im,lm), VFALLRN_CN_X(im,lm), QRAIN_CN(im,lm),&
1064! & QSNOW_CN(im,lm), QRAIN_CN(im,lm),&
1065! & REV_CN_X(im,lm), RSU_CN_X(im,lm), DLPDF_X(im,lm), &
1066! & DIPDF_X(im,lm), ALPHT_X(im,lm), PFRZ(im,lm), &
1067 & alpht_x(im,lm), pfrz(im,lm))
1068! & ACLL_CN_X(im,lm), ACIL_CN_X(im,lm), DQRL_X(im,lm)
1069! & ACLL_CN_X(im,lm), ACIL_CN_X(im,lm), DQRL_X(im,lm), &
1070! & DZET(im,lm))
1071! & DZET(im,lm), qst3(im,lm))
1072! allocate (PFI_CN_X(im,0:lm), PFL_CN_X(im,0:lm))
1073
1074! do L=LM,1,-1
1075! do i=1,im
1076! DZET(i,L) = ZET(i,L) - ZET(i,L+1)
1077! tx1 = plo(i,l)*100.0
1078! est3 = min(tx1, fpvs(temp(i,l)))
1079! qst3(i,l) = min(eps*est3/max(tx1+epsm1*est3,1.0e-10),1.0)
1080!! MASS(i,l) = (ple(i,l) - ple(i,l-1)) * (100.0/grav)
1081! enddo
1082! enddo
1083! do k=1,lm
1084! do i=1,im
1085! REV_CN_X(i,k) = 0.0
1086! RSU_CN_X(i,k) = 0.0
1087! enddo
1088! enddo
1089! do k=0,lm
1090! do i=1,im
1091! PFI_CN_X(i,k) = 0.0
1092! PFL_CN_X(i,k) = 0.0
1093! enddo
1094! enddo
1095
1096! call macro_cloud (IM, LM, DT_MOIST, PLO, PLE, PK, FRLAND, &
1097! call macro_cloud (IM, LM, DT_MOIST, PLO, PLE, FRLAND, &
1099 call macro_cloud (im, lm, dt_moist, alf_fac, plo, ple, &
1100 & cnv_dqldt, &
1101! & CNV_MFD, CNV_DQLDT, &
1102! & CNV_MFD, CNV_DQLDT, CNV_PRC3, CNV_UPDF, &
1103! & U1, V1, temp, Q1, QLLS, QLCN, QILS, QICN, &
1104 & temp, q1, qlls, qlcn, qils, qicn, &
1105! & U1, V1, TH1, Q1, QLLS, QLCN, QILS, QICN, &
1106 & clcn, clls, &
1107! & CLCN, CLLS, CN_PRC2, CN_ARFX, CN_SNR, &
1108 & cloudparams, sclmfdfr, &
1109! & CLOUDPARAMS, SCLMFDFR, QST3, DZET, QDDF3, &
1110! & RHX_X, REV_CN_X, RSU_CN_X, &
1111! & ACLL_CN_X, ACIL_CN_X, PFL_CN_X, &
1112! & PFI_CN_X, DLPDF_X, DIPDF_X, &
1113! & ALPHT_X, CFPDF_X, DQRL_X, VFALLSN_CN_X, &
1114! & VFALLRN_CN_X, CNV_FICE, CNV_NDROP, CNV_NICE, &
1115 & alpht_x, &
1116 & cnv_fice, cnv_ndrop, cnv_nice, &
1117 & sc_ice, ncpl, ncpi, pfrz, &
1118 & lprnt, ipr, rhc, pdfflag, qc_min)
1119! & QRAIN_CN, QSNOW_CN, KCBL, lprnt, ipr, rhc)
1120
1121
1122! if (lprnt) write(0,*) ' in micro qicn3=',qicn(ipr,25)
1123! if(lprnt) write(0,*)' aft macro_cloud clcn=',clcn(ipr,:)
1124! if(lprnt) write(0,*)' aft macro_cloud clls=',clls(ipr,:)
1125! if(lprnt) write(0,*)' aft macro_cloud q1=',q1(ipr,:)
1126! if(lprnt) write(0,*)' aft macro_cloud qils=',qils(ipr,:)
1127
1128 do k=1,lm
1129 do i=1,im
1130 if (cnv_mfd(i,k) > 1.0e-6_kp) then
1131 tx1 = one / cnv_mfd(i,k)
1132 cnv_ndrop(i,k) = cnv_ndrop(i,k) * tx1
1133 cnv_nice(i,k) = cnv_nice(i,k) * tx1
1134 else
1135 cnv_ndrop(i,k) = zero
1136 cnv_nice(i,k) = zero
1137 endif
1138! temp(i,k) = th1(i,k) * PK(i,k)
1139 rad_cf(i,k) = min(clls(i,k)+clcn(i,k), one)
1140!
1141 if (iccn /= 1) then
1142 if (pfrz(i,k) > zero) then
1143 inc_nuc(i,k) = inc_nuc(i,k) * pfrz(i,k)
1144 nhet_nuc(i,k) = nhet_nuc(i,k) * pfrz(i,k)
1145 else
1146 inc_nuc(i,k) = zero
1147 nhet_nuc(i,k) = zero
1148 endif
1149 endif
1150
1151 enddo
1152 enddo
1153
1154
1155!make sure QI , NI stay within T limits
1156! call meltfrz_inst(IM, LM, TEMP, QLLS, QLCN, QILS, QICN, NCPL, NCPI)
1157
1158!============ a little treatment of cloud before micorphysics
1159! call update_cld(im,lm,DT_MOIST, ALPHT_X, qc_min &
1160! &, pdfflag, PLO , Q1, QLLS &
1161! &, QLCN, QILS, QICN, TEMP &
1162! &, CLLS, CLCN, SC_ICE, NCPI &
1163! &, NCPL)
1164!! &, NCPL, INC_NUC)
1165!============ Put cloud fraction back in contact with the PDF (Barahona et al., GMD, 2014)============
1166
1167!make sure QI , NI stay within T limits
1169 call meltfrz_inst(im, lm, temp, qlls, qlcn, qils, qicn, ncpl, ncpi)
1170
1171
1172! deallocate(RHX_X, CFPDF_X, VFALLSN_CN_X, &
1173 deallocate( &
1174! & QSNOW_CN, VFALLRN_CN_X, QRAIN_CN, REV_CN_X, RSU_CN_X,&
1175! & QSNOW_CN, QRAIN_CN, &
1176 & pfrz)
1177! & DLPDF_X, DIPDF_X, PFRZ, ACLL_CN_X, ACIL_CN_X, DQRL_X,&
1178! & PFI_CN_X, PFL_CN_X)
1179! & PFI_CN_X, PFL_CN_X, DZET, qst3, qddf3)
1180
1181 else
1182! do i=1,im
1183! CN_PRC2(i) = 0.0
1184! CN_SNR(i) = 0.0
1185! enddo
1186
1187
1188 endif ! .not. skip_macro
1189
1190
1191!===========================End of Cloud Macrophysics ========================
1192! --------------------------
1193!
1194
1195
1196!TVQX1 = SUM( ( Q1 + QLCN + QICN )*DM, 3)
1197
1198 do k=1,lm
1199 do i=1,im
1200 qcntot = qlcn(i,k) + qicn(i,k)
1201 ql_tot(i,k) = qlcn(i,k) + qlls(i,k)
1202 qi_tot(i,k) = qicn(i,k) + qils(i,k)
1203! Anning if negative, borrow water and ice from vapor 11/23/2016
1204 if (ql_tot(i,k) < zero) then
1205 q1(i,k) = q1(i,k) + ql_tot(i,k)
1206 temp(i,k) = temp(i,k) - lvbcp*ql_tot(i,k)
1207 ql_tot(i,k) = zero
1208 endif
1209 if (qi_tot(i,k) < zero) then
1210 q1(i,k) = q1(i,k) + qi_tot(i,k)
1211 temp(i,k) = temp(i,k) - lsbcp*qi_tot(i,k)
1212 qi_tot(i,k) = zero
1213 endif
1214 qtot = ql_tot(i,k) + qi_tot(i,k)
1215 if (qtot > zero) then
1216 fqa(i,k) = min(max(qcntot / qtot, zero), one)
1217 else
1218 fqa(i,k) = zero
1219 endif
1220 enddo
1221 enddo
1222
1223!=============================================================================================
1224!===========================Two-moment stratiform microphysics ===============================
1225!===========This is the implementation of the Morrison and Gettelman (2008) microphysics =====
1226!=============================================================================================
1229
1230 do i=1,im
1231 ls_snr(i) = zero
1232 ls_prc2(i) = zero
1233
1234 nbincontactdust = 1
1235
1236 do l=1,10
1237 do k=1,lm
1238 naconr8(k,l) = zero
1239 rndstr8(k,l) = 2.0e-7_kp
1240 enddo
1241 enddo
1242 do k=1,lm
1243 npccninr8(k) = zero
1244 naair8(k) = zero
1245 omegr8(k) = zero
1246
1247! tx1 = MIN(CLLS(I,k) + CLCN(I,k), 0.99)
1248 tx1 = min(clls(i,k) + clcn(i,k), one)
1249 if (tx1 > zero) then
1250 cldfr8(k) = min(max(tx1, 0.00001_kp), one)
1251 else
1252 cldfr8(k) = zero
1253 endif
1254
1255 if (temp(i,k) > tice) then
1256 liqcldfr8(k) = cldfr8(k)
1257 icecldfr8(k) = zero
1258 elseif (temp(i,k) <= t_ice_all) then
1259 liqcldfr8(k) = zero
1260 icecldfr8(k) = cldfr8(k)
1261 else
1262 icecldfr8(k) = cldfr8(k) * (tice - temp(i,k))/(tice-t_ice_all)
1263 liqcldfr8(k) = cldfr8(k) - icecldfr8(k)
1264 endif
1265
1266
1267 cldor8(k) = cldfr8(k)
1268 ter8(k) = temp(i,k)
1269 qvr8(k) = q1(i,k)
1270
1271 qcr8(k) = ql_tot(i,k)
1272 qir8(k) = qi_tot(i,k)
1273 ncr8(k) = max(ncpl(i,k), zero)
1274 nir8(k) = max(ncpi(i,k), zero)
1275 qrr8(k) = rnw(i,k)
1276 qsr8(k) = snw(i,k)
1277 qgr8(k) = qgl(i,k)
1278 nrr8(k) = max(ncpr(i,k), zero)
1279 nsr8(k) = max(ncps(i,k), zero)
1280 ngr8(k) = max(ncgl(i,k), zero)
1281
1282
1283 naair8(k) = inc_nuc(i,k)
1284 npccninr8(k) = cdnc_nuc(i,k)
1285
1286 if (cldfr8(k) >= 0.001_kp) then
1287 nimmr8(k) = min(dnhet_imm(i,k),ncr8(k)/(cldfr8(k)*dt_moist))
1288 else
1289 nimmr8(k) = zero
1290 endif
1291
1292
1293! if(aero_in) then
1294 aeroaux = aeroprops(i, k)
1295! else
1296! call init_Aer(AeroAux)
1297! end if
1299 call getinsubset(1, aeroaux, aeroaux_b)
1300 naux = aeroaux_b%nmods
1301 if (nbincontactdust < naux) then
1302 nbincontactdust = naux
1303 endif
1304 naconr8(k, 1:naux) = aeroaux_b%num(1:naux)
1305 rndstr8(k, 1:naux) = aeroaux_b%dpg(1:naux) * half
1306
1307! The following moved inside of if(fprcp <= 0) then loop
1308! Get black carbon properties for contact ice nucleation
1309! call getINsubset(2, AeroAux, AeroAux_b)
1310! nsootr8 (K) = sum(AeroAux_b%num)
1311! naux = AeroAux_b%nmods
1312! rnsootr8 (K) = sum(AeroAux_b%dpg(1:naux))/naux
1313
1314 pdelr8(k) = (ple(i,k) - ple(i,k-1)) * 100.0_kp
1315 rpdelr8(k) = one / pdelr8(k)
1316 plevr8(k) = 100.0_kp * plo(i,k)
1317 zmr8(k) = zlo(i,k)
1318 ficer8(k) = qir8(k) / (qcr8(k)+qir8(k) + 1.0e-10_kp)
1319 omegr8(k) = wsub(i,k)
1320! alphar8(k) = max(alpht_x(i,k)/maxval(alpht_x(i,:))*8.,0.5)
1321! alphar8(k) = qcvar2
1322 rhr8(k) = rhc(i,k)
1323
1324 END DO
1325 do k=1,lm+1
1326 pintr8(k) = ple(i,k-1) * 100.0_kp
1327 kkvhr8(k) = kh(i,k-1)
1328 END DO
1329
1330 lev_sed_strt = 0
1331 tx1 = one / pintr8(lm+1)
1332 do k=1,lm
1333 if (plevr8(k)*tx1 < sig_sed_strt) then
1334 lev_sed_strt(1) = k
1335 endif
1336 enddo
1337 lev_sed_strt(1) = max(lm/6, min(lm/3,lev_sed_strt(1)))
1338! if (kdt == 1) &
1339! write(0,*)' lev_sed_strt=',lev_sed_strt,' plevr8=',plevr8(lev_sed_strt), &
1340! ' pintr8=',pintr8(lm+1),' sig_sed_strt=',sig_sed_strt
1341!
1342! do k=1,lm
1343! if (cldfr8(k) <= 0.2 ) then
1344! alphar8(k) = 0.5
1345! elseif (cldfr8(k) <= 0.999) then
1346!! tx1 = 0.0284 * exp(4.4*cldfr8(k))
1347!! alphar8(k) = tx1 / (cldfr8(k) - tx1*(one-cldfr8(k)))
1348!! alphar8(k) = 0.5 + (7.5/0.799)*(cldfr8(k)-0.2)
1349! alphar8(k) = 0.5 + (7.5/0.799)*(cldfr8(k)-0.2)
1350! else
1351! alphar8(k) = 8.0
1352! endif
1353! alphar8(k) = min(8.0, max(alphar8(k), 0.5))
1354! enddo
1355
1356 kbmin = kcbl(i)
1357
1358!!!Call to MG microphysics. Lives in cldwat2m_micro.f
1359! ttendr8, qtendr8,cwtendr8, not used so far Anning noted August 2015
1360
1361 if (fprcp <= 0) then ! if fprcp = -1, then Anning's code for MG2 will be used
1362 ! if fprcp = 0, then MG1 is used
1363
1364! Get black carbon properties for contact ice nucleation
1365 do k=1,lm
1366 call getinsubset(2, aeroaux, aeroaux_b)
1367 nsootr8(k) = sum(aeroaux_b%num)
1368 naux = aeroaux_b%nmods
1369 rnsootr8(k) = sum(aeroaux_b%dpg(1:naux))/naux
1370 enddo
1371
1372 call mmicro_pcond ( ncolmicro, ncolmicro, &
1373 & dt_kp, ter8, ttendr8, &
1374 & ncolmicro, lm , qvr8, &
1375 & qtendr8, cwtendr8, qcr8, qir8, ncr8, nir8, &
1376 & abs(fprcp), qrr8, qsr8, nrr8, nsr8, &
1377 & plevr8, pdelr8, cldfr8, liqcldfr8, &
1378 & icecldfr8, cldor8, pintr8, &
1379 & rpdelr8, zmr8, rate1ord_cw2pr, &
1380 & naair8, npccninr8, &
1381 & rndstr8, naconr8, rhdfdar8, rhu00r8, ficer8, &
1382 & tlatr8, qvlatr8, qctendr8, &
1383 & qitendr8, nctendr8, nitendr8, effcr8, &
1384 & effc_fnr8, effir8, prectr8, precir8, &
1385 & nevaprr8, evapsnowr8, &
1386 & prainr8, prodsnowr8, cmeoutr8, &
1387 & deffir8, pgamradr8, lamcradr8, &
1388 & qsoutr8, dsoutr8, qroutr8, droutr8, &
1389 & qcsevapr8, qisevapr8, qvresr8, &
1390 & cmeioutr8, vtrmcr8, vtrmir8, &
1391 & qcsedtenr8, qisedtenr8, praor8, prcor8, &
1392 & mnucccor8, mnucctor8, msacwior8, &
1393 & psacwsor8, bergsor8, bergor8, &
1394 & meltor8, homoor8, qcresor8, prcior8, &
1395 & praior8, qiresor8, mnuccror8, &
1396 & pracsor8, meltsdtr8, frzrdtr8, ncalr8, &
1397 & ncair8, mnuccdor8, &
1398 & nnucctor8, nsoutr8, nroutr8, &
1399 & ncnstr8, ninstr8, nimmr8, disp_liu, &
1400 & nsootr8, rnsootr8, ui_scale, dcrit, &
1401 & nnuccdor8, nnucccor8, &
1402 & nsacwior8, nsubior8, nprcior8, &
1403 & npraior8, npccnor8, npsacwsor8, &
1404 & nsubcor8, npraor8, nprc1or8, &
1405 & tlatauxr8, nbincontactdust, &
1406 & lprnt, xlat(i), xlon(i), rhr8)
1407
1408! if (lprint) write(0,*)' prectr8=',prectr8(1), &
1409! & ' precir8=',precir8(1)
1410
1411 ls_prc2(i) = max(1000.0_kp*(prectr8(1)-precir8(1)), zero)
1412 ls_snr(i) = max(1000.0_kp*precir8(1), zero)
1413
1414
1415 do k=1,lm
1416 ql_tot(i,k) = ql_tot(i,k) + qctendr8(k)*dt_kp
1417 qi_tot(i,k) = qi_tot(i,k) + qitendr8(k)*dt_kp
1418 q1(i,k) = q1(i,k) + qvlatr8(k)*dt_kp
1419! if(lprnt .and. i == ipr) write(0,*)' k=',k,' q1aftm=',q1(i,k) &
1420! &,' qvlatr8=',qvlatr8(k)
1421 temp(i,k) = temp(i,k) + tlatr8(k)*dt_kp*onebcp
1422
1423 ncpl(i,k) = max(ncpl(i,k) + nctendr8(k) * dt_kp, zero)
1424 ncpi(i,k) = max(ncpi(i,k) + nitendr8(k) * dt_kp, zero)
1425 rnw(i,k) = qrr8(k)
1426 snw(i,k) = qsr8(k)
1427 ncpr(i,k) = nrr8(k)
1428 ncps(i,k) = nsr8(k)
1429
1430 cldreffl(i,k) = min(max(effcr8(k), 10.0_kp), 150.0_kp)
1431 cldreffi(i,k) = min(max(effir8(k), 20.0_kp), 150.0_kp)
1432 cldreffr(i,k) = max(droutr8(k)*0.5_kp*1.0e6_kp, 150.0_kp)
1433 cldreffs(i,k) = max(0.192_kp*dsoutr8(k)*0.5_kp*1.0e6_kp, 250.0_kp)
1434
1435 enddo ! K loop
1436
1437 elseif (fprcp == 1) then ! Call MG2
1438! --------
1439! if (lprnt .and. i == ipr) then
1440! write(0,*)' bef micro_mg_tend ter8= ', ter8(:)
1441! write(0,*)' bef micro_mg_tend qvr8= ', qvr8(:),'dt_kp=',dt_kp
1442! write(0,*)' bef micro_mg_tend rhr8= ', rhr8(:)
1443! endif
1444 lprint = lprnt .and. i == ipr
1445 ltrue = any(qcr8 >= qsmall) .or. any(qir8 >= qsmall) &
1446 .or. any(qsr8 >= qsmall) .or. any(qrr8 >= qsmall)
1447 if (ltrue) then
1448 alphar8(:) = qcvar2
1449
1450! if(lprint) then
1451! write(0,*)' calling micro_mg_tend2_0 qcvar2=',qcvar2
1452! write(0,*)' qcr8=',qcr8(:)
1453! write(0,*)' ncr8=',ncr8(:)
1454! write(0,*)' npccninr8=',npccninr8(:)
1455! write(0,*)' plevr8=',plevr8(:)
1456! write(0,*)' ter8=',ter8(:)
1457! endif
1458
1459 call micro_mg_tend2_0 ( &
1460 & ncolmicro, lm, dt_kp, &
1461 & ter8, qvr8, &
1462 & qcr8, qir8, &
1463 & ncr8, nir8, &
1464 & qrr8, qsr8, &
1465 & nrr8, nsr8, &
1466 & alphar8, 1., &
1467 & plevr8, pdelr8, &
1468! & cldfr8, liqcldfr8, icecldfr8, rhc, &
1469 & cldfr8, liqcldfr8, icecldfr8, rhr8, &
1470 & qcsinksum_rate1ord, &
1471 & naair8, npccninr8, &
1472 & rndstr8, naconr8, &
1473 & tlatr8, qvlatr8, &
1474 & qctendr8, qitendr8, &
1475 & nctendr8, nitendr8, &
1476 & qrtend, qstend, &
1477 & nrtend, nstend, &
1478 & effcr8, effc_fnr8, effir8, &
1479 & sadice, sadsnow, &
1480 & prectr8, precir8, &
1481 & nevaprr8, evapsnowr8, &
1482 & am_evp_st, &
1483 & prainr8, prodsnowr8, &
1484 & cmeoutr8, deffir8, &
1485 & pgamradr8, lamcradr8, &
1486 & qsoutr8, dsoutr8, &
1487 & lflx, iflx, &
1488 & rflx, sflx, qroutr8, &
1489 & reff_rain, reff_snow, &
1490 & qcsevapr8, qisevapr8, qvresr8, &
1491 & cmeioutr8, vtrmcr8, vtrmir8, &
1492 & umr, ums, &
1493 & qcsedtenr8, qisedtenr8, &
1494 & qrsedten, qssedten, &
1495 & praor8, prcor8, &
1496 & mnucccor8, mnucctor8, msacwior8, &
1497 & psacwsor8, bergsor8, bergor8, &
1498 & meltor8, homoor8, &
1499 & qcresor8, prcior8, praior8, &
1500 & qiresor8, mnuccror8, pracsor8, &
1501 & meltsdtr8, frzrdtr8, mnuccdor8, &
1502 & nroutr8, nsoutr8, &
1503 & refl, arefl, areflz, &
1504 & frefl, csrfl, acsrfl, &
1505 & fcsrfl, rercld, &
1506 & ncair8, ncalr8, &
1507 & qrout2, qsout2, &
1508 & nrout2, nsout2, &
1509 & drout2, dsout2, &
1510 & freqs, freqr, &
1511 & nfice, qcrat, &
1512 & prer_evap, xlat(i), xlon(i), lprint, iccn, &
1513 & lev_sed_strt)
1514!
1515 ls_prc2(i) = max(1000.0_kp*(prectr8(1)-precir8(1)), zero)
1516 ls_snr(i) = max(1000.0_kp*precir8(1), zero)
1517 do k=1,lm
1518 ql_tot(i,k) = ql_tot(i,k) + qctendr8(k)*dt_kp
1519 qi_tot(i,k) = qi_tot(i,k) + qitendr8(k)*dt_kp
1520 q1(i,k) = q1(i,k) + qvlatr8(k)*dt_kp
1521 temp(i,k) = temp(i,k) + tlatr8(k)*dt_kp*onebcp
1522 rnw(i,k) = rnw(i,k) + qrtend(k)*dt_kp
1523 snw(i,k) = snw(i,k) + qstend(k)*dt_kp
1524
1525 ncpl(i,k) = max(ncpl(i,k) + nctendr8(k)*dt_kp, zero)
1526 ncpi(i,k) = max(ncpi(i,k) + nitendr8(k)*dt_kp, zero)
1527 ncpr(i,k) = max(ncpr(i,k) + nrtend(k)*dt_kp, zero)
1528 ncps(i,k) = max(ncps(i,k) + nstend(k)*dt_kp, zero)
1529
1530 cldreffl(i,k) = min(max(effcr8(k), 10.0_kp), 150.0_kp)
1531 cldreffi(i,k) = min(max(effir8(k), 20.0_kp), 150.0_kp)
1532 cldreffr(i,k) = max(reff_rain(k), 150.0_kp)
1533 cldreffs(i,k) = max(reff_snow(k), 250.0_kp)
1534 enddo ! K loop
1535! if (lprint) then
1536! write(0,*)' aft micro_mg_tend temp= ', temp(i,:)
1537! write(0,*)' aft micro_mg_tend q1= ', q1(i,:)
1538! write(0,*)' aft micro_mg_tend LS_PRC2= ', LS_PRC2(i),' ls_snr=',ls_snr(i)
1539! endif
1540 else
1541 ls_prc2(i) = zero
1542 ls_snr(i) = zero
1543 do k=1,lm
1544 cldreffl(i,k) = 10.0_kp
1545 cldreffi(i,k) = 50.0_kp
1546 cldreffr(i,k) = 1000.0_kp
1547 cldreffs(i,k) = 250.0_kp
1548 enddo ! K loop
1549 endif
1550!
1551 else ! Call MG3
1552! --------
1553 ltrue = any(qcr8 >= qsmall) .or. any(qir8 >= qsmall) &
1554 .or. any(qsr8 >= qsmall) .or. any(qrr8 >= qsmall) &
1555 .or. any(qgr8 >= qsmall)
1556 lprint = lprnt .and. i == ipr
1557 if (ltrue) then
1558 alphar8(:) = qcvar3
1559
1560! if(lprint) then
1561! write(0,*)' calling micro_mg_tend3_0 qcvar3=',qcvar3,' i=',i
1562! write(0,*)' qcr8=',qcr8(:)
1563! write(0,*)' qir8=',qir8(:)
1564! write(0,*)' ncr8=',ncr8(:)
1565! write(0,*)' nir8=',nir8(:)
1566! write(0,*)' npccninr8=',npccninr8(:)
1567! write(0,*)' plevr8=',plevr8(:)
1568! write(0,*)' ter8=',ter8(:)
1569! endif
1572 call micro_mg_tend3_0 ( &
1573 & ncolmicro, lm, dt_kp, &
1574 & ter8, qvr8, &
1575 & qcr8, qir8, &
1576 & ncr8, nir8, &
1577 & qrr8, qsr8, &
1578 & nrr8, nsr8, &
1579 & qgr8, ngr8, &
1580 & alphar8, 1., &
1581 & plevr8, pdelr8, &
1582! & cldfr8, liqcldfr8, icecldfr8, rhc, &
1583 & cldfr8, liqcldfr8, icecldfr8, rhr8, &
1584 & qcsinksum_rate1ord, &
1585 & naair8, npccninr8, &
1586 & rndstr8, naconr8, &
1587 & tlatr8, qvlatr8, &
1588 & qctendr8, qitendr8, &
1589 & nctendr8, nitendr8, &
1590 & qrtend, qstend, &
1591 & nrtend, nstend, &
1592!
1593 & qgtend, ngtend, &
1594!
1595 & effcr8, effc_fnr8, effir8, &
1596 & sadice, sadsnow, &
1597 & prectr8, precir8, &
1598 & nevaprr8, evapsnowr8, &
1599 & am_evp_st, &
1600 & prainr8, prodsnowr8, &
1601 & cmeoutr8, deffir8, &
1602 & pgamradr8, lamcradr8, &
1603 & qsoutr8, dsoutr8, &
1604!
1605 & qgoutr8, ngoutr8, dgoutr8, &
1606!
1607 & lflx, iflx, gflx, &
1608!
1609 & rflx, sflx, qroutr8, &
1610!
1611 & reff_rain, reff_snow, reff_grau, &
1612!
1613 & qcsevapr8, qisevapr8, qvresr8, &
1614 & cmeioutr8, vtrmcr8, vtrmir8, &
1615 & umr, ums, &
1616!
1617 & umg, qgsedtenr8, &
1618!
1619 & qcsedtenr8, qisedtenr8, &
1620 & qrsedten, qssedten, &
1621 & praor8, prcor8, &
1622 & mnucccor8, mnucctor8, msacwior8, &
1623 & psacwsor8, bergsor8, bergor8, &
1624 & meltor8, homoor8, &
1625 & qcresor8, prcior8, praior8, &
1626!
1627 & qiresor8, mnuccror8, mnuccrior8, pracsor8, &
1628!
1629 & meltsdtr8, frzrdtr8, mnuccdor8, &
1630!
1631 & pracgr8, psacwgr8, pgsacwr8, &
1632 & pgracsr8, prdgr8, &
1633 & qmultgr8, qmultrgr8, psacrr8, &
1634 & npracgr8, nscngr8, ngracsr8, &
1635 & nmultgr8, nmultrgr8, npsacwgr8, &
1636!
1637 & nroutr8, nsoutr8, &
1638 & refl, arefl, areflz, &
1639 & frefl, csrfl, acsrfl, &
1640 & fcsrfl, rercld, &
1641 & ncair8, ncalr8, &
1642 & qrout2, qsout2, &
1643 & nrout2, nsout2, &
1644 & drout2, dsout2, &
1645!
1646 & qgout2, ngout2, dgout2, freqg, &
1647 & freqs, freqr, &
1648 & nfice, qcrat, &
1649 & prer_evap, xlat(i), xlon(i), lprint, iccn, &
1650 & lev_sed_strt)
1651
1652 ls_prc2(i) = max(1000.0_kp*(prectr8(1)-precir8(1)), zero)
1653 ls_snr(i) = max(1000.0_kp*precir8(1), zero)
1654 do k=1,lm
1655 ql_tot(i,k) = ql_tot(i,k) + qctendr8(k)*dt_kp
1656 qi_tot(i,k) = qi_tot(i,k) + qitendr8(k)*dt_kp
1657 q1(i,k) = q1(i,k) + qvlatr8(k)*dt_kp
1658 temp(i,k) = temp(i,k) + tlatr8(k)*dt_kp*onebcp
1659 rnw(i,k) = rnw(i,k) + qrtend(k)*dt_kp
1660 snw(i,k) = snw(i,k) + qstend(k)*dt_kp
1661 qgl(i,k) = qgl(i,k) + qgtend(k)*dt_kp
1662
1663 ncpl(i,k) = max(ncpl(i,k) + nctendr8(k)*dt_kp, zero)
1664 ncpi(i,k) = max(ncpi(i,k) + nitendr8(k)*dt_kp, zero)
1665 ncpr(i,k) = max(ncpr(i,k) + nrtend(k)*dt_kp, zero)
1666 ncps(i,k) = max(ncps(i,k) + nstend(k)*dt_kp, zero)
1667 ncgl(i,k) = max(ncgl(i,k) + ngtend(k)*dt_kp, zero)
1668
1669 cldreffl(i,k) = min(max(effcr8(k), 10.0_kp), 150.0_kp)
1670 cldreffi(i,k) = min(max(effir8(k), 20.0_kp), 150.0_kp)
1671 cldreffr(i,k) = max(reff_rain(k), 150.0_kp)
1672 cldreffs(i,k) = max(reff_snow(k), 250.0_kp)
1673 cldreffg(i,k) = max(reff_grau(k), 250.0_kp)
1674 enddo ! K loop
1675! if (lprint) then
1676! write(0,*)' aft micro_mg_tend temp= ', temp(i,:)
1677! write(0,*)' aft micro_mg_tend q1= ', q1(i,:)
1678! write(0,*)' aft micro_mg_tend LS_PRC2= ', LS_PRC2(i),' ls_snr=',ls_snr(i)
1679! endif
1680 else
1681 ls_prc2(i) = zero
1682 ls_snr(i) = zero
1683 do k=1,lm
1684 cldreffl(i,k) = 10.0_kp
1685 cldreffi(i,k) = 50.0_kp
1686 cldreffr(i,k) = 1000.0_kp
1687 cldreffs(i,k) = 250.0_kp
1688 cldreffg(i,k) = 250.0_kp
1689 enddo ! K loop
1690 endif
1691 endif
1692
1693 enddo ! I loop
1694!============================================Finish 2-moment micro implementation===========================
1695
1696!TVQX1 = SUM( ( Q1 + QL_TOT + QI_TOT(1:im,:,:))*DM, 3) &
1697
1698 if (skip_macro) then
1699 do k=1,lm
1700 do i=1,im
1701 qlcn(i,k) = ql_tot(i,k) * fqa(i,k)
1702 qlls(i,k) = ql_tot(i,k) - qlcn(i,k)
1703 qicn(i,k) = qi_tot(i,k) * fqa(i,k)
1704 qils(i,k) = qi_tot(i,k) - qicn(i,k)
1705
1706 CALL fix_up_clouds_2m(q1(i,k), temp(i,k), qlls(i,k), &
1707 & qils(i,k), clls(i,k), qlcn(i,k), &
1708 & qicn(i,k), clcn(i,k), ncpl(i,k), &
1709 & ncpi(i,k), qc_min)
1710
1711 ql_tot(i,k) = qlls(i,k) + qlcn(i,k)
1712 qi_tot(i,k) = qils(i,k) + qicn(i,k)
1713 if (rnw(i,k) <= qc_min(1)) then
1714 ncpr(i,k) = zero
1715 rnw(i,k) = zero
1716 elseif (ncpr(i,k) <= nmin) then ! make sure NL > 0 if Q >0
1717 ncpr(i,k) = max(rnw(i,k) / (fourb3 * pi *rl_cub*997.0_kp), nmin)
1718 endif
1719 if (snw(i,k) <= qc_min(2)) then
1720 ncps(i,k) = zero
1721 snw(i,k) = zero
1722 elseif (ncps(i,k) <= nmin) then
1723 ncps(i,k) = max(snw(i,k) / (fourb3 * pi *rl_cub*500.0_kp), nmin)
1724 endif
1725 if (qgl(i,k) <= qc_min(2)) then
1726 ncgl(i,k) = zero
1727 qgl(i,k) = zero
1728 elseif (ncgl(i,k) <= nmin) then
1729 ncgl(i,k) = max(qgl(i,k) / (fourb3 * pi *rl_cub*500.0_kp), nmin)
1730 endif
1731 enddo
1732 enddo
1733 else
1734 do k=1,lm
1735 do i=1,im
1736 qlcn(i,k) = ql_tot(i,k) * fqa(i,k)
1737 qlls(i,k) = ql_tot(i,k) - qlcn(i,k)
1738 qicn(i,k) = qi_tot(i,k) * fqa(i,k)
1739 qils(i,k) = qi_tot(i,k) - qicn(i,k)
1740 enddo
1741 enddo
1742
1744 call update_cld(im, lm, dt_moist, alpht_x, qc_min &
1745 &, pdfflag, plo, q1, qlls, qlcn &
1746 &, qils, qicn, temp, clls, clcn &
1747 &, sc_ice, ncpi, ncpl)
1748
1749! if(lprnt) write(0,*)' aft update_cloud clls=',clls(ipr,:)
1750
1751 do k=1,lm
1752 do i=1,im
1753 ql_tot(i,k) = qlls(i,k) + qlcn(i,k)
1754 qi_tot(i,k) = qils(i,k) + qicn(i,k)
1755!
1756 if (rnw(i,k) <= qc_min(1)) then
1757 ncpr(i,k) = zero
1758 rnw(i,k) = zero
1759 elseif (ncpr(i,k) <= nmin) then ! make sure NL > 0 if Q >0
1760 ncpr(i,k) = max(rnw(i,k) / (fourb3 * pi *rl_cub*997.0_kp), nmin)
1761 endif
1762 if (snw(i,k) <= qc_min(2)) then
1763 ncps(i,k) = zero
1764 snw(i,k) = zero
1765 elseif (ncps(i,k) <= nmin) then
1766 ncps(i,k) = max(snw(i,k) / (fourb3 * pi *rl_cub*500.0_kp), nmin)
1767 endif
1768 if (qgl(i,k) <= qc_min(2)) then
1769 ncgl(i,k) = zero
1770 qgl(i,k) = zero
1771 elseif (ncgl(i,k) <= nmin) then
1772 ncgl(i,k) = max(qgl(i,k) / (fourb3 * pi *rl_cub*500.0_kp), nmin)
1773 endif
1774 enddo
1775 enddo
1776 deallocate(cnv_mfd,cnv_fice,cnv_ndrop,cnv_nice)
1777! deallocate(CNV_MFD,CNV_PRC3,CNV_FICE,CNV_NDROP,CNV_NICE)
1778 endif
1779
1780! do I=1,IM
1781! TPREC(i) = CN_PRC2(i) + LS_PRC2(i) + CN_SNR(i) + LS_SNR(i)
1782! enddo
1783
1784 do k= 1, lm
1785 do i=1,im
1786 if (qi_tot(i,k) <= zero) ncpi(i,k) = zero
1787 if (ql_tot(i,k) <= zero) ncpl(i,k) = zero
1788 end do
1789 end do
1790
1791
1792!=============================================End Stratiform cloud processes==========================================
1793!======================================================================================================================
1794!===========================Clean stuff and send it to radiation ======================================================
1795!======================================================================================================================
1796! outputs
1797 if(flipv) then
1798 DO k=1, lm
1799 ll = lm-k+1
1800 DO i = 1,im
1801 t_io(i,k) = temp(i,ll)
1802 q_io(i,k) = q1(i,ll)
1803 ncpi_io(i,k) = ncpi(i,ll)
1804 ncpl_io(i,k) = ncpl(i,ll)
1805 rnw_io(i,k) = rnw(i,ll)
1806 snw_io(i,k) = snw(i,ll)
1807 qgl_io(i,k) = qgl(i,ll)
1808 ncpr_io(i,k) = ncpr(i,ll)
1809 ncps_io(i,k) = ncps(i,ll)
1810 ncgl_io(i,k) = ncgl(i,ll)
1811 lwm_o(i,k) = ql_tot(i,ll)
1812 qi_o(i,k) = qi_tot(i,ll)
1813 END DO
1814 END DO
1815 if (skip_macro) then
1816 DO k=1, lm
1817 ll = lm-k+1
1818 DO i = 1,im
1819 clls_io(i,k) = max(zero, min(clls(i,ll)+clcn(i,ll),one))
1820 enddo
1821 enddo
1822 else
1823 DO k=1, lm
1824 ll = lm-k+1
1825 DO i = 1,im
1826 clls_io(i,k) = clls(i,ll)
1827 enddo
1828 enddo
1829 endif
1830 else
1831 DO k=1, lm
1832 DO i = 1,im
1833 t_io(i,k) = temp(i,k)
1834 q_io(i,k) = q1(i,k)
1835 ncpi_io(i,k) = ncpi(i,k)
1836 ncpl_io(i,k) = ncpl(i,k)
1837 rnw_io(i,k) = rnw(i,k)
1838 snw_io(i,k) = snw(i,k)
1839 qgl_io(i,k) = qgl(i,k)
1840 ncpr_io(i,k) = ncpr(i,k)
1841 ncps_io(i,k) = ncps(i,k)
1842 ncgl_io(i,k) = ncgl(i,k)
1843 lwm_o(i,k) = ql_tot(i,k)
1844 qi_o(i,k) = qi_tot(i,k)
1845 END DO
1846 END DO
1847 if (skip_macro) then
1848 DO k=1, lm
1849 DO i = 1,im
1850 clls_io(i,k) = max(zero, min(clls(i,k)+clcn(i,k),one))
1851 enddo
1852 enddo
1853 else
1854 DO k=1, lm
1855 DO i = 1,im
1856 clls_io(i,k) = clls(i,k)
1857 enddo
1858 enddo
1859 endif
1860 endif ! end of flipv if
1861
1862 DO i = 1,im
1863 tx1 = ls_prc2(i) + ls_snr(i)
1864 rn_o(i) = tx1 * dt_i * 0.001_kp
1865
1866 if (rn_o(i) < rainmin) then
1867 sr_o(i) = zero
1868 else
1869 sr_o(i) = max(zero, min(one, ls_snr(i)/tx1))
1870 endif
1871 ENDDO
1872
1873 if (allocated(alpht_x)) deallocate (alpht_x)
1874
1875! if (lprnt) then
1876! write(0,*)' rn_o=',rn_o(ipr),' ls_prc2=',ls_prc2(ipr),' ls_snr=',ls_snr(ipr),' kdt=',kdt
1877! write(0,*)' end micro_mg_tend t_io= ', t_io(ipr,:)
1878! write(0,*)' end micro_mg_tend clls_io= ', clls_io(ipr,:)
1879! endif
1880! do k=1,lm
1881! do i=1,im
1882! dum(i,k) = clls_io(i,k)
1883! enddo
1884! enddo
1885! do k=2,lm-1
1886! do i=1,im
1887! clls_io(i,k) = 0.25*dum(i,k-1) + 0.5*dum(i,k)+0.25*dum(i,k+1)
1888! enddo
1889! enddo
1890! do i=1,im
1891! clls_io(i,lm) = 0.5 * (dum(i,lm-1) + dum(i,lm))
1892! enddo
1893
1894
1895
1896!=======================================================================
1897
1898 end subroutine m_micro_run
1900
1901!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1902
1903!DONIF Calculate the Brunt_Vaisala frequency
1904
1905!===============================================================================
1911 subroutine gw_prof (pcols, pver, ncol, t, pm, pi, rhoi, ni, ti, &
1912 nm, sph)
1913 use machine , only : kind_phys
1914 use physcons, grav => con_g, cp => con_cp, rgas => con_rd, &
1915 fv => con_fvirt
1916 implicit none
1917 integer, parameter :: kp = kind_phys
1918!-----------------------------------------------------------------------
1919! Compute profiles of background state quantities for the multiple
1920! gravity wave drag parameterization.
1921!
1922! The parameterization is assumed to operate only where water vapor
1923! concentrations are negligible in determining the density.
1924!-----------------------------------------------------------------------
1925!------------------------------Arguments--------------------------------
1926 integer, intent(in) :: ncol, pcols, pver
1927
1928
1929
1930 real(kind=kind_phys), intent(in) :: t(pcols,pver)
1931 real(kind=kind_phys), intent(in) :: pm(pcols,pver)
1932 real(kind=kind_phys), intent(in) :: pi(pcols,0:pver)
1933 real(kind=kind_phys), intent(in) :: sph(pcols,pver)
1934
1935 real(kind=kind_phys), intent(out) :: rhoi(pcols,0:pver)
1936 real(kind=kind_phys), intent(out) :: ni(pcols,0:pver)
1937 real(kind=kind_phys), intent(out) :: ti(pcols,0:pver)
1938 real(kind=kind_phys), intent(out) :: nm(pcols,pver)
1939
1940 real(kind=kind_phys), parameter :: r=rgas, cpair=cp, g=grav, &
1941 oneocp=1.0_kp/cp, n2min=1.0e-8_kp
1942
1943!---------------------------Local storage-------------------------------
1944 integer :: ix,kx
1945
1946 real :: dtdp, n2
1947
1948!-----------------------------------------------------------------------------
1950!-----------------------------------------------------------------------------
1951
1952! The top interface values are calculated assuming an isothermal atmosphere
1953! above the top level.
1954 kx = 0
1955 do ix = 1, ncol
1956 ti(ix,kx) = t(ix,kx+1)
1957 rhoi(ix,kx) = pi(ix,kx) / (r*(ti(ix,kx)*(1.0_kp+fv*sph(ix,kx+1))))
1958 ni(ix,kx) = sqrt(g*g / (cpair*ti(ix,kx)))
1959 end do
1960
1961! Interior points use centered differences
1962 do kx = 1, pver-1
1963 do ix = 1, ncol
1964 ti(ix,kx) = 0.5_kp * (t(ix,kx) + t(ix,kx+1))
1965 rhoi(ix,kx) = pi(ix,kx) / (r*ti(ix,kx)*(1.0_kp+0.5_kp*fv*(sph(ix,kx)+sph(ix,kx+1))))
1966 dtdp = (t(ix,kx+1)-t(ix,kx)) / (pm(ix,kx+1)-pm(ix,kx))
1967 n2 = g*g/ti(ix,kx) * (oneocp - rhoi(ix,kx)*dtdp)
1968 ni(ix,kx) = sqrt(max(n2min, n2))
1969 end do
1970 end do
1971
1972! Bottom interface uses bottom level temperature, density; next interface
1973! B-V frequency.
1974 kx = pver
1975 do ix = 1, ncol
1976 ti(ix,kx) = t(ix,kx)
1977 rhoi(ix,kx) = pi(ix,kx) / (r*ti(ix,kx)*(1.0_kp+fv*sph(ix,kx)))
1978 ni(ix,kx) = ni(ix,kx-1)
1979 end do
1980
1981!-----------------------------------------------------------------------------
1983!-----------------------------------------------------------------------------
1984 do kx=1,pver
1985 do ix=1,ncol
1986 nm(ix,kx) = 0.5_kp * (ni(ix,kx-1) + ni(ix,kx))
1987 end do
1988 end do
1989
1990 return
1991 end subroutine gw_prof
1993
1996 subroutine find_cldtop(ncol, pver, cf, kcldtop)
1997 implicit none
1998
1999 integer, intent(in) :: pver , ncol
2000 real, intent(in) :: cf(ncol,pver)
2001 integer, intent(out) :: kcldtop
2002 integer :: kuppest, ibot, k
2003 real :: stab, cfcrit, cf00, cfp1
2004
2005
2006 ibot = pver-1
2007 kcldtop = ibot+1
2008 kuppest = 20
2009 cfcrit = 1.0d-2
2010
2011
2012 do k = kuppest , ibot
2013 cfp1 = cf(ncol, k+1)
2014
2015 if ( ( cfp1 >= cfcrit ) ) then
2016 kcldtop = k +1
2017 exit
2018 end if
2019 end do
2020
2021 if (kcldtop >= ibot) then
2022 kcldtop = pver
2023 return
2024 endif
2025
2026
2027 end subroutine find_cldtop
2028
2029end module m_micro
subroutine, public aer_cloud_init()
This subroutine calculates.
Definition aer_cloud.F:116
subroutine, public init_aer(aerout)
This subroutine initialize aerosol properties in MG sheme.
Definition aer_cloud.F:1182
subroutine, public aerconversion1(aer_mass, aerpr)
This subroutine sets the properties of the aerosol distributions.
Definition aer_cloud.F:930
subroutine, public getinsubset(typ, aerin, aerout)
This subroutine extracts aerosol props with INactive = typ.
Definition aer_cloud.F:1102
subroutine, public aerosol_activate(tparc_in, pparc_in, sigwparc_in, wparc_ls, aer_props, npre_in, dpre_in, ccn_diagr8, ndropr8, cdncr8, smaxliqr8, incr8, smaxicer8, nheticer8, inimmr8, dinimmr8, ncdepr8, ncdhfr8, sc_icer8, fdust_immr8, fdust_depr8, fdust_dhfr8, nlimr8, use_average_v, ccn_param, in_param, fd_dust, fd_soot, pfrz_inc_r8, sigma_nuc, rhi_cell, nccn)
This subroutine sets the variables needed for the activation subroutines and return the activated dro...
Definition aer_cloud.F:191
subroutine, public fix_up_clouds_2m(qv, te, qlc, qic, cf, qla, qia, af, nl, ni, qc_min)
Definition cldmacro.F:623
subroutine, public update_cld(irun, lm, dt, alpha, qc_min, pdfshape, pl, qv, qcl, qal, qci, qai, te, cf, af, scice, ni, nl)
This subroutine calculates the cloud fraction that would correspond to the current condensate.
Definition cldmacro.F:713
subroutine, public macro_cloud(irun, lm, dt, alf_fac, pp_dev, ppe_dev, qlwdtr_dev, th_dev, q_dev, qlw_ls_dev, qlw_an_dev, qiw_ls_dev, qiw_an_dev, anvfrc_dev, cldfrc_dev, physparams, sclmfdfr, alpht_dev, cnv_fice_dev, cnv_ndrop_dev, cnv_nice_dev, scice_dev, ncpl_dev, ncpi_dev, pfrz_dev, lprnt, ipr, rhc, pdfflag, qc_min)
This subroutine is the cloud macrophysics scheme in MG micriphysics.
Definition cldmacro.F:123
subroutine, public meltfrz_inst(im, lm, te, qcl, qal, qci, qai, nl, ni)
This subroutine calculates instantaneous freezing of condensate.
Definition cldmacro.F:2098
real(r8), parameter rh2o
real(r8), parameter tmelt
subroutine, public mmicro_pcond(lchnk, ncol, deltatin, tn, ttend, pcols, pver, qn, qtend, cwtend, qc, qi, nc, ni, fprcp, qrn, qsnw, nrn, nsnw, p, pdel, cldn, liqcldf, icecldf, cldo, pint, rpdel, zm, rate1ord_cw2pr_st, naai, npccnin, rndst, nacon, rhdfda, rhu00, fice, tlat, qvlat, qctend, qitend, nctend, nitend, effc, effc_fn, effi, prect, preci, nevapr, evapsnow, prain, prodsnow, cmeout, deffi, pgamrad, lamcrad, qsout2, dsout2, qrout2, drout2, qcsevap, qisevap, qvres, cmeiout, vtrmc, vtrmi, qcsedten, qisedten, prao, prco, mnuccco, mnuccto, msacwio, psacwso, bergso, bergo, melto, homoo, qcreso, prcio, praio, qireso, mnuccro, pracso, meltsdt, frzrdt, ncal, ncai, mnuccdo, nnuccto, nsout2, nrout2, ncnst, ninst, nimm, miu_disp, nsoot, rnsoot, ui_scale, dcrit, nnuccdo, nnuccco, nsacwio, nsubio, nprcio, npraio, npccno, npsacwso, nsubco, nprao, nprc1o, tlataux, nbincontactdust, lprint, xlat, xlon, rhc)
This subroutine is the microphysics routine for each timestep goes here...
real(r8), parameter oneb3
real(r8), parameter cpair
real(r8), parameter latice
real(r8), parameter gravit
real(r8), parameter half
real(r8), parameter one
real(r8), private qsmall
real(r8), parameter latvap
real(r8), parameter rair
subroutine, public ini_micro(dcs_, qcvar_, ts_auto_ice_)
This subroutine initializes constants for MG microphysics.
real(r8), private qcvar
real(r8), parameter onebcp
subroutine, public micro_mg_init(kind, gravit, rair, rh2o, cpair, eps, tmelt_in, latvap, latice, rhmini_in, micro_mg_dcs, ts_auto, mg_qcvar, microp_uniform_in, do_cldice_in, use_hetfrz_classnuc_in, micro_mg_precip_frac_method_in, micro_mg_berg_eff_factor_in, allow_sed_supersat_in, do_sb_physics_in, do_ice_gmao_in, do_liq_liu_in, nccons_in, nicons_in, ncnst_in, ninst_in)
This subroutine calculates.
subroutine find_cldtop(ncol, pver, cf, kcldtop)
This subroutine is to find cloud top based on cloud fraction.
Definition m_micro.F90:1997
subroutine gw_prof(pcols, pver, ncol, t, pm, pi, rhoi, ni, ti, nm, sph)
This subroutine computes profiles of background state quantities for the multiple gravity wave drag p...
Definition m_micro.F90:1913
subroutine, public m_micro_run(im, lm, flipv, dt_i, prsl_i, prsi_i, phil, phii, omega_i, qlls_i, qlcn_i, qils_i, qicn_i, lwheat_i, swheat_i, w_upi, cf_upi, frland, zpbl, cnv_mfd_i, cnv_dqldt_i, clcn_i, u_i, v_i, taugwx, taugwy, tauorox, tauoroy, cnv_fice_i, cnv_ndrop_i, cnv_nice_i, q_io, lwm_o, qi_o, t_io, rn_o, sr_o, ncpl_io, ncpi_io, fprcp, rnw_io, snw_io, qgl_io, ncpr_io, ncps_io, ncgl_io, clls_io, kcbl, rainmin, cldreffl, cldreffi, cldreffr, cldreffs, cldreffg, ntrcaer, aerfld_i, naai_i, npccn_i, iccn, skip_macro, alf_fac, qc_min, pdfflag, kdt, xlat, xlon, rhc_i, errmsg, errflg)
Definition m_micro.F90:165
This module contains the CCPP-compliant Morrison-Gettelman microphysics (MG1, MG2 and MG3) scheme.
Definition m_micro.F90:8
MG microphysics version 3.0 - Update of MG microphysics with prognostic hail OR graupel.
This module contains some of the most frequently used math and physics constants for GCM models.
Definition physcons.F90:36
Definition sflx.f:3