ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
Functions/Subroutines | Variables
module_mp_morr_two_moment Module Reference

Functions/Subroutines

subroutine morr_two_moment_init (morr_rimed_ice, morr_noice)
 
subroutine, public set_morrison_ndcnst (ndcnst_in)
 
subroutine, public mp_morr_two_moment (ITIMESTEP, TH, QV, QC, QR, QI, QS, QG, NI, NS, NR, NG, RHO, PII, P, DT_IN, DZ, W, RAINNC, RAINNCV, SR, SNOWNC, SNOWNCV, GRAUPELNC, GRAUPELNCV, refl_10cm, diagflag, do_radar_ref, qrcuten, qscuten, qicuten, F_QNDROP, qndrop, IDS, IDE, JDS, JDE, KDS, KDE, IMS, IME, JMS, JME, KMS, KME, ITS, ITE, JTS, JTE, KTS, KTE, wetscav_on, rainprod, evapprod, QLSINK, PRECR, PRECI, PRECS, PRECG)
 
subroutine, private morr_two_moment_micro (i, j, istep, kts, kte, QC3DTEN, QI3DTEN, QNI3DTEN, QR3DTEN, NI3DTEN, NS3DTEN, NR3DTEN, QC3D, QI3D, QNI3D, QR3D, NI3D, NS3D, NR3D, T3DTEN, QV3DTEN, T3D, QV3D, PRES, DZQ, W3D, PRECRT, SNOWRT, SNOWPRT, GRPLPRT, EFFC, EFFI, EFFS, EFFR, DT, QG3DTEN, NG3DTEN, QG3D, NG3D, EFFG, qrcu1d, qscu1d, qicu1d, QGSTEN, QRSTEN, QISTEN, QNISTEN, QCSTEN, nc3d, nc3dten, iinum, c2prec, CSED, ISED, SSED, GSED, RSED)
 
real(c_double) function, public polysvp (T, TYPE)
 
real(c_double) function, private gamma (X)
 

Variables

real(c_double), parameter, private pi = 3.1415926535897932384626434D0
 
real(c_double), parameter, private xxx = 0.9189385332046727417803297D0
 
integer, private iact
 
integer, private inum
 
real(c_double), private ndcnst
 
integer, private iliq
 
integer, private inuc
 
integer, private ibase
 
integer, private isub
 
integer, private igraup
 
integer, private ihail
 
real(c_double), private ai
 
real(c_double), private ac
 
real(c_double), private as
 
real(c_double), private ar
 
real(c_double), private ag
 
real(c_double), private bi
 
real(c_double), private bc
 
real(c_double), private bs
 
real(c_double), private br
 
real(c_double), private bg
 
real(c_double), private rhosu
 
real(c_double), private rhow
 
real(c_double), private rhoi
 
real(c_double), private rhosn
 
real(c_double), private rhog
 
real(c_double), private aimm
 
real(c_double), private bimm
 
real(c_double), private ecr
 
real(c_double), private dcs
 
real(c_double), private mi0
 
real(c_double), private mg0
 
real(c_double), private f1s
 
real(c_double), private f2s
 
real(c_double), private f1r
 
real(c_double), private f2r
 
real(c_double), private qsmall
 
real(c_double), private ci
 
real(c_double), private di
 
real(c_double), private cs
 
real(c_double), private ds
 
real(c_double), private cg
 
real(c_double), private dg
 
real(c_double), private eii
 
real(c_double), private eci
 
real(c_double), private rin
 
real(c_double), private cpw
 
real(c_double), private c1
 
real(c_double), private k1
 
real(c_double), private mw
 
real(c_double), private osm
 
real(c_double), private vi
 
real(c_double), private epsm
 
real(c_double), private rhoa
 
real(c_double), private map
 
real(c_double), private ma
 
real(c_double), private rr
 
real(c_double), private bact
 
real(c_double), private rm1
 
real(c_double), private rm2
 
real(c_double), private nanew1
 
real(c_double), private nanew2
 
real(c_double), private sig1
 
real(c_double), private sig2
 
real(c_double), private f11
 
real(c_double), private f12
 
real(c_double), private f21
 
real(c_double), private f22
 
real(c_double), private mmult
 
real(c_double), private lammaxi
 
real(c_double), private lammini
 
real(c_double), private lammaxr
 
real(c_double), private lamminr
 
real(c_double), private lammaxs
 
real(c_double), private lammins
 
real(c_double), private lammaxg
 
real(c_double), private lamming
 
real(c_double), private cons1
 
real(c_double), private cons2
 
real(c_double), private cons3
 
real(c_double), private cons4
 
real(c_double), private cons5
 
real(c_double), private cons6
 
real(c_double), private cons7
 
real(c_double), private cons8
 
real(c_double), private cons9
 
real(c_double), private cons10
 
real(c_double), private cons11
 
real(c_double), private cons12
 
real(c_double), private cons13
 
real(c_double), private cons14
 
real(c_double), private cons15
 
real(c_double), private cons16
 
real(c_double), private cons17
 
real(c_double), private cons18
 
real(c_double), private cons19
 
real(c_double), private cons20
 
real(c_double), private cons21
 
real(c_double), private cons22
 
real(c_double), private cons23
 
real(c_double), private cons24
 
real(c_double), private cons25
 
real(c_double), private cons26
 
real(c_double), private cons27
 
real(c_double), private cons28
 
real(c_double), private cons29
 
real(c_double), private cons30
 
real(c_double), private cons31
 
real(c_double), private cons32
 
real(c_double), private cons33
 
real(c_double), private cons34
 
real(c_double), private cons35
 
real(c_double), private cons36
 
real(c_double), private cons37
 
real(c_double), private cons38
 
real(c_double), private cons39
 
real(c_double), private cons40
 
real(c_double), private cons41
 

Function/Subroutine Documentation

◆ gamma()

real(c_double) function, private module_mp_morr_two_moment::gamma ( real(c_double)  X)
private
4108 !----------------------------------------------------------------------
4109 !
4110 ! THIS ROUTINE CALCULATES THE GAMMA FUNCTION FOR A REAL(C_DOUBLE) ARGUMENT X.
4111 ! COMPUTATION IS BASED ON AN ALGORITHM OUTLINED IN REFERENCE 1.
4112 ! THE PROGRAM USES RATIONAL FUNCTIONS THAT APPROXIMATE THE GAMMA
4113 ! FUNCTION TO AT LEAST 20 SIGNIFICANT DECIMAL DIGITS. COEFFICIENTS
4114 ! FOR THE APPROXIMATION OVER THE INTERVAL (1,2) ARE UNPUBLISHED.
4115 ! THOSE FOR THE APPROXIMATION FOR X .GE. 12 ARE FROM REFERENCE 2.
4116 ! THE ACCURACY ACHIEVED DEPENDS ON THE ARITHMETIC SYSTEM, THE
4117 ! COMPILER, THE INTRINSIC FUNCTIONS, AND PROPER SELECTION OF THE
4118 ! MACHINE-DEPENDENT CONSTANTS.
4119 !
4120 !
4121 !*******************************************************************
4122 !*******************************************************************
4123 !
4124 ! EXPLANATION OF MACHINE-DEPENDENT CONSTANTS
4125 !
4126 ! BETA - RADIX FOR THE FLOATING-POINT REPRESENTATION
4127 ! MAXEXP - THE SMALLEST POSITIVE POWER OF BETA THAT OVERFLOWS
4128 ! XBIG - THE LARGEST ARGUMENT FOR WHICH GAMMA(X) IS REPRESENTABLE
4129 ! IN THE MACHINE, I.E., THE SOLUTION TO THE EQUATION
4130 ! GAMMA(XBIG) = BETA**MAXEXP
4131 ! XINF - THE LARGEST MACHINE REPRESENTABLE FLOATING-POINT NUMBER;
4132 ! APPROXIMATELY BETA**MAXEXP
4133 ! EPS - THE SMALLEST POSITIVE FLOATING-POINT NUMBER SUCH THAT
4134 ! 1.0+EPS .GT. 1.0
4135 ! XMININ - THE SMALLEST POSITIVE FLOATING-POINT NUMBER SUCH THAT
4136 ! 1/XMININ IS MACHINE REPRESENTABLE
4137 !
4138 ! APPROXIMATE VALUES FOR SOME IMPORTANT MACHINES ARE:
4139 !
4140 ! BETA MAXEXP XBIG
4141 !
4142 ! CRAY-1 (S.P.) 2 8191 966.961
4143 ! CYBER 180/855
4144 ! UNDER NOS (S.P.) 2 1070 177.803
4145 ! IEEE (IBM/XT,
4146 ! SUN, ETC.) (S.P.) 2 128 35.040
4147 ! IEEE (IBM/XT,
4148 ! SUN, ETC.) (D.P.) 2 1024 171.624
4149 ! IBM 3033 (D.P.) 16 63 57.574
4150 ! VAX D-FORMAT (D.P.) 2 127 34.844
4151 ! VAX G-FORMAT (D.P.) 2 1023 171.489
4152 !
4153 ! XINF EPS XMININ
4154 !
4155 ! CRAY-1 (S.P.) 5.45E+2465 7.11E-15 1.84E-2466
4156 ! CYBER 180/855
4157 ! UNDER NOS (S.P.) 1.26E+322 3.55E-15 3.14E-294
4158 ! IEEE (IBM/XT,
4159 ! SUN, ETC.) (S.P.) 3.40E+38 1.19E-7 1.18E-38
4160 ! IEEE (IBM/XT,
4161 ! SUN, ETC.) (D.P.) 1.79D+308 2.22D-16 2.23D-308
4162 ! IBM 3033 (D.P.) 7.23D+75 2.22D-16 1.39D-76
4163 ! VAX D-FORMAT (D.P.) 1.70D+38 1.39D-17 5.88D-39
4164 ! VAX G-FORMAT (D.P.) 8.98D+307 1.11D-16 1.12D-308
4165 !
4166 !*******************************************************************
4167 !*******************************************************************
4168 !
4169 ! ERROR RETURNS
4170 !
4171 ! THE PROGRAM RETURNS THE VALUE XINF FOR SINGULARITIES OR
4172 ! WHEN OVERFLOW WOULD OCCUR. THE COMPUTATION IS BELIEVED
4173 ! TO BE FREE OF UNDERFLOW AND OVERFLOW.
4174 !
4175 !
4176 ! INTRINSIC FUNCTIONS REQUIRED ARE:
4177 !
4178 ! INT, DBLE, EXP, LOG, REAL(C_DOUBLE), SIN
4179 !
4180 !
4181 ! REFERENCES: AN OVERVIEW OF SOFTWARE DEVELOPMENT FOR SPECIAL
4182 ! FUNCTIONS W. J. CODY, LECTURE NOTES IN MATHEMATICS,
4183 ! 506, NUMERICAL ANALYSIS DUNDEE, 1975, G. A. WATSON
4184 ! (ED.), SPRINGER VERLAG, BERLIN, 1976.
4185 !
4186 ! COMPUTER APPROXIMATIONS, HART, ET. AL., WILEY AND
4187 ! SONS, NEW YORK, 1968.
4188 !
4189 ! LATEST MODIFICATION: OCTOBER 12, 1989
4190 !
4191 ! AUTHORS: W. J. CODY AND L. STOLTZ
4192 ! APPLIED MATHEMATICS DIVISION
4193 ! ARGONNE NATIONAL LABORATORY
4194 ! ARGONNE, IL 60439
4195 !
4196 !----------------------------------------------------------------------
4197  implicit none
4198  INTEGER I,N
4199  LOGICAL PARITY
4200  REAL(C_DOUBLE) &
4201  conv,eps,fact,half,one,res,sum,twelve, &
4202  two,x,xbig,xden,xinf,xminin,xnum,y,y1,ysq,z,zero
4203  REAL(C_DOUBLE), DIMENSION(7) :: C
4204  REAL(C_DOUBLE), DIMENSION(8) :: P
4205  REAL(C_DOUBLE), DIMENSION(8) :: Q
4206 !----------------------------------------------------------------------
4207 ! MATHEMATICAL CONSTANTS
4208 !----------------------------------------------------------------------
4209  DATA one,half,twelve,two,zero/1.0d0,0.5d0,12.0d0,2.0d0,0.0d0/
4210 
4211 
4212 !----------------------------------------------------------------------
4213 ! MACHINE DEPENDENT PARAMETERS
4214 !----------------------------------------------------------------------
4215  DATA xbig,xminin,eps/35.040d0,1.18d-38,1.19d-7/,xinf/3.4d38/
4216 !----------------------------------------------------------------------
4217 ! NUMERATOR AND DENOMINATOR COEFFICIENTS FOR RATIONAL MINIMAX
4218 ! APPROXIMATION OVER (1,2).
4219 !----------------------------------------------------------------------
4220  DATA p/-1.71618513886549492533811d+0,2.47656508055759199108314d+1, &
4221  -3.79804256470945635097577d+2,6.29331155312818442661052d+2, &
4222  8.66966202790413211295064d+2,-3.14512729688483675254357d+4, &
4223  -3.61444134186911729807069d+4,6.64561438202405440627855d+4/
4224  DATA q/-3.08402300119738975254353d+1,3.15350626979604161529144d+2, &
4225  -1.01515636749021914166146d+3,-3.10777167157231109440444d+3, &
4226  2.25381184209801510330112d+4,4.75584627752788110767815d+3, &
4227  -1.34659959864969306392456d+5,-1.15132259675553483497211d+5/
4228 !----------------------------------------------------------------------
4229 ! COEFFICIENTS FOR MINIMAX APPROXIMATION OVER (12, INF).
4230 !----------------------------------------------------------------------
4231  DATA c/-1.910444077728d-03,8.4171387781295d-04, &
4232  -5.952379913043012d-04,7.93650793500350248d-04, &
4233  -2.777777777777681622553d-03,8.333333333333333331554247d-02, &
4234  5.7083835261d-03/
4235 !----------------------------------------------------------------------
4236 ! STATEMENT FUNCTIONS FOR CONVERSION BETWEEN INTEGER AND FLOAT
4237 !----------------------------------------------------------------------
4238  conv(i) = real(i)
4239  parity=.false.
4240  fact=one
4241  n=0
4242  y=x
4243  IF(y.LE.zero)THEN
4244 !----------------------------------------------------------------------
4245 ! ARGUMENT IS NEGATIVE
4246 !----------------------------------------------------------------------
4247  y=-x
4248  y1=aint(y)
4249  res=y-y1
4250  IF(res.NE.zero)THEN
4251  IF(y1.NE.aint(y1*half)*two)parity=.true.
4252  fact=-pi/sin(pi*res)
4253  y=y+one
4254  ELSE
4255  res=xinf
4256  GOTO 900
4257  ENDIF
4258  ENDIF
4259 !----------------------------------------------------------------------
4260 ! ARGUMENT IS POSITIVE
4261 !----------------------------------------------------------------------
4262  IF(y.LT.eps)THEN
4263 !----------------------------------------------------------------------
4264 ! ARGUMENT .LT. EPS
4265 !----------------------------------------------------------------------
4266  IF(y.GE.xminin)THEN
4267  res=one/y
4268  ELSE
4269  res=xinf
4270  GOTO 900
4271  ENDIF
4272  ELSEIF(y.LT.twelve)THEN
4273  y1=y
4274  IF(y.LT.one)THEN
4275 !----------------------------------------------------------------------
4276 ! 0.0 .LT. ARGUMENT .LT. 1.0
4277 !----------------------------------------------------------------------
4278  z=y
4279  y=y+one
4280  ELSE
4281 !----------------------------------------------------------------------
4282 ! 1.0 .LT. ARGUMENT .LT. 12.0, REDUCE ARGUMENT IF NECESSARY
4283 !----------------------------------------------------------------------
4284  n=int(y)-1
4285  y=y-conv(n)
4286  z=y-one
4287  ENDIF
4288 !----------------------------------------------------------------------
4289 ! EVALUATE APPROXIMATION FOR 1.0 .LT. ARGUMENT .LT. 2.0
4290 !----------------------------------------------------------------------
4291  xnum=zero
4292  xden=one
4293  DO i=1,8
4294  xnum=(xnum+p(i))*z
4295  xden=xden*z+q(i)
4296  END DO
4297  res=xnum/xden+one
4298  IF(y1.LT.y)THEN
4299 !----------------------------------------------------------------------
4300 ! ADJUST RESULT FOR CASE 0.0 .LT. ARGUMENT .LT. 1.0
4301 !----------------------------------------------------------------------
4302  res=res/y1
4303  ELSEIF(y1.GT.y)THEN
4304 !----------------------------------------------------------------------
4305 ! ADJUST RESULT FOR CASE 2.0 .LT. ARGUMENT .LT. 12.0
4306 !----------------------------------------------------------------------
4307  DO i=1,n
4308  res=res*y
4309  y=y+one
4310  END DO
4311  ENDIF
4312  ELSE
4313 !----------------------------------------------------------------------
4314 ! EVALUATE FOR ARGUMENT .GE. 12.0,
4315 !----------------------------------------------------------------------
4316  IF(y.LE.xbig)THEN
4317  ysq=y*y
4318  sum=c(7)
4319  DO i=1,6
4320  sum=sum/ysq+c(i)
4321  END DO
4322  sum=sum/y-y+xxx
4323  sum=sum+(y-half)*log(y)
4324  res=exp(sum)
4325  ELSE
4326  res=xinf
4327  GOTO 900
4328  ENDIF
4329  ENDIF
4330 !----------------------------------------------------------------------
4331 ! FINAL ADJUSTMENTS AND RETURN
4332 !----------------------------------------------------------------------
4333  IF(parity)res=-res
4334  IF(fact.NE.one)res=fact/res
4335  900 gamma=res
4336  RETURN
4337 ! ---------- LAST LINE OF GAMMA ----------
real
double real
Definition: ERF_OrbCosZenith.H:9

◆ morr_two_moment_init()

subroutine module_mp_morr_two_moment::morr_two_moment_init ( integer, intent(in)  morr_rimed_ice,
integer, intent(in)  morr_noice 
)
252 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
253 ! THIS SUBROUTINE INITIALIZES ALL PHYSICAL CONSTANTS AMND PARAMETERS
254 ! NEEDED BY THE MICROPHYSICS SCHEME.
255 ! NEEDS TO BE CALLED AT FIRST TIME STEP, PRIOR TO CALL TO MAIN MICROPHYSICS INTERFACE
256 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
257 
258  IMPLICIT NONE
259 
260  INTEGER, INTENT(IN):: morr_rimed_ice ! RAS
261  INTEGER, INTENT(IN):: morr_noice
262 
263  integer n,i
264 
265  print *,'IN MORR_TWO_MOMENT_INIT ',"with noice = ",morr_noice
266 
267 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
268 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
269 
270 ! THE FOLLOWING PARAMETERS ARE USER-DEFINED SWITCHES AND NEED TO BE
271 ! SET PRIOR TO CODE COMPILATION
272 
273 ! INUM IS AUTOMATICALLY SET TO 0 FOR WRF-CHEM BELOW,
274 ! ALLOWING PREDICTION OF DROPLET CONCENTRATION
275 ! THUS, THIS PARAMETER SHOULD NOT BE CHANGED HERE
276 ! AND SHOULD BE LEFT TO 1
277 
278  inum = 1
279 
280 ! SET CONSTANT DROPLET CONCENTRATION (UNITS OF CM-3)
281 ! IF NO COUPLING WITH WRF-CHEM
282 
283  ndcnst = 250.
284 
285 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
286 ! NOTE, THE FOLLOWING OPTIONS RELATED TO DROPLET ACTIVATION
287 ! (IACT, IBASE, ISUB) ARE NOT AVAILABLE IN CURRENT VERSION
288 ! FOR WRF-CHEM, DROPLET ACTIVATION IS PERFORMED
289 ! IN 'MIX_ACTIVATE', NOT IN MICROPHYSICS SCHEME
290 
291 
292 ! IACT = 1, USE POWER-LAW CCN SPECTRA, NCCN = CS^K
293 ! IACT = 2, USE LOGNORMAL AEROSOL SIZE DIST TO DERIVE CCN SPECTRA
294 
295  iact = 2
296 
297 ! IBASE = 1, NEGLECT DROPLET ACTIVATION AT LATERAL CLOUD EDGES DUE TO
298 ! UNRESOLVED ENTRAINMENT AND MIXING, ACTIVATE
299 ! AT CLOUD BASE OR IN REGION WITH LITTLE CLOUD WATER USING
300 ! NON-EQULIBRIUM SUPERSATURATION ASSUMING NO INITIAL CLOUD WATER,
301 ! IN CLOUD INTERIOR ACTIVATE USING EQUILIBRIUM SUPERSATURATION
302 ! IBASE = 2, ASSUME DROPLET ACTIVATION AT LATERAL CLOUD EDGES DUE TO
303 ! UNRESOLVED ENTRAINMENT AND MIXING DOMINATES,
304 ! ACTIVATE DROPLETS EVERYWHERE IN THE CLOUD USING NON-EQUILIBRIUM
305 ! SUPERSATURATION ASSUMING NO INITIAL CLOUD WATER, BASED ON THE
306 ! LOCAL SUB-GRID AND/OR GRID-SCALE VERTICAL VELOCITY
307 ! AT THE GRID POINT
308 
309 ! NOTE: ONLY USED FOR PREDICTED DROPLET CONCENTRATION (INUM = 0)
310 
311  ibase = 2
312 
313 ! INCLUDE SUB-GRID VERTICAL VELOCITY (standard deviation of w) IN DROPLET ACTIVATION
314 ! ISUB = 0, INCLUDE SUB-GRID W (RECOMMENDED FOR LOWER RESOLUTION)
315 ! currently, sub-grid w is constant of 0.5 m/s (not coupled with PBL/turbulence scheme)
316 ! ISUB = 1, EXCLUDE SUB-GRID W, ONLY USE GRID-SCALE W
317 
318 ! NOTE: ONLY USED FOR PREDICTED DROPLET CONCENTRATION (INUM = 0)
319 
320  isub = 0
321 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
322  if(morr_noice .eq. 1) then
323 
324 ! SWITCH FOR LIQUID-ONLY RUN
325 ! ILIQ = 0, INCLUDE ICE
326 ! ILIQ = 1, LIQUID ONLY, NO ICE
327 
328  iliq = 1
329 
330 ! SWITCH FOR ICE NUCLEATION
331 ! INUC = 0, USE FORMULA FROM RASMUSSEN ET AL. 2002 (MID-LATITUDE)
332 ! = 1, USE MPACE OBSERVATIONS (ARCTIC ONLY)
333 
334  inuc = 0
335 
336 ! SWITCH FOR GRAUPEL/HAIL NO GRAUPEL/HAIL
337 ! IGRAUP = 0, INCLUDE GRAUPEL/HAIL
338 ! IGRAUP = 1, NO GRAUPEL/HAIL
339 
340  igraup = 1
341 
342 ! HM ADDED 11/7/07
343 ! SWITCH FOR HAIL/GRAUPEL
344 ! IHAIL = 0, DENSE PRECIPITATING ICE IS GRAUPEL
345 ! IHAIL = 1, DENSE PRECIPITATING ICE IS HAIL
346 ! NOTE ---> RECOMMEND IHAIL = 1 FOR CONTINENTAL DEEP CONVECTION
347 
348  !IHAIL = 0 !changed to namelist option (morr_rimed_ice) by RAS
349  ! Check if namelist option is feasible, otherwise default to graupel - RAS
350  IF (morr_rimed_ice .eq. 1) THEN
351  ihail = 1
352  ELSE
353  ihail = 0
354  ENDIF
355 else
356 
357 ! SWITCH FOR LIQUID-ONLY RUN
358 ! ILIQ = 0, INCLUDE ICE
359 ! ILIQ = 1, LIQUID ONLY, NO ICE
360 
361  iliq = 0
362 
363 ! SWITCH FOR ICE NUCLEATION
364 ! INUC = 0, USE FORMULA FROM RASMUSSEN ET AL. 2002 (MID-LATITUDE)
365 ! = 1, USE MPACE OBSERVATIONS (ARCTIC ONLY)
366 
367  inuc = 0
368 
369 ! SWITCH FOR GRAUPEL/HAIL NO GRAUPEL/HAIL
370 ! IGRAUP = 0, INCLUDE GRAUPEL/HAIL
371 ! IGRAUP = 1, NO GRAUPEL/HAIL
372 
373  igraup = 0
374 
375 ! HM ADDED 11/7/07
376 ! SWITCH FOR HAIL/GRAUPEL
377 ! IHAIL = 0, DENSE PRECIPITATING ICE IS GRAUPEL
378 ! IHAIL = 1, DENSE PRECIPITATING ICE IS HAIL
379 ! NOTE ---> RECOMMEND IHAIL = 1 FOR CONTINENTAL DEEP CONVECTION
380 
381  !IHAIL = 0 !changed to namelist option (morr_rimed_ice) by RAS
382  ! Check if namelist option is feasible, otherwise default to graupel - RAS
383  IF (morr_rimed_ice .eq. 1) THEN
384  ihail = 1
385  ELSE
386  ihail = 0
387  ENDIF
388 endif
389 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
390 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
391 ! SET PHYSICAL CONSTANTS
392 
393 ! FALLSPEED PARAMETERS (V=AD^B)
394  ai = 700.
395  ac = 3.e7
396  as = 11.72
397  ar = 841.99667
398  bi = 1.
399  bc = 2.
400  bs = 0.41
401  br = 0.8
402  IF (ihail.EQ.0) THEN
403  ag = 19.3
404  bg = 0.37
405  ELSE ! (MATSUN AND HUGGINS 1980)
406  ag = 114.5
407  bg = 0.5
408  END IF
409 
410 ! CONSTANTS AND PARAMETERS
411 ! R = 287.15
412 ! RV = 461.5
413 ! CP = 1005.
414  rhosu = 85000./(287.15*273.15)
415  rhow = 997.
416  rhoi = 500.
417  rhosn = 100.
418  IF (ihail.EQ.0) THEN
419  rhog = 400.
420  ELSE
421  rhog = 900.
422  END IF
423  aimm = 0.66
424  bimm = 100.
425  ecr = 1.
426  dcs = 125.e-6
427  mi0 = 4./3.*pi*rhoi*(10.e-6)**3
428  mg0 = 1.6e-10
429  f1s = 0.86
430  f2s = 0.28
431  f1r = 0.78
432 ! F2R = 0.32
433 ! fix 053011
434  f2r = 0.308
435 ! G = 9.806
436 
437  qsmall = 1.d-14
438  print *,'SETTING SMALL TO ',qsmall
439 
440  eii = 0.1
441  eci = 0.7
442 ! HM, ADD FOR V3.2
443 ! hm, 7/23/13
444 ! CPW = 4218.
445  cpw = 4187.
446 
447 ! SIZE DISTRIBUTION PARAMETERS
448 
449  ci = rhoi*pi/6.
450  di = 3.
451  cs = rhosn*pi/6.
452  ds = 3.
453  cg = rhog*pi/6.
454  dg = 3.
455 
456 ! RADIUS OF CONTACT NUCLEI
457  rin = 0.1e-6
458 
459  mmult = 4./3.*pi*rhoi*(5.e-6)**3
460 
461 ! SIZE LIMITS FOR LAMBDA
462 
463  lammaxi = 1./1.e-6
464  lammini = 1./(2.*dcs+100.e-6)
465  lammaxr = 1./20.e-6
466 ! LAMMINR = 1./500.E-6
467  lamminr = 1./2800.e-6
468 
469  lammaxs = 1./10.e-6
470  lammins = 1./2000.e-6
471  lammaxg = 1./20.e-6
472  lamming = 1./2000.e-6
473 
474 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
475 ! note: these parameters only used by the non-wrf-chem version of the
476 ! scheme with predicted droplet number
477 
478 ! CCN SPECTRA FOR IACT = 1
479 
480 ! MARITIME
481 ! MODIFIED FROM RASMUSSEN ET AL. 2002
482 ! NCCN = C*S^K, NCCN IS IN CM-3, S IS SUPERSATURATION RATIO IN %
483 
484  k1 = 0.4
485  c1 = 120.
486 
487 ! CONTINENTAL
488 
489 ! K1 = 0.5
490 ! C1 = 1000.
491 
492 ! AEROSOL ACTIVATION PARAMETERS FOR IACT = 2
493 ! PARAMETERS CURRENTLY SET FOR AMMONIUM SULFATE
494 
495  mw = 0.018
496  osm = 1.
497  vi = 3.
498  epsm = 0.7
499  rhoa = 1777.
500  map = 0.132
501  ma = 0.0284
502 ! hm fix 6/23/16
503 ! RR = 8.3187
504  rr = 8.3145
505  bact = vi*osm*epsm*mw*rhoa/(map*rhow)
506 
507 ! AEROSOL SIZE DISTRIBUTION PARAMETERS CURRENTLY SET FOR MPACE
508 ! (see morrison et al. 2007, JGR)
509 ! MODE 1
510 
511  rm1 = 0.052e-6
512  sig1 = 2.04
513  nanew1 = 72.2e6
514  f11 = 0.5*exp(2.5*(log(sig1))**2)
515  f21 = 1.+0.25*log(sig1)
516 
517 ! MODE 2
518 
519  rm2 = 1.3e-6
520  sig2 = 2.5
521  nanew2 = 1.8e6
522  f12 = 0.5*exp(2.5*(log(sig2))**2)
523  f22 = 1.+0.25*log(sig2)
524 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
525 
526 ! CONSTANTS FOR EFFICIENCY
527 
528  cons1=gamma(1.+ds)*cs
529  cons2=gamma(1.+dg)*cg
530  cons3=gamma(4.+bs)/6.
531  cons4=gamma(4.+br)/6.
532  cons5=gamma(1.+bs)
533  cons6=gamma(1.+br)
534  cons7=gamma(4.+bg)/6.
535  cons8=gamma(1.+bg)
536  cons9=gamma(5./2.+br/2.)
537  cons10=gamma(5./2.+bs/2.)
538  cons11=gamma(5./2.+bg/2.)
539  cons12=gamma(1.+di)*ci
540  cons13=gamma(bs+3.)*pi/4.*eci
541  cons14=gamma(bg+3.)*pi/4.*eci
542  cons15=-1108.*eii*pi**((1.-bs)/3.)*rhosn**((-2.-bs)/3.)/(4.*720.)
543  cons16=gamma(bi+3.)*pi/4.*eci
544  cons17=4.*2.*3.*rhosu*pi*eci*eci*gamma(2.*bs+2.)/(8.*(rhog-rhosn))
545  cons18=rhosn*rhosn
546  cons19=rhow*rhow
547  cons20=20.*pi*pi*rhow*bimm
548  cons21=4./(dcs*rhoi)
549  cons22=pi*rhoi*dcs**3/6.
550  cons23=pi/4.*eii*gamma(bs+3.)
551  cons24=pi/4.*ecr*gamma(br+3.)
552  cons25=pi*pi/24.*rhow*ecr*gamma(br+6.)
553  cons26=pi/6.*rhow
554  cons27=gamma(1.+bi)
555  cons28=gamma(4.+bi)/6.
556  cons29=4./3.*pi*rhow*(25.e-6)**3
557  cons30=4./3.*pi*rhow
558  cons31=pi*pi*ecr*rhosn
559  cons32=pi/2.*ecr
560  cons33=pi*pi*ecr*rhog
561  cons34=5./2.+br/2.
562  cons35=5./2.+bs/2.
563  cons36=5./2.+bg/2.
564  cons37=4.*pi*1.38e-23/(6.*pi*rin)
565  cons38=pi*pi/3.*rhow
566  cons39=pi*pi/36.*rhow*bimm
567  cons40=pi/6.*bimm
568  cons41=pi*pi*ecr*rhow
569 

Referenced by mp_morr_two_moment_isohelper::morr_two_moment_init_c().

Here is the caller graph for this function:

◆ morr_two_moment_micro()

subroutine, private module_mp_morr_two_moment::morr_two_moment_micro ( integer, intent(in)  i,
integer, intent(in)  j,
integer, intent(in)  istep,
integer, intent(in)  kts,
integer, intent(in)  kte,
real(c_double), dimension(kts:kte)  QC3DTEN,
real(c_double), dimension(kts:kte)  QI3DTEN,
real(c_double), dimension(kts:kte)  QNI3DTEN,
real(c_double), dimension(kts:kte)  QR3DTEN,
real(c_double), dimension(kts:kte)  NI3DTEN,
real(c_double), dimension(kts:kte)  NS3DTEN,
real(c_double), dimension(kts:kte)  NR3DTEN,
real(c_double), dimension(kts:kte)  QC3D,
real(c_double), dimension(kts:kte)  QI3D,
real(c_double), dimension(kts:kte)  QNI3D,
real(c_double), dimension(kts:kte)  QR3D,
real(c_double), dimension(kts:kte)  NI3D,
real(c_double), dimension(kts:kte)  NS3D,
real(c_double), dimension(kts:kte)  NR3D,
real(c_double), dimension(kts:kte)  T3DTEN,
real(c_double), dimension(kts:kte)  QV3DTEN,
real(c_double), dimension(kts:kte)  T3D,
real(c_double), dimension(kts:kte)  QV3D,
real(c_double), dimension(kts:kte)  PRES,
real(c_double), dimension(kts:kte)  DZQ,
real(c_double), dimension(kts:kte)  W3D,
real(c_double)  PRECRT,
real(c_double)  SNOWRT,
real(c_double)  SNOWPRT,
real(c_double)  GRPLPRT,
real(c_double), dimension(kts:kte)  EFFC,
real(c_double), dimension(kts:kte)  EFFI,
real(c_double), dimension(kts:kte)  EFFS,
real(c_double), dimension(kts:kte)  EFFR,
real(c_double)  DT,
real(c_double), dimension(kts:kte)  QG3DTEN,
real(c_double), dimension(kts:kte)  NG3DTEN,
real(c_double), dimension(kts:kte)  QG3D,
real(c_double), dimension(kts:kte)  NG3D,
real(c_double), dimension(kts:kte)  EFFG,
real(c_double), dimension(kts:kte)  qrcu1d,
real(c_double), dimension(kts:kte)  qscu1d,
real(c_double), dimension(kts:kte)  qicu1d,
real(c_double), dimension(kts:kte)  QGSTEN,
real(c_double), dimension(kts:kte)  QRSTEN,
real(c_double), dimension(kts:kte)  QISTEN,
real(c_double), dimension(kts:kte)  QNISTEN,
real(c_double), dimension(kts:kte)  QCSTEN,
real(c_double), dimension(kts:kte)  nc3d,
real(c_double), dimension(kts:kte)  nc3dten,
integer, intent(in)  iinum,
real(c_double), dimension(kts:kte)  c2prec,
real(c_double), dimension(kts:kte)  CSED,
real(c_double), dimension(kts:kte)  ISED,
real(c_double), dimension(kts:kte)  SSED,
real(c_double), dimension(kts:kte)  GSED,
real(c_double), dimension(kts:kte)  RSED 
)
private
933 
934 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
935 ! THIS PROGRAM IS THE MAIN TWO-MOMENT MICROPHYSICS SUBROUTINE DESCRIBED BY
936 ! MORRISON ET AL. 2005 JAS AND MORRISON ET AL. 2009 MWR
937 
938 ! THIS SCHEME IS A BULK DOUBLE-MOMENT SCHEME THAT PREDICTS MIXING
939 ! RATIOS AND NUMBER CONCENTRATIONS OF FIVE HYDROMETEOR SPECIES:
940 ! CLOUD DROPLETS, CLOUD (SMALL) ICE, RAIN, SNOW, AND GRAUPEL/HAIL.
941 
942 ! CODE STRUCTURE: MAIN SUBROUTINE IS 'MORR_TWO_MOMENT'. ALSO INCLUDED IN THIS FILE IS
943 ! 'FUNCTION POLYSVP', 'FUNCTION DERF1', AND
944 ! 'FUNCTION GAMMA'.
945 
946 ! NOTE: THIS SUBROUTINE USES 1D ARRAY IN VERTICAL (COLUMN), EVEN THOUGH VARIABLES ARE CALLED '3D'......
947 
948 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
949 
950 ! DECLARATIONS
951 
952  IMPLICIT NONE
953 
954 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
955 ! THESE VARIABLES BELOW MUST BE LINKED WITH THE MAIN MODEL.
956 ! DEFINE ARRAY SIZES
957 
958 ! INPUT NUMBER OF GRID CELLS
959 
960 ! INPUT/OUTPUT PARAMETERS ! DESCRIPTION (UNITS)
961  INTEGER, INTENT( IN) :: i,j,istep,kts,kte
962 
963  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QC3DTEN ! CLOUD WATER MIXING RATIO TENDENCY (KG/KG/S)
964  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QI3DTEN ! CLOUD ICE MIXING RATIO TENDENCY (KG/KG/S)
965  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QNI3DTEN ! SNOW MIXING RATIO TENDENCY (KG/KG/S)
966  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QR3DTEN ! RAIN MIXING RATIO TENDENCY (KG/KG/S)
967  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NI3DTEN ! CLOUD ICE NUMBER CONCENTRATION (1/KG/S)
968  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NS3DTEN ! SNOW NUMBER CONCENTRATION (1/KG/S)
969  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NR3DTEN ! RAIN NUMBER CONCENTRATION (1/KG/S)
970  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QC3D ! CLOUD WATER MIXING RATIO (KG/KG)
971  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QI3D ! CLOUD ICE MIXING RATIO (KG/KG)
972  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QNI3D ! SNOW MIXING RATIO (KG/KG)
973  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QR3D ! RAIN MIXING RATIO (KG/KG)
974  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NI3D ! CLOUD ICE NUMBER CONCENTRATION (1/KG)
975  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NS3D ! SNOW NUMBER CONCENTRATION (1/KG)
976  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NR3D ! RAIN NUMBER CONCENTRATION (1/KG)
977  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: T3DTEN ! TEMPERATURE TENDENCY (K/S)
978  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QV3DTEN ! WATER VAPOR MIXING RATIO TENDENCY (KG/KG/S)
979  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: T3D ! TEMPERATURE (K)
980  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QV3D ! WATER VAPOR MIXING RATIO (KG/KG)
981  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PRES ! ATMOSPHERIC PRESSURE (PA)
982  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: DZQ ! DIFFERENCE IN HEIGHT ACROSS LEVEL (m)
983  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: W3D ! GRID-SCALE VERTICAL VELOCITY (M/S)
984 ! below for wrf-chem
985  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: nc3d
986  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: nc3dten
987  integer, intent(in) :: iinum
988 
989 ! HM ADDED GRAUPEL VARIABLES
990  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QG3DTEN ! GRAUPEL MIX RATIO TENDENCY (KG/KG/S)
991  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NG3DTEN ! GRAUPEL NUMB CONC TENDENCY (1/KG/S)
992  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QG3D ! GRAUPEL MIX RATIO (KG/KG)
993  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NG3D ! GRAUPEL NUMBER CONC (1/KG)
994 
995 ! HM, ADD 1/16/07, SEDIMENTATION TENDENCIES FOR MIXING RATIO
996 
997  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QGSTEN ! GRAUPEL SED TEND (KG/KG/S)
998  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QRSTEN ! RAIN SED TEND (KG/KG/S)
999  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QISTEN ! CLOUD ICE SED TEND (KG/KG/S)
1000  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QNISTEN ! SNOW SED TEND (KG/KG/S)
1001  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QCSTEN ! CLOUD WAT SED TEND (KG/KG/S)
1002 
1003 ! hm add cumulus tendencies for precip
1004  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: qrcu1d
1005  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: qscu1d
1006  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: qicu1d
1007 
1008 ! OUTPUT VARIABLES
1009 
1010  REAL(C_DOUBLE) PRECRT ! TOTAL PRECIP PER TIME STEP (mm)
1011  REAL(C_DOUBLE) SNOWRT ! SNOW PER TIME STEP (mm)
1012 ! hm added 7/13/13
1013  REAL(C_DOUBLE) SNOWPRT ! TOTAL CLOUD ICE PLUS SNOW PER TIME STEP (mm)
1014  REAL(C_DOUBLE) GRPLPRT ! TOTAL GRAUPEL PER TIME STEP (mm)
1015 
1016  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: EFFC ! DROPLET EFFECTIVE RADIUS (MICRON)
1017  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: EFFI ! CLOUD ICE EFFECTIVE RADIUS (MICRON)
1018  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: EFFS ! SNOW EFFECTIVE RADIUS (MICRON)
1019  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: EFFR ! RAIN EFFECTIVE RADIUS (MICRON)
1020  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: EFFG ! GRAUPEL EFFECTIVE RADIUS (MICRON)
1021 
1022 ! MODEL INPUT PARAMETERS (FORMERLY IN COMMON BLOCKS)
1023 
1024  REAL(C_DOUBLE) DT ! MODEL TIME STEP (SEC)
1025 
1026 
1027 !.....................................................................................................
1028 ! LOCAL VARIABLES: ALL PARAMETERS BELOW ARE LOCAL TO SCHEME AND DON'T NEED TO COMMUNICATE WITH THE
1029 ! REST OF THE MODEL.
1030 
1031 ! SIZE PARAMETER VARIABLES
1032 
1033  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: LAMC ! SLOPE PARAMETER FOR DROPLETS (M-1)
1034  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: LAMI ! SLOPE PARAMETER FOR CLOUD ICE (M-1)
1035  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: LAMS ! SLOPE PARAMETER FOR SNOW (M-1)
1036  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: LAMR ! SLOPE PARAMETER FOR RAIN (M-1)
1037  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: LAMG ! SLOPE PARAMETER FOR GRAUPEL (M-1)
1038  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: CDIST1 ! PSD PARAMETER FOR DROPLETS
1039  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: N0I ! INTERCEPT PARAMETER FOR CLOUD ICE (KG-1 M-1)
1040  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: N0S ! INTERCEPT PARAMETER FOR SNOW (KG-1 M-1)
1041  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: N0RR ! INTERCEPT PARAMETER FOR RAIN (KG-1 M-1)
1042  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: N0G ! INTERCEPT PARAMETER FOR GRAUPEL (KG-1 M-1)
1043  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PGAM ! SPECTRAL SHAPE PARAMETER FOR DROPLETS
1044 
1045 ! MICROPHYSICAL PROCESSES
1046 
1047  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NSUBC ! LOSS OF NC DURING EVAP
1048  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NSUBI ! LOSS OF NI DURING SUB.
1049  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NSUBS ! LOSS OF NS DURING SUB.
1050  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NSUBR ! LOSS OF NR DURING EVAP
1051  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PRD ! DEP CLOUD ICE
1052  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PRE ! EVAP OF RAIN
1053  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PRDS ! DEP SNOW
1054  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NNUCCC ! CHANGE N DUE TO CONTACT FREEZ DROPLETS
1055  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: MNUCCC ! CHANGE Q DUE TO CONTACT FREEZ DROPLETS
1056  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PRA ! ACCRETION DROPLETS BY RAIN
1057  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PRC ! AUTOCONVERSION DROPLETS
1058  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PCC ! COND/EVAP DROPLETS
1059  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NNUCCD ! CHANGE N FREEZING AEROSOL (PRIM ICE NUCLEATION)
1060  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: MNUCCD ! CHANGE Q FREEZING AEROSOL (PRIM ICE NUCLEATION)
1061  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: MNUCCR ! CHANGE Q DUE TO CONTACT FREEZ RAIN
1062  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NNUCCR ! CHANGE N DUE TO CONTACT FREEZ RAIN
1063  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NPRA ! CHANGE IN N DUE TO DROPLET ACC BY RAIN
1064  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NRAGG ! SELF-COLLECTION/BREAKUP OF RAIN
1065  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NSAGG ! SELF-COLLECTION OF SNOW
1066  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NPRC ! CHANGE NC AUTOCONVERSION DROPLETS
1067  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NPRC1 ! CHANGE NR AUTOCONVERSION DROPLETS
1068  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PRAI ! CHANGE Q ACCRETION CLOUD ICE BY SNOW
1069  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PRCI ! CHANGE Q AUTOCONVERSIN CLOUD ICE TO SNOW
1070  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PSACWS ! CHANGE Q DROPLET ACCRETION BY SNOW
1071  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NPSACWS ! CHANGE N DROPLET ACCRETION BY SNOW
1072  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PSACWI ! CHANGE Q DROPLET ACCRETION BY CLOUD ICE
1073  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NPSACWI ! CHANGE N DROPLET ACCRETION BY CLOUD ICE
1074  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NPRCI ! CHANGE N AUTOCONVERSION CLOUD ICE BY SNOW
1075  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NPRAI ! CHANGE N ACCRETION CLOUD ICE
1076  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NMULTS ! ICE MULT DUE TO RIMING DROPLETS BY SNOW
1077  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NMULTR ! ICE MULT DUE TO RIMING RAIN BY SNOW
1078  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QMULTS ! CHANGE Q DUE TO ICE MULT DROPLETS/SNOW
1079  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QMULTR ! CHANGE Q DUE TO ICE RAIN/SNOW
1080  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PRACS ! CHANGE Q RAIN-SNOW COLLECTION
1081  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NPRACS ! CHANGE N RAIN-SNOW COLLECTION
1082  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PCCN ! CHANGE Q DROPLET ACTIVATION
1083  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PSMLT ! CHANGE Q MELTING SNOW TO RAIN
1084  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: EVPMS ! CHNAGE Q MELTING SNOW EVAPORATING
1085  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NSMLTS ! CHANGE N MELTING SNOW
1086  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NSMLTR ! CHANGE N MELTING SNOW TO RAIN
1087 ! HM ADDED 12/13/06
1088  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PIACR ! CHANGE QR, ICE-RAIN COLLECTION
1089  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NIACR ! CHANGE N, ICE-RAIN COLLECTION
1090  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PRACI ! CHANGE QI, ICE-RAIN COLLECTION
1091  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PIACRS ! CHANGE QR, ICE RAIN COLLISION, ADDED TO SNOW
1092  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NIACRS ! CHANGE N, ICE RAIN COLLISION, ADDED TO SNOW
1093  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PRACIS ! CHANGE QI, ICE RAIN COLLISION, ADDED TO SNOW
1094  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: EPRD ! SUBLIMATION CLOUD ICE
1095  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: EPRDS ! SUBLIMATION SNOW
1096 ! HM ADDED GRAUPEL PROCESSES
1097  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PRACG ! CHANGE IN Q COLLECTION RAIN BY GRAUPEL
1098  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PSACWG ! CHANGE IN Q COLLECTION DROPLETS BY GRAUPEL
1099  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PGSACW ! CONVERSION Q TO GRAUPEL DUE TO COLLECTION DROPLETS BY SNOW
1100  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PGRACS ! CONVERSION Q TO GRAUPEL DUE TO COLLECTION RAIN BY SNOW
1101  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PRDG ! DEP OF GRAUPEL
1102  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: EPRDG ! SUB OF GRAUPEL
1103  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: EVPMG ! CHANGE Q MELTING OF GRAUPEL AND EVAPORATION
1104  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PGMLT ! CHANGE Q MELTING OF GRAUPEL
1105  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NPRACG ! CHANGE N COLLECTION RAIN BY GRAUPEL
1106  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NPSACWG ! CHANGE N COLLECTION DROPLETS BY GRAUPEL
1107  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NSCNG ! CHANGE N CONVERSION TO GRAUPEL DUE TO COLLECTION DROPLETS BY SNOW
1108  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NGRACS ! CHANGE N CONVERSION TO GRAUPEL DUE TO COLLECTION RAIN BY SNOW
1109  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NGMLTG ! CHANGE N MELTING GRAUPEL
1110  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NGMLTR ! CHANGE N MELTING GRAUPEL TO RAIN
1111  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NSUBG ! CHANGE N SUB/DEP OF GRAUPEL
1112  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: PSACR ! CONVERSION DUE TO COLL OF SNOW BY RAIN
1113  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NMULTG ! ICE MULT DUE TO ACC DROPLETS BY GRAUPEL
1114  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NMULTRG ! ICE MULT DUE TO ACC RAIN BY GRAUPEL
1115  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QMULTG ! CHANGE Q DUE TO ICE MULT DROPLETS/GRAUPEL
1116  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QMULTRG ! CHANGE Q DUE TO ICE MULT RAIN/GRAUPEL
1117 
1118 ! TIME-VARYING ATMOSPHERIC PARAMETERS
1119 
1120  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: KAP ! THERMAL CONDUCTIVITY OF AIR
1121  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: EVS ! SATURATION VAPOR PRESSURE
1122  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: EIS ! ICE SATURATION VAPOR PRESSURE
1123  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QVS ! SATURATION MIXING RATIO
1124  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QVI ! ICE SATURATION MIXING RATIO
1125  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QVQVS ! SAUTRATION RATIO
1126  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: QVQVSI! ICE SATURAION RATIO
1127  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: DV ! DIFFUSIVITY OF WATER VAPOR IN AIR
1128  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: XXLS ! LATENT HEAT OF SUBLIMATION
1129  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: XXLV ! LATENT HEAT OF VAPORIZATION
1130  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: CPM ! SPECIFIC HEAT AT CONST PRESSURE FOR MOIST AIR
1131  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: MU ! VISCOCITY OF AIR
1132  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: SC ! SCHMIDT NUMBER
1133  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: XLF ! LATENT HEAT OF FREEZING
1134  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: RHO ! AIR DENSITY
1135  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: AB ! CORRECTION TO CONDENSATION RATE DUE TO LATENT HEATING
1136  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: ABI ! CORRECTION TO DEPOSITION RATE DUE TO LATENT HEATING
1137 
1138 ! TIME-VARYING MICROPHYSICS PARAMETERS
1139 
1140  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: DAP ! DIFFUSIVITY OF AEROSOL
1141  REAL(C_DOUBLE) NACNT ! NUMBER OF CONTACT IN
1142  REAL(C_DOUBLE) FMULT ! TEMP.-DEP. PARAMETER FOR RIME-SPLINTERING
1143  REAL(C_DOUBLE) COFFI ! ICE AUTOCONVERSION PARAMETER
1144 
1145 ! FALL SPEED WORKING VARIABLES (DEFINED IN CODE)
1146 
1147  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: DUMI,DUMR,DUMFNI,DUMG,DUMFNG
1148  REAL(C_DOUBLE) UNI, UMI,UMR
1149  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: FR, FI, FNI,FG,FNG
1150  REAL(C_DOUBLE) RGVM
1151  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: FALOUTR,FALOUTI,FALOUTNI
1152  REAL(C_DOUBLE) FALTNDR,FALTNDI,FALTNDNI,RHO2
1153  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: DUMQS,DUMFNS
1154  REAL(C_DOUBLE) UMS,UNS
1155  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: FS,FNS, FALOUTS,FALOUTNS,FALOUTG,FALOUTNG
1156  REAL(C_DOUBLE) FALTNDS,FALTNDNS,UNR,FALTNDG,FALTNDNG
1157  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: DUMC,DUMFNC
1158  REAL(C_DOUBLE) UNC,UMC,UNG,UMG
1159  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: FC,FALOUTC,FALOUTNC
1160  REAL(C_DOUBLE) FALTNDC,FALTNDNC
1161  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: FNC,DUMFNR,FALOUTNR
1162  REAL(C_DOUBLE) FALTNDNR
1163  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: FNR
1164 
1165 ! FALL-SPEED PARAMETER 'A' WITH AIR DENSITY CORRECTION
1166 
1167  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: AIN,ARN,ASN,ACN,AGN
1168 
1169 ! EXTERNAL FUNCTION CALL RETURN VARIABLES
1170 
1171 ! REAL(C_DOUBLE) GAMMA, ! EULER GAMMA FUNCTION
1172 ! REAL(C_DOUBLE) POLYSVP, ! SAT. PRESSURE FUNCTION
1173 ! REAL(C_DOUBLE) DERF1 ! ERROR FUNCTION
1174 
1175 ! DUMMY VARIABLES
1176 
1177  REAL(C_DOUBLE) DUM,DUM1,DUM2,DUMT,DUMQV,DUMQSS,DUMQSI,DUMS
1178 
1179 ! PROGNOSTIC SUPERSATURATION
1180 
1181  REAL(C_DOUBLE) DQSDT ! CHANGE OF SAT. MIX. RAT. WITH TEMPERATURE
1182  REAL(C_DOUBLE) DQSIDT ! CHANGE IN ICE SAT. MIXING RAT. WITH T
1183  REAL(C_DOUBLE) EPSI ! 1/PHASE REL. TIME (SEE M2005), ICE
1184  REAL(C_DOUBLE) EPSS ! 1/PHASE REL. TIME (SEE M2005), SNOW
1185  REAL(C_DOUBLE) EPSR ! 1/PHASE REL. TIME (SEE M2005), RAIN
1186  REAL(C_DOUBLE) EPSG ! 1/PHASE REL. TIME (SEE M2005), GRAUPEL
1187 
1188 ! NEW DROPLET ACTIVATION VARIABLES
1189  REAL(C_DOUBLE) TAUC ! PHASE REL. TIME (SEE M2005), DROPLETS
1190  REAL(C_DOUBLE) TAUR ! PHASE REL. TIME (SEE M2005), RAIN
1191  REAL(C_DOUBLE) TAUI ! PHASE REL. TIME (SEE M2005), CLOUD ICE
1192  REAL(C_DOUBLE) TAUS ! PHASE REL. TIME (SEE M2005), SNOW
1193  REAL(C_DOUBLE) TAUG ! PHASE REL. TIME (SEE M2005), GRAUPEL
1194  REAL(C_DOUBLE) DUMACT,DUM3
1195 
1196 ! COUNTING/INDEX VARIABLES
1197 
1198  INTEGER K,NSTEP,N ! ,I
1199 
1200 ! LTRUE IS ONLY USED TO SPEED UP THE CODE !!
1201 ! LTRUE, SWITCH = 0, NO HYDROMETEORS IN COLUMN,
1202 ! = 1, HYDROMETEORS IN COLUMN
1203 
1204  INTEGER LTRUE
1205 
1206 ! DROPLET ACTIVATION/FREEZING AEROSOL
1207 
1208 
1209  REAL(C_DOUBLE) CT ! DROPLET ACTIVATION PARAMETER
1210  REAL(C_DOUBLE) TEMP1 ! DUMMY TEMPERATURE
1211  REAL(C_DOUBLE) SAT1 ! DUMMY SATURATION
1212  REAL(C_DOUBLE) SIGVL ! SURFACE TENSION LIQ/VAPOR
1213  REAL(C_DOUBLE) KEL ! KELVIN PARAMETER
1214  REAL(C_DOUBLE) KC2 ! TOTAL ICE NUCLEATION RATE
1215 
1216  REAL(C_DOUBLE) CRY,KRY ! AEROSOL ACTIVATION PARAMETERS
1217 
1218 ! MORE WORKING/DUMMY VARIABLES
1219 
1220  REAL(C_DOUBLE) DUMQI,DUMNI,DC0,DS0,DG0
1221  REAL(C_DOUBLE) DUMQC,DUMQR,RATIO,SUM_DEP,FUDGEF
1222 
1223 ! EFFECTIVE VERTICAL VELOCITY (M/S)
1224  REAL(C_DOUBLE) WEF
1225 
1226 ! WORKING PARAMETERS FOR ICE NUCLEATION
1227 
1228  REAL(C_DOUBLE) ANUC,BNUC
1229 
1230 ! WORKING PARAMETERS FOR AEROSOL ACTIVATION
1231 
1232  REAL(C_DOUBLE) AACT,GAMM,GG,PSI,ETA1,ETA2,SM1,SM2,SMAX,UU1,UU2,ALPHA
1233 
1234 ! DUMMY SIZE DISTRIBUTION PARAMETERS
1235 
1236  REAL(C_DOUBLE) DLAMS,DLAMR,DLAMI,DLAMC,DLAMG,LAMMAX,LAMMIN
1237 
1238  INTEGER IDROP
1239 
1240 ! FOR WRF-CHEM
1241  REAL(C_DOUBLE), DIMENSION(KTS:KTE)::C2PREC,CSED,ISED,SSED,GSED,RSED
1242  REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: tqimelt ! melting of cloud ice (tendency)
1243 
1244 ! comment lines for wrf-chem since these are intent(in) in that case
1245 ! REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NC3DTEN ! CLOUD DROPLET NUMBER CONCENTRATION (1/KG/S)
1246 ! REAL(C_DOUBLE), DIMENSION(KTS:KTE) :: NC3D ! CLOUD DROPLET NUMBER CONCENTRATION (1/KG)
1247 
1248 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
1249 
1250  ! print *,'IN MICRO KTS:KTE ', i,j,kts, kte
1251 
1252 ! SET LTRUE INITIALLY TO 0
1253 
1254  ltrue = 0
1255 
1256 ! ATMOSPHERIC PARAMETERS THAT VARY IN TIME AND HEIGHT
1257  DO k = kts,kte
1258 
1259 ! NC3DTEN LOCAL ARRAY INITIALIZED
1260  nc3dten(k) = 0.
1261 ! INITIALIZE VARIABLES FOR WRF-CHEM OUTPUT TO ZERO
1262 
1263  c2prec(k)=0.
1264  csed(k)=0.
1265  ised(k)=0.
1266  ssed(k)=0.
1267  gsed(k)=0.
1268  rsed(k)=0.
1269 
1270 ! LATENT HEAT OF VAPORATION
1271 
1272  xxlv(k) = 3.1484e6-2370.*t3d(k)
1273 
1274 ! LATENT HEAT OF SUBLIMATION
1275 
1276  xxls(k) = 3.15e6-2370.*t3d(k)+0.3337e6
1277 
1278  cpm(k) = cp*(1.+0.887*qv3d(k))
1279 
1280 ! SATURATION VAPOR PRESSURE AND MIXING RATIO
1281 
1282 ! hm, add fix for low pressure, 5/12/10
1283 
1284  evs(k) = min(0.99*pres(k),polysvp(t3d(k),0)) ! PA
1285  eis(k) = min(0.99*pres(k),polysvp(t3d(k),1)) ! PA
1286 ! Fortran version
1287 ! MAKE SURE ICE SATURATION DOESN'T EXCEED WATER SAT. NEAR FREEZING
1288 
1289  IF (eis(k).GT.evs(k)) eis(k) = evs(k)
1290 
1291  qvs(k) = ep_2*evs(k)/(pres(k)-evs(k))
1292  qvi(k) = ep_2*eis(k)/(pres(k)-eis(k))
1293 
1294  qvqvs(k) = qv3d(k)/qvs(k)
1295  qvqvsi(k) = qv3d(k)/qvi(k)
1296 
1297 ! AIR DENSITY
1298 
1299  rho(k) = pres(k)/(r*t3d(k))
1300 
1301 ! ADD NUMBER CONCENTRATION DUE TO CUMULUS TENDENCY
1302 ! ASSUME N0 ASSOCIATED WITH CUMULUS PARAM RAIN IS 10^7 M^-4
1303 ! ASSUME N0 ASSOCIATED WITH CUMULUS PARAM SNOW IS 2 X 10^7 M^-4
1304 ! FOR DETRAINED CLOUD ICE, ASSUME MEAN VOLUME DIAM OF 80 MICRON
1305 
1306  IF (qrcu1d(k).GE.1.e-10) THEN
1307  dum=1.8e5*(qrcu1d(k)*dt/(pi*rhow*rho(k)**3))**0.25
1308  nr3d(k)=nr3d(k)+dum
1309  END IF
1310  IF (qscu1d(k).GE.1.e-10) THEN
1311  dum=3.e5*(qscu1d(k)*dt/(cons1*rho(k)**3))**(1./(ds+1.))
1312  ns3d(k)=ns3d(k)+dum
1313  END IF
1314  IF (qicu1d(k).GE.1.e-10) THEN
1315  dum=qicu1d(k)*dt/(ci*(80.e-6)**di)
1316  ni3d(k)=ni3d(k)+dum
1317  END IF
1318 
1319 ! AT SUBSATURATION, REMOVE SMALL AMOUNTS OF CLOUD/PRECIP WATER
1320 ! hm modify 7/0/09 change limit to 1.e-8
1321 
1322  IF (qvqvs(k).LT.0.9) THEN
1323  IF (qr3d(k).LT.1.e-8) THEN
1324  qv3d(k)=qv3d(k)+qr3d(k)
1325  t3d(k)=t3d(k)-qr3d(k)*xxlv(k)/cpm(k)
1326  qr3d(k)=0.
1327  END IF
1328  IF (qc3d(k).LT.1.e-8) THEN
1329  qv3d(k)=qv3d(k)+qc3d(k)
1330  t3d(k)=t3d(k)-qc3d(k)*xxlv(k)/cpm(k)
1331  qc3d(k)=0.
1332  END IF
1333  END IF
1334 ! Fortran version
1335  IF (qvqvsi(k).LT.0.9) THEN
1336  IF (qi3d(k).LT.1.e-8) THEN
1337  qv3d(k)=qv3d(k)+qi3d(k)
1338  t3d(k)=t3d(k)-qi3d(k)*xxls(k)/cpm(k)
1339  qi3d(k)=0.
1340  END IF
1341  IF (qni3d(k).LT.1.e-8) THEN
1342  qv3d(k)=qv3d(k)+qni3d(k)
1343  t3d(k)=t3d(k)-qni3d(k)*xxls(k)/cpm(k)
1344  qni3d(k)=0.
1345  END IF
1346  IF (qg3d(k).LT.1.e-8) THEN
1347  qv3d(k)=qv3d(k)+qg3d(k)
1348  t3d(k)=t3d(k)-qg3d(k)*xxls(k)/cpm(k)
1349  qg3d(k)=0.
1350  END IF
1351  END IF
1352  ! Fortran version
1353 ! HEAT OF FUSION
1354 
1355  xlf(k) = xxls(k)-xxlv(k)
1356 
1357 !..................................................................
1358 ! IF MIXING RATIO < QSMALL SET MIXING RATIO AND NUMBER CONC TO ZERO
1359 
1360  IF (qc3d(k).LT.qsmall) THEN
1361  qc3d(k) = 0.
1362  nc3d(k) = 0.
1363  effc(k) = 0.
1364  END IF
1365  IF (qr3d(k).LT.qsmall) THEN
1366  qr3d(k) = 0.
1367  nr3d(k) = 0.
1368  effr(k) = 0.
1369  END IF
1370  IF (qi3d(k).LT.qsmall) THEN
1371  qi3d(k) = 0.
1372  ni3d(k) = 0.
1373  effi(k) = 0.
1374  END IF
1375  IF (qni3d(k).LT.qsmall) THEN
1376  qni3d(k) = 0.
1377  ns3d(k) = 0.
1378  effs(k) = 0.
1379  END IF
1380  IF (qg3d(k).LT.qsmall) THEN
1381  qg3d(k) = 0.
1382  ng3d(k) = 0.
1383  effg(k) = 0.
1384  END IF
1385 ! Fortran version
1386 ! INITIALIZE SEDIMENTATION TENDENCIES FOR MIXING RATIO
1387 
1388  qrsten(k) = 0.
1389  qisten(k) = 0.
1390  qnisten(k) = 0.
1391  qcsten(k) = 0.
1392  qgsten(k) = 0.
1393 
1394 !..................................................................
1395 ! MICROPHYSICS PARAMETERS VARYING IN TIME/HEIGHT
1396 
1397 ! fix 053011
1398  mu(k) = 1.496e-6*t3d(k)**1.5/(t3d(k)+120.)
1399 
1400 ! FALL SPEED WITH DENSITY CORRECTION (HEYMSFIELD AND BENSSEMER 2006)
1401 
1402  dum = (rhosu/rho(k))**0.54
1403 
1404 ! fix 053011
1405 ! AIN(K) = DUM*AI
1406 ! AA revision 4/1/11: Ikawa and Saito 1991 air-density correction
1407  ain(k) = (rhosu/rho(k))**0.35*ai
1408  arn(k) = dum*ar
1409  asn(k) = dum*as
1410 ! ACN(K) = DUM*AC
1411 ! AA revision 4/1/11: temperature-dependent Stokes fall speed
1412  acn(k) = g*rhow/(18.*mu(k))
1413 ! HM ADD GRAUPEL 8/28/06
1414  agn(k) = dum*ag
1415 !hm 4/7/09 bug fix, initialize lami to prevent later division by zero
1416  lami(k)=0.
1417 
1418 !..................................
1419 ! IF THERE IS NO CLOUD/PRECIP WATER, AND IF SUBSATURATED, THEN SKIP MICROPHYSICS
1420 ! FOR THIS LEVEL
1421 
1422  IF ( qc3d(k).LT.qsmall.AND. &
1423  qi3d(k).LT.qsmall.AND. &
1424  qni3d(k).LT.qsmall.AND. &
1425  qr3d(k).LT.qsmall.AND. &
1426  qg3d(k).LT.qsmall) THEN
1427  IF (t3d(k).LT.273.15.AND.qvqvsi(k).LT.0.999) then
1428  GOTO 200
1429  endif
1430  IF (t3d(k).GE.273.15.AND.qvqvs(k).LT.0.999) then
1431  GOTO 200
1432  endif
1433  END IF
1434 
1435 ! THERMAL CONDUCTIVITY FOR AIR
1436 
1437 ! fix 053011
1438  kap(k) = 1.414e3*mu(k)
1439 
1440 ! DIFFUSIVITY OF WATER VAPOR
1441 
1442  dv(k) = 8.794e-5*t3d(k)**1.81/pres(k)
1443 
1444 ! SCHMIT NUMBER
1445 
1446 ! fix 053011
1447  sc(k) = mu(k)/(rho(k)*dv(k))
1448 
1449 ! PSYCHOMETIC CORRECTIONS
1450 
1451 ! RATE OF CHANGE SAT. MIX. RATIO WITH TEMPERATURE
1452 
1453  dum = (rv*t3d(k)**2)
1454 
1455  dqsdt = xxlv(k)*qvs(k)/dum
1456  dqsidt = xxls(k)*qvi(k)/dum
1457 
1458  abi(k) = 1.+dqsidt*xxls(k)/cpm(k)
1459  ab(k) = 1.+dqsdt*xxlv(k)/cpm(k)
1460 !
1461 !.....................................................................
1462 !.....................................................................
1463 ! CASE FOR TEMPERATURE ABOVE FREEZING
1464 
1465  IF (t3d(k).GE.273.15) THEN
1466 
1467 !......................................................................
1468 !HM ADD, ALLOW FOR CONSTANT DROPLET NUMBER
1469 ! INUM = 0, PREDICT DROPLET NUMBER
1470 ! INUM = 1, SET CONSTANT DROPLET NUMBER
1471 
1472  IF (iinum.EQ.1) THEN
1473 ! CONVERT NDCNST FROM CM-3 TO KG-1
1474  nc3d(k)=ndcnst*1.e6/rho(k)
1475  END IF
1476 
1477 ! GET SIZE DISTRIBUTION PARAMETERS
1478 
1479 ! MELT VERY SMALL SNOW AND GRAUPEL MIXING RATIOS, ADD TO RAIN
1480  IF (qni3d(k).LT.1.e-6) THEN
1481  qr3d(k)=qr3d(k)+qni3d(k)
1482  nr3d(k)=nr3d(k)+ns3d(k)
1483  t3d(k)=t3d(k)-qni3d(k)*xlf(k)/cpm(k)
1484  qni3d(k) = 0.
1485  ns3d(k) = 0.
1486  END IF
1487  IF (qg3d(k).LT.1.e-6) THEN
1488  qr3d(k)=qr3d(k)+qg3d(k)
1489  nr3d(k)=nr3d(k)+ng3d(k)
1490  t3d(k)=t3d(k)-qg3d(k)*xlf(k)/cpm(k)
1491  qg3d(k) = 0.
1492  ng3d(k) = 0.
1493  END IF
1494  IF (qc3d(k).LT.qsmall.AND.qni3d(k).LT.1.e-8.AND.qr3d(k).LT.qsmall.AND.qg3d(k).LT.1.e-8) THEN
1495  GOTO 300
1496  ENDIF
1497 
1498 ! MAKE SURE NUMBER CONCENTRATIONS AREN'T NEGATIVE
1499 
1500  ns3d(k) = max(0.,ns3d(k))
1501  nc3d(k) = max(0.,nc3d(k))
1502  nr3d(k) = max(0.,nr3d(k))
1503  ng3d(k) = max(0.,ng3d(k))
1504 
1505 !......................................................................
1506 ! RAIN
1507 
1508  IF (qr3d(k).GE.qsmall) THEN
1509  lamr(k) = (pi*rhow*nr3d(k)/qr3d(k))**(1./3.)
1510  n0rr(k) = nr3d(k)*lamr(k)
1511 
1512 ! CHECK FOR SLOPE
1513 
1514 ! ADJUST VARS
1515 
1516  IF (lamr(k).LT.lamminr) THEN
1517 
1518  lamr(k) = lamminr
1519 
1520  n0rr(k) = lamr(k)**4*qr3d(k)/(pi*rhow)
1521 
1522  nr3d(k) = n0rr(k)/lamr(k)
1523  ELSE IF (lamr(k).GT.lammaxr) THEN
1524  lamr(k) = lammaxr
1525  n0rr(k) = lamr(k)**4*qr3d(k)/(pi*rhow)
1526 
1527  nr3d(k) = n0rr(k)/lamr(k)
1528  END IF
1529  END IF
1530 
1531 !......................................................................
1532 ! CLOUD DROPLETS
1533 
1534 ! MARTIN ET AL. (1994) FORMULA FOR PGAM
1535 
1536  IF (qc3d(k).GE.qsmall) THEN
1537 
1538  dum = pres(k)/(287.15*t3d(k))
1539  pgam(k)=0.0005714*(nc3d(k)/1.e6*dum)+0.2714
1540  pgam(k)=1./(pgam(k)**2)-1.
1541  pgam(k)=max(pgam(k),2.)
1542  pgam(k)=min(pgam(k),10.)
1543 
1544 ! CALCULATE LAMC
1545 
1546  lamc(k) = (cons26*nc3d(k)*gamma(pgam(k)+4.)/ &
1547  (qc3d(k)*gamma(pgam(k)+1.)))**(1./3.)
1548 
1549 ! LAMMIN, 60 MICRON DIAMETER
1550 ! LAMMAX, 1 MICRON
1551 
1552  lammin = (pgam(k)+1.)/60.e-6
1553  lammax = (pgam(k)+1.)/1.e-6
1554 
1555  IF (lamc(k).LT.lammin) THEN
1556  lamc(k) = lammin
1557 
1558  nc3d(k) = exp(3.*log(lamc(k))+log(qc3d(k))+ &
1559  log(gamma(pgam(k)+1.))-log(gamma(pgam(k)+4.)))/cons26
1560  ELSE IF (lamc(k).GT.lammax) THEN
1561  lamc(k) = lammax
1562 
1563  nc3d(k) = exp(3.*log(lamc(k))+log(qc3d(k))+ &
1564  log(gamma(pgam(k)+1.))-log(gamma(pgam(k)+4.)))/cons26
1565 
1566  END IF
1567 
1568  END IF
1569 
1570 !......................................................................
1571 ! SNOW
1572 
1573  IF (qni3d(k).GE.qsmall) THEN
1574  lams(k) = (cons1*ns3d(k)/qni3d(k))**(1./ds)
1575  n0s(k) = ns3d(k)*lams(k)
1576 
1577 ! CHECK FOR SLOPE
1578 
1579 ! ADJUST VARS
1580 
1581  IF (lams(k).LT.lammins) THEN
1582  lams(k) = lammins
1583  n0s(k) = lams(k)**4*qni3d(k)/cons1
1584 
1585  ns3d(k) = n0s(k)/lams(k)
1586 
1587  ELSE IF (lams(k).GT.lammaxs) THEN
1588 
1589  lams(k) = lammaxs
1590  n0s(k) = lams(k)**4*qni3d(k)/cons1
1591 
1592  ns3d(k) = n0s(k)/lams(k)
1593  END IF
1594  END IF
1595 
1596 !......................................................................
1597 ! GRAUPEL
1598 
1599  IF (qg3d(k).GE.qsmall) THEN
1600  lamg(k) = (cons2*ng3d(k)/qg3d(k))**(1./dg)
1601  n0g(k) = ng3d(k)*lamg(k)
1602 
1603 ! ADJUST VARS
1604 
1605  IF (lamg(k).LT.lamming) THEN
1606  lamg(k) = lamming
1607  n0g(k) = lamg(k)**4*qg3d(k)/cons2
1608 
1609  ng3d(k) = n0g(k)/lamg(k)
1610 
1611  ELSE IF (lamg(k).GT.lammaxg) THEN
1612 
1613  lamg(k) = lammaxg
1614  n0g(k) = lamg(k)**4*qg3d(k)/cons2
1615 
1616  ng3d(k) = n0g(k)/lamg(k)
1617  END IF
1618  END IF
1619 !.....................................................................
1620 ! ZERO OUT PROCESS RATES
1621 
1622  prc(k) = 0.
1623  nprc(k) = 0.
1624  nprc1(k) = 0.
1625  pra(k) = 0.
1626  npra(k) = 0.
1627  nragg(k) = 0.
1628  nsmlts(k) = 0.
1629  nsmltr(k) = 0.
1630  evpms(k) = 0.
1631  pcc(k) = 0.
1632  pre(k) = 0.
1633  nsubc(k) = 0.
1634  nsubr(k) = 0.
1635  pracg(k) = 0.
1636  npracg(k) = 0.
1637  psmlt(k) = 0.
1638  pgmlt(k) = 0.
1639  evpmg(k) = 0.
1640  pracs(k) = 0.
1641  npracs(k) = 0.
1642  ngmltg(k) = 0.
1643  ngmltr(k) = 0.
1644 
1645 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
1646 ! CALCULATION OF MICROPHYSICAL PROCESS RATES, T > 273.15 K
1647 
1648 !.................................................................
1649 !.......................................................................
1650 ! AUTOCONVERSION OF CLOUD LIQUID WATER TO RAIN
1651 ! FORMULA FROM BEHENG (1994)
1652 ! USING NUMERICAL SIMULATION OF STOCHASTIC COLLECTION EQUATION
1653 ! AND INITIAL CLOUD DROPLET SIZE DISTRIBUTION SPECIFIED
1654 ! AS A GAMMA DISTRIBUTION
1655 
1656 ! USE MINIMUM VALUE OF 1.E-6 TO PREVENT FLOATING POINT ERROR
1657 
1658  IF (qc3d(k).GE.1.e-6) THEN
1659 
1660 ! HM ADD 12/13/06, REPLACE WITH NEWER FORMULA
1661 ! FROM KHAIROUTDINOV AND KOGAN 2000, MWR
1662 
1663  prc(k)=1350.*qc3d(k)**2.47* &
1664  (nc3d(k)/1.e6*rho(k))**(-1.79)
1665 
1666 ! note: nprc1 is change in Nr,
1667 ! nprc is change in Nc
1668 
1669  nprc1(k) = prc(k)/cons29
1670  nprc(k) = prc(k)/(qc3d(k)/nc3d(k))
1671 
1672 ! hm bug fix 3/20/12
1673  nprc(k) = min(nprc(k),nc3d(k)/dt)
1674  nprc1(k) = min(nprc1(k),nprc(k))
1675 
1676  END IF
1677 
1678 !.......................................................................
1679 ! HM ADD 12/13/06, COLLECTION OF SNOW BY RAIN ABOVE FREEZING
1680 ! FORMULA FROM IKAWA AND SAITO (1991)
1681 
1682  IF (qr3d(k).GE.1.e-8.AND.qni3d(k).GE.1.e-8) THEN
1683 
1684  ums = asn(k)*cons3/(lams(k)**bs)
1685  umr = arn(k)*cons4/(lamr(k)**br)
1686  uns = asn(k)*cons5/lams(k)**bs
1687  unr = arn(k)*cons6/lamr(k)**br
1688 
1689 ! SET REASLISTIC LIMITS ON FALLSPEEDS
1690 
1691 ! bug fix, 10/08/09
1692  dum=(rhosu/rho(k))**0.54
1693  ums=min(ums,1.2*dum)
1694  uns=min(uns,1.2*dum)
1695  umr=min(umr,9.1*dum)
1696  unr=min(unr,9.1*dum)
1697 
1698 ! hm fix, 2/12/13
1699 ! for above freezing conditions to get accelerated melting of snow,
1700 ! we need collection of rain by snow (following Lin et al. 1983)
1701 ! PRACS(K) = CONS31*(((1.2*UMR-0.95*UMS)**2+ &
1702 ! 0.08*UMS*UMR)**0.5*RHO(K)* &
1703 ! N0RR(K)*N0S(K)/LAMS(K)**3* &
1704 ! (5./(LAMS(K)**3*LAMR(K))+ &
1705 ! 2./(LAMS(K)**2*LAMR(K)**2)+ &
1706 ! 0.5/(LAMS(K)*LAMR(K)**3)))
1707 
1708  pracs(k) = cons41*(((1.2*umr-0.95*ums)**2+ &
1709  0.08*ums*umr)**0.5*rho(k)* &
1710  n0rr(k)*n0s(k)/lamr(k)**3* &
1711  (5./(lamr(k)**3*lams(k))+ &
1712  2./(lamr(k)**2*lams(k)**2)+ &
1713  0.5/(lamr(k)*lams(k)**3)))
1714 
1715 ! fix 053011, npracs no longer subtracted from snow
1716 ! NPRACS(K) = CONS32*RHO(K)*(1.7*(UNR-UNS)**2+ &
1717 ! 0.3*UNR*UNS)**0.5*N0RR(K)*N0S(K)* &
1718 ! (1./(LAMR(K)**3*LAMS(K))+ &
1719 ! 1./(LAMR(K)**2*LAMS(K)**2)+ &
1720 ! 1./(LAMR(K)*LAMS(K)**3))
1721 
1722  END IF
1723 
1724 ! ADD COLLECTION OF GRAUPEL BY RAIN ABOVE FREEZING
1725 ! ASSUME ALL RAIN COLLECTION BY GRAUPEL ABOVE FREEZING IS SHED
1726 ! ASSUME SHED DROPS ARE 1 MM IN SIZE
1727 
1728  IF (qr3d(k).GE.1.e-8.AND.qg3d(k).GE.1.e-8) THEN
1729 
1730  umg = agn(k)*cons7/(lamg(k)**bg)
1731  umr = arn(k)*cons4/(lamr(k)**br)
1732  ung = agn(k)*cons8/lamg(k)**bg
1733  unr = arn(k)*cons6/lamr(k)**br
1734 
1735 ! SET REASLISTIC LIMITS ON FALLSPEEDS
1736 ! bug fix, 10/08/09
1737  dum=(rhosu/rho(k))**0.54
1738  umg=min(umg,20.*dum)
1739  ung=min(ung,20.*dum)
1740  umr=min(umr,9.1*dum)
1741  unr=min(unr,9.1*dum)
1742 
1743 ! PRACG IS MIXING RATIO OF RAIN PER SEC COLLECTED BY GRAUPEL/HAIL
1744  pracg(k) = cons41*(((1.2*umr-0.95*umg)**2+ &
1745  0.08*umg*umr)**0.5*rho(k)* &
1746  n0rr(k)*n0g(k)/lamr(k)**3* &
1747  (5./(lamr(k)**3*lamg(k))+ &
1748  2./(lamr(k)**2*lamg(k)**2)+ &
1749  0.5/(lamr(k)*lamg(k)**3)))
1750 
1751 ! ASSUME 1 MM DROPS ARE SHED, GET NUMBER SHED PER SEC
1752 
1753  dum = pracg(k)/5.2e-7
1754 
1755  npracg(k) = cons32*rho(k)*(1.7*(unr-ung)**2+ &
1756  0.3*unr*ung)**0.5*n0rr(k)*n0g(k)* &
1757  (1./(lamr(k)**3*lamg(k))+ &
1758  1./(lamr(k)**2*lamg(k)**2)+ &
1759  1./(lamr(k)*lamg(k)**3))
1760 
1761 ! hm 7/15/13, remove limit so that the number of collected drops can smaller than
1762 ! number of shed drops
1763 ! NPRACG(K)=MAX(NPRACG(K)-DUM,0.)
1764  npracg(k)=npracg(k)-dum
1765 
1766  END IF
1767 
1768 !.......................................................................
1769 ! ACCRETION OF CLOUD LIQUID WATER BY RAIN
1770 ! CONTINUOUS COLLECTION EQUATION WITH
1771 ! GRAVITATIONAL COLLECTION KERNEL, DROPLET FALL SPEED NEGLECTED
1772 
1773  IF (qr3d(k).GE.1.e-8 .AND. qc3d(k).GE.1.e-8) THEN
1774 
1775 ! 12/13/06 HM ADD, REPLACE WITH NEWER FORMULA FROM
1776 ! KHAIROUTDINOV AND KOGAN 2000, MWR
1777 
1778  dum=(qc3d(k)*qr3d(k))
1779  pra(k) = 67.*(dum)**1.15
1780  npra(k) = pra(k)/(qc3d(k)/nc3d(k))
1781 
1782  END IF
1783 !.......................................................................
1784 ! SELF-COLLECTION OF RAIN DROPS
1785 ! FROM BEHENG(1994)
1786 ! FROM NUMERICAL SIMULATION OF THE STOCHASTIC COLLECTION EQUATION
1787 ! AS DESCRINED ABOVE FOR AUTOCONVERSION
1788 
1789  IF (qr3d(k).GE.1.e-8) THEN
1790 ! include breakup add 10/09/09
1791  dum1=300.e-6
1792  if (1./lamr(k).lt.dum1) then
1793  dum=1.
1794  else if (1./lamr(k).ge.dum1) then
1795  dum=2.-exp(2300.*(1./lamr(k)-dum1))
1796  end if
1797 ! NRAGG(K) = -8.*NR3D(K)*QR3D(K)*RHO(K)
1798  nragg(k) = -5.78*dum*nr3d(k)*qr3d(k)*rho(k)
1799  END IF
1800 
1801 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
1802 ! CALCULATE EVAP OF RAIN (RUTLEDGE AND HOBBS 1983)
1803 
1804  IF (qr3d(k).GE.qsmall) THEN
1805  epsr = 2.*pi*n0rr(k)*rho(k)*dv(k)* &
1806  (f1r/(lamr(k)*lamr(k))+ &
1807  f2r*(arn(k)*rho(k)/mu(k))**0.5* &
1808  sc(k)**(1./3.)*cons9/ &
1809  (lamr(k)**cons34))
1810  ELSE
1811  epsr = 0.
1812  END IF
1813 ! NO CONDENSATION ONTO RAIN, ONLY EVAP ALLOWED
1814 
1815  IF (qv3d(k).LT.qvs(k)) THEN
1816  pre(k) = epsr*(qv3d(k)-qvs(k))/ab(k)
1817  pre(k) = min(pre(k),0.)
1818  ELSE
1819  pre(k) = 0.
1820  END IF
1821 !.......................................................................
1822 ! MELTING OF SNOW
1823 
1824 ! SNOW MAY PERSITS ABOVE FREEZING, FORMULA FROM RUTLEDGE AND HOBBS, 1984
1825 ! IF WATER SUPERSATURATION, SNOW MELTS TO FORM RAIN
1826 
1827  IF (qni3d(k).GE.1.e-8) THEN
1828 
1829 ! fix 053011
1830 ! HM, MODIFY FOR V3.2, ADD ACCELERATED MELTING DUE TO COLLISION WITH RAIN
1831 ! DUM = -CPW/XLF(K)*T3D(K)*PRACS(K)
1832  dum = -cpw/xlf(k)*(t3d(k)-273.15)*pracs(k)
1833 
1834 ! hm fix 1/20/15
1835 ! PSMLT(K)=2.*PI*N0S(K)*KAP(K)*(273.15-T3D(K))/ &
1836 ! XLF(K)*RHO(K)*(F1S/(LAMS(K)*LAMS(K))+ &
1837 ! F2S*(ASN(K)*RHO(K)/MU(K))**0.5* &
1838 ! SC(K)**(1./3.)*CONS10/ &
1839 ! (LAMS(K)**CONS35))+DUM
1840  psmlt(k)=2.*pi*n0s(k)*kap(k)*(273.15-t3d(k))/ &
1841  xlf(k)*(f1s/(lams(k)*lams(k))+ &
1842  f2s*(asn(k)*rho(k)/mu(k))**0.5* &
1843  sc(k)**(1./3.)*cons10/ &
1844  (lams(k)**cons35))+dum
1845 
1846 ! IN WATER SUBSATURATION, SNOW MELTS AND EVAPORATES
1847 
1848  IF (qvqvs(k).LT.1.) THEN
1849  epss = 2.*pi*n0s(k)*rho(k)*dv(k)* &
1850  (f1s/(lams(k)*lams(k))+ &
1851  f2s*(asn(k)*rho(k)/mu(k))**0.5* &
1852  sc(k)**(1./3.)*cons10/ &
1853  (lams(k)**cons35))
1854 ! hm fix 8/4/08
1855  evpms(k) = (qv3d(k)-qvs(k))*epss/ab(k)
1856  evpms(k) = max(evpms(k),psmlt(k))
1857  psmlt(k) = psmlt(k)-evpms(k)
1858  END IF
1859  END IF
1860 
1861 !.......................................................................
1862 ! MELTING OF GRAUPEL
1863 
1864 ! GRAUPEL MAY PERSITS ABOVE FREEZING, FORMULA FROM RUTLEDGE AND HOBBS, 1984
1865 ! IF WATER SUPERSATURATION, GRAUPEL MELTS TO FORM RAIN
1866 
1867  IF (qg3d(k).GE.1.e-8) THEN
1868 
1869 ! fix 053011
1870 ! HM, MODIFY FOR V3.2, ADD ACCELERATED MELTING DUE TO COLLISION WITH RAIN
1871 ! DUM = -CPW/XLF(K)*T3D(K)*PRACG(K)
1872  dum = -cpw/xlf(k)*(t3d(k)-273.15)*pracg(k)
1873 
1874 ! hm fix 1/20/15
1875 ! PGMLT(K)=2.*PI*N0G(K)*KAP(K)*(273.15-T3D(K))/ &
1876 ! XLF(K)*RHO(K)*(F1S/(LAMG(K)*LAMG(K))+ &
1877 ! F2S*(AGN(K)*RHO(K)/MU(K))**0.5* &
1878 ! SC(K)**(1./3.)*CONS11/ &
1879 ! (LAMG(K)**CONS36))+DUM
1880  pgmlt(k)=2.*pi*n0g(k)*kap(k)*(273.15-t3d(k))/ &
1881  xlf(k)*(f1s/(lamg(k)*lamg(k))+ &
1882  f2s*(agn(k)*rho(k)/mu(k))**0.5* &
1883  sc(k)**(1./3.)*cons11/ &
1884  (lamg(k)**cons36))+dum
1885 
1886 ! IN WATER SUBSATURATION, GRAUPEL MELTS AND EVAPORATES
1887 
1888  IF (qvqvs(k).LT.1.) THEN
1889  epsg = 2.*pi*n0g(k)*rho(k)*dv(k)* &
1890  (f1s/(lamg(k)*lamg(k))+ &
1891  f2s*(agn(k)*rho(k)/mu(k))**0.5* &
1892  sc(k)**(1./3.)*cons11/ &
1893  (lamg(k)**cons36))
1894 ! hm fix 8/4/08
1895  evpmg(k) = (qv3d(k)-qvs(k))*epsg/ab(k)
1896  evpmg(k) = max(evpmg(k),pgmlt(k))
1897  pgmlt(k) = pgmlt(k)-evpmg(k)
1898  END IF
1899  END IF
1900 
1901 ! HM, V3.2
1902 ! RESET PRACG AND PRACS TO ZERO, THIS IS DONE BECAUSE THERE IS NO
1903 ! TRANSFER OF MASS FROM SNOW AND GRAUPEL TO RAIN DIRECTLY FROM COLLECTION
1904 ! ABOVE FREEZING, IT IS ONLY USED FOR ENHANCEMENT OF MELTING AND SHEDDING
1905 
1906  pracg(k) = 0.
1907  pracs(k) = 0.
1908 
1909 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
1910 
1911 ! FOR CLOUD ICE, ONLY PROCESSES OPERATING AT T > 273.15 IS
1912 ! MELTING, WHICH IS ALREADY CONSERVED DURING PROCESS
1913 ! CALCULATION
1914 
1915 ! CONSERVATION OF QC
1916 
1917  dum = (prc(k)+pra(k))*dt
1918 
1919  IF (dum.GT.qc3d(k).AND.qc3d(k).GE.qsmall) THEN
1920 
1921  ratio = qc3d(k)/dum
1922 
1923  prc(k) = prc(k)*ratio
1924  pra(k) = pra(k)*ratio
1925 
1926  END IF
1927 
1928 ! CONSERVATION OF SNOW
1929 
1930  dum = (-psmlt(k)-evpms(k)+pracs(k))*dt
1931 
1932  IF (dum.GT.qni3d(k).AND.qni3d(k).GE.qsmall) THEN
1933 
1934 ! NO SOURCE TERMS FOR SNOW AT T > FREEZING
1935  ratio = qni3d(k)/dum
1936 
1937  psmlt(k) = psmlt(k)*ratio
1938  evpms(k) = evpms(k)*ratio
1939  pracs(k) = pracs(k)*ratio
1940 
1941  END IF
1942 
1943 ! CONSERVATION OF GRAUPEL
1944 
1945  dum = (-pgmlt(k)-evpmg(k)+pracg(k))*dt
1946 
1947  IF (dum.GT.qg3d(k).AND.qg3d(k).GE.qsmall) THEN
1948 
1949 ! NO SOURCE TERM FOR GRAUPEL ABOVE FREEZING
1950  ratio = qg3d(k)/dum
1951 
1952  pgmlt(k) = pgmlt(k)*ratio
1953  evpmg(k) = evpmg(k)*ratio
1954  pracg(k) = pracg(k)*ratio
1955 
1956  END IF
1957 
1958 ! CONSERVATION OF QR
1959 ! HM 12/13/06, ADDED CONSERVATION OF RAIN SINCE PRE IS NEGATIVE
1960 
1961  dum = (-pracs(k)-pracg(k)-pre(k)-pra(k)-prc(k)+psmlt(k)+pgmlt(k))*dt
1962 
1963  IF (dum.GT.qr3d(k).AND.qr3d(k).GE.qsmall) THEN
1964 
1965  ratio = (qr3d(k)/dt+pracs(k)+pracg(k)+pra(k)+prc(k)-psmlt(k)-pgmlt(k))/ &
1966  (-pre(k))
1967  pre(k) = pre(k)*ratio
1968  END IF
1969 
1970 !....................................
1971 
1972  qv3dten(k) = qv3dten(k)+(-pre(k)-evpms(k)-evpmg(k))
1973 
1974  t3dten(k) = t3dten(k)+(pre(k)*xxlv(k)+(evpms(k)+evpmg(k))*xxls(k)+&
1975  (psmlt(k)+pgmlt(k)-pracs(k)-pracg(k))*xlf(k))/cpm(k)
1976 
1977  qc3dten(k) = qc3dten(k)+(-pra(k)-prc(k))
1978  qr3dten(k) = qr3dten(k)+(pre(k)+pra(k)+prc(k)-psmlt(k)-pgmlt(k)+pracs(k)+pracg(k))
1979  qni3dten(k) = qni3dten(k)+(psmlt(k)+evpms(k)-pracs(k))
1980  qg3dten(k) = qg3dten(k)+(pgmlt(k)+evpmg(k)-pracg(k))
1981 ! fix 053011
1982 ! NS3DTEN(K) = NS3DTEN(K)-NPRACS(K)
1983 ! HM, bug fix 5/12/08, npracg is subtracted from nr not ng
1984 ! NG3DTEN(K) = NG3DTEN(K)
1985  nc3dten(k) = nc3dten(k)+ (-npra(k)-nprc(k))
1986  nr3dten(k) = nr3dten(k)+ (nprc1(k)+nragg(k)-npracg(k))
1987 ! Fortran version
1988 ! HM ADD, WRF-CHEM, ADD TENDENCIES FOR C2PREC
1989 
1990  c2prec(k) = pra(k)+prc(k)
1991  IF (pre(k).LT.0.) THEN
1992  dum = pre(k)*dt/qr3d(k)
1993  dum = max(-1.,dum)
1994  nsubr(k) = dum*nr3d(k)/dt
1995  END IF
1996 
1997  IF (evpms(k)+psmlt(k).LT.0.) THEN
1998  dum = (evpms(k)+psmlt(k))*dt/qni3d(k)
1999  dum = max(-1.,dum)
2000  nsmlts(k) = dum*ns3d(k)/dt
2001  END IF
2002  IF (psmlt(k).LT.0.) THEN
2003  dum = psmlt(k)*dt/qni3d(k)
2004  dum = max(-1.0,dum)
2005  nsmltr(k) = dum*ns3d(k)/dt
2006  END IF
2007  IF (evpmg(k)+pgmlt(k).LT.0.) THEN
2008  dum = (evpmg(k)+pgmlt(k))*dt/qg3d(k)
2009  dum = max(-1.,dum)
2010  ngmltg(k) = dum*ng3d(k)/dt
2011  END IF
2012  IF (pgmlt(k).LT.0.) THEN
2013  dum = pgmlt(k)*dt/qg3d(k)
2014  dum = max(-1.0,dum)
2015  ngmltr(k) = dum*ng3d(k)/dt
2016  END IF
2017 
2018  ns3dten(k) = ns3dten(k)+(nsmlts(k))
2019  ng3dten(k) = ng3dten(k)+(ngmltg(k))
2020  nr3dten(k) = nr3dten(k)+(nsubr(k)-nsmltr(k)-ngmltr(k))
2021 
2022  300 CONTINUE
2023 
2024 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
2025 ! NOW CALCULATE SATURATION ADJUSTMENT TO CONDENSE EXTRA VAPOR ABOVE
2026 ! WATER SATURATION
2027 
2028  dumt = t3d(k)+dt*t3dten(k)
2029  dumqv = qv3d(k)+dt*qv3dten(k)
2030 ! hm, add fix for low pressure, 5/12/10
2031  dum=min(0.99*pres(k),polysvp(dumt,0))
2032  dumqss = ep_2*dum/(pres(k)-dum)
2033  dumqc = qc3d(k)+dt*qc3dten(k)
2034  dumqc = max(dumqc,0.)
2035 
2036 ! SATURATION ADJUSTMENT FOR LIQUID
2037 
2038  dums = dumqv-dumqss
2039  pcc(k) = dums/(1.+xxlv(k)**2*dumqss/(cpm(k)*rv*dumt**2))/dt
2040  IF (pcc(k)*dt+dumqc.LT.0.) THEN
2041  pcc(k) = -dumqc/dt
2042  END IF
2043 
2044  qv3dten(k) = qv3dten(k)-pcc(k)
2045  t3dten(k) = t3dten(k)+pcc(k)*xxlv(k)/cpm(k)
2046  qc3dten(k) = qc3dten(k)+pcc(k)
2047 !.......................................................................
2048 ! ACTIVATION OF CLOUD DROPLETS
2049 ! ACTIVATION OF DROPLET CURRENTLY NOT CALCULATED
2050 ! DROPLET CONCENTRATION IS SPECIFIED !!!!!
2051 
2052 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
2053 ! SUBLIMATE, MELT, OR EVAPORATE NUMBER CONCENTRATION
2054 ! THIS FORMULATION ASSUMES 1:1 RATIO BETWEEN MASS LOSS AND
2055 ! LOSS OF NUMBER CONCENTRATION
2056 
2057 ! IF (PCC(K).LT.0.) THEN
2058 ! DUM = PCC(K)*DT/QC3D(K)
2059 ! DUM = MAX(-1.,DUM)
2060 ! NSUBC(K) = DUM*NC3D(K)/DT
2061 ! END IF
2062 
2063 ! UPDATE TENDENCIES
2064 
2065 ! NC3DTEN(K) = NC3DTEN(K)+NSUBC(K)
2066 
2067 !.....................................................................
2068 !.....................................................................
2069  ELSE ! TEMPERATURE < 273.15
2070 
2071 !......................................................................
2072 !HM ADD, ALLOW FOR CONSTANT DROPLET NUMBER
2073 ! INUM = 0, PREDICT DROPLET NUMBER
2074 ! INUM = 1, SET CONSTANT DROPLET NUMBER
2075 
2076  IF (iinum.EQ.1) THEN
2077 ! CONVERT NDCNST FROM CM-3 TO KG-1
2078  nc3d(k)=ndcnst*1.e6/rho(k)
2079  END IF
2080 
2081 ! CALCULATE SIZE DISTRIBUTION PARAMETERS
2082 ! MAKE SURE NUMBER CONCENTRATIONS AREN'T NEGATIVE
2083 
2084  ni3d(k) = max(0.,ni3d(k))
2085  ns3d(k) = max(0.,ns3d(k))
2086  nc3d(k) = max(0.,nc3d(k))
2087  nr3d(k) = max(0.,nr3d(k))
2088  ng3d(k) = max(0.,ng3d(k))
2089 
2090 !......................................................................
2091 ! CLOUD ICE
2092 
2093  IF (qi3d(k).GE.qsmall) THEN
2094  lami(k) = (cons12* &
2095  ni3d(k)/qi3d(k))**(1./di)
2096  n0i(k) = ni3d(k)*lami(k)
2097 
2098 ! CHECK FOR SLOPE
2099 
2100 ! ADJUST VARS
2101 
2102  IF (lami(k).LT.lammini) THEN
2103 
2104  lami(k) = lammini
2105 
2106  n0i(k) = lami(k)**4*qi3d(k)/cons12
2107 
2108  ni3d(k) = n0i(k)/lami(k)
2109  ELSE IF (lami(k).GT.lammaxi) THEN
2110  lami(k) = lammaxi
2111  n0i(k) = lami(k)**4*qi3d(k)/cons12
2112 
2113  ni3d(k) = n0i(k)/lami(k)
2114  END IF
2115  END IF
2116 
2117 !......................................................................
2118 ! RAIN
2119 
2120  IF (qr3d(k).GE.qsmall) THEN
2121  lamr(k) = (pi*rhow*nr3d(k)/qr3d(k))**(1./3.)
2122  n0rr(k) = nr3d(k)*lamr(k)
2123 
2124 ! CHECK FOR SLOPE
2125 
2126 ! ADJUST VARS
2127 
2128  IF (lamr(k).LT.lamminr) THEN
2129 
2130  lamr(k) = lamminr
2131 
2132  n0rr(k) = lamr(k)**4*qr3d(k)/(pi*rhow)
2133 
2134  nr3d(k) = n0rr(k)/lamr(k)
2135  ELSE IF (lamr(k).GT.lammaxr) THEN
2136  lamr(k) = lammaxr
2137  n0rr(k) = lamr(k)**4*qr3d(k)/(pi*rhow)
2138 
2139  nr3d(k) = n0rr(k)/lamr(k)
2140  END IF
2141  END IF
2142 !......................................................................
2143 ! CLOUD DROPLETS
2144 
2145 ! MARTIN ET AL. (1994) FORMULA FOR PGAM
2146 
2147  IF (qc3d(k).GE.qsmall) THEN
2148 
2149  dum = pres(k)/(287.15*t3d(k))
2150  pgam(k)=0.0005714*(nc3d(k)/1.e6*dum)+0.2714
2151  pgam(k)=1./(pgam(k)**2)-1.
2152  pgam(k)=max(pgam(k),2.)
2153  pgam(k)=min(pgam(k),10.)
2154 
2155 ! CALCULATE LAMC
2156 
2157  lamc(k) = (cons26*nc3d(k)*gamma(pgam(k)+4.)/ &
2158  (qc3d(k)*gamma(pgam(k)+1.)))**(1./3.)
2159 
2160 ! LAMMIN, 60 MICRON DIAMETER
2161 ! LAMMAX, 1 MICRON
2162 
2163  lammin = (pgam(k)+1.)/60.e-6
2164  lammax = (pgam(k)+1.)/1.e-6
2165 
2166  IF (lamc(k).LT.lammin) THEN
2167  lamc(k) = lammin
2168 
2169  nc3d(k) = exp(3.*log(lamc(k))+log(qc3d(k))+ &
2170  log(gamma(pgam(k)+1.))-log(gamma(pgam(k)+4.)))/cons26
2171  ELSE IF (lamc(k).GT.lammax) THEN
2172  lamc(k) = lammax
2173  nc3d(k) = exp(3.*log(lamc(k))+log(qc3d(k))+ &
2174  log(gamma(pgam(k)+1.))-log(gamma(pgam(k)+4.)))/cons26
2175 
2176  END IF
2177 
2178 ! TO CALCULATE DROPLET FREEZING
2179 
2180  cdist1(k) = nc3d(k)/gamma(pgam(k)+1.)
2181 
2182  END IF
2183 
2184 !......................................................................
2185 ! SNOW
2186 
2187  IF (qni3d(k).GE.qsmall) THEN
2188  lams(k) = (cons1*ns3d(k)/qni3d(k))**(1./ds)
2189  n0s(k) = ns3d(k)*lams(k)
2190 
2191 ! CHECK FOR SLOPE
2192 
2193 ! ADJUST VARS
2194 
2195  IF (lams(k).LT.lammins) THEN
2196  lams(k) = lammins
2197  n0s(k) = lams(k)**4*qni3d(k)/cons1
2198 
2199  ns3d(k) = n0s(k)/lams(k)
2200 
2201  ELSE IF (lams(k).GT.lammaxs) THEN
2202 
2203  lams(k) = lammaxs
2204  n0s(k) = lams(k)**4*qni3d(k)/cons1
2205 
2206  ns3d(k) = n0s(k)/lams(k)
2207  END IF
2208  END IF
2209 
2210 !......................................................................
2211 ! GRAUPEL
2212 
2213  IF (qg3d(k).GE.qsmall) THEN
2214  lamg(k) = (cons2*ng3d(k)/qg3d(k))**(1./dg)
2215  n0g(k) = ng3d(k)*lamg(k)
2216 
2217 ! CHECK FOR SLOPE
2218 
2219 ! ADJUST VARS
2220 
2221  IF (lamg(k).LT.lamming) THEN
2222  lamg(k) = lamming
2223  n0g(k) = lamg(k)**4*qg3d(k)/cons2
2224 
2225  ng3d(k) = n0g(k)/lamg(k)
2226 
2227  ELSE IF (lamg(k).GT.lammaxg) THEN
2228 
2229  lamg(k) = lammaxg
2230  n0g(k) = lamg(k)**4*qg3d(k)/cons2
2231 
2232  ng3d(k) = n0g(k)/lamg(k)
2233  END IF
2234  END IF
2235 !.....................................................................
2236 ! ZERO OUT PROCESS RATES
2237 
2238  mnuccc(k) = 0.
2239  nnuccc(k) = 0.
2240  prc(k) = 0.
2241  nprc(k) = 0.
2242  nprc1(k) = 0.
2243  nsagg(k) = 0.
2244  psacws(k) = 0.
2245  npsacws(k) = 0.
2246  psacwi(k) = 0.
2247  npsacwi(k) = 0.
2248  pracs(k) = 0.
2249  npracs(k) = 0.
2250  nmults(k) = 0.
2251  qmults(k) = 0.
2252  nmultr(k) = 0.
2253  qmultr(k) = 0.
2254  nmultg(k) = 0.
2255  qmultg(k) = 0.
2256  nmultrg(k) = 0.
2257  qmultrg(k) = 0.
2258  mnuccr(k) = 0.
2259  nnuccr(k) = 0.
2260  pra(k) = 0.
2261  npra(k) = 0.
2262  nragg(k) = 0.
2263  prci(k) = 0.
2264  nprci(k) = 0.
2265  prai(k) = 0.
2266  nprai(k) = 0.
2267  nnuccd(k) = 0.
2268  mnuccd(k) = 0.
2269  pcc(k) = 0.
2270  pre(k) = 0.
2271  prd(k) = 0.
2272  prds(k) = 0.
2273  eprd(k) = 0.
2274  eprds(k) = 0.
2275  nsubc(k) = 0.
2276  nsubi(k) = 0.
2277  nsubs(k) = 0.
2278  nsubr(k) = 0.
2279  piacr(k) = 0.
2280  niacr(k) = 0.
2281  praci(k) = 0.
2282  piacrs(k) = 0.
2283  niacrs(k) = 0.
2284  pracis(k) = 0.
2285 ! HM: ADD GRAUPEL PROCESSES
2286  pracg(k) = 0.
2287  psacr(k) = 0.
2288  psacwg(k) = 0.
2289  pgsacw(k) = 0.
2290  pgracs(k) = 0.
2291  prdg(k) = 0.
2292  eprdg(k) = 0.
2293  npracg(k) = 0.
2294  npsacwg(k) = 0.
2295  nscng(k) = 0.
2296  ngracs(k) = 0.
2297  nsubg(k) = 0.
2298 
2299 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
2300 ! CALCULATION OF MICROPHYSICAL PROCESS RATES
2301 ! ACCRETION/AUTOCONVERSION/FREEZING/MELTING/COAG.
2302 !.......................................................................
2303 ! FREEZING OF CLOUD DROPLETS
2304 ! ONLY ALLOWED BELOW -4 C
2305  IF (qc3d(k).GE.qsmall .AND. t3d(k).LT.269.15) THEN
2306 
2307 ! NUMBER OF CONTACT NUCLEI (M^-3) FROM MEYERS ET AL., 1992
2308 ! FACTOR OF 1000 IS TO CONVERT FROM L^-1 TO M^-3
2309 
2310 ! MEYERS CURVE
2311 
2312  nacnt = exp(-2.80+0.262*(273.15-t3d(k)))*1000.
2313 
2314 ! COOPER CURVE
2315 ! NACNT = 5.*EXP(0.304*(273.15-T3D(K)))
2316 
2317 ! FLECTHER
2318 ! NACNT = 0.01*EXP(0.6*(273.15-T3D(K)))
2319 
2320 ! CONTACT FREEZING
2321 
2322 ! MEAN FREE PATH
2323 
2324  dum = 7.37*t3d(k)/(288.*10.*pres(k))/100.
2325 
2326 ! EFFECTIVE DIFFUSIVITY OF CONTACT NUCLEI
2327 ! BASED ON BROWNIAN DIFFUSION
2328 
2329  dap(k) = cons37*t3d(k)*(1.+dum/rin)/mu(k)
2330 
2331  mnuccc(k) = cons38*dap(k)*nacnt*exp(log(cdist1(k))+ &
2332  log(gamma(pgam(k)+5.))-4.*log(lamc(k)))
2333  nnuccc(k) = 2.*pi*dap(k)*nacnt*cdist1(k)* &
2334  gamma(pgam(k)+2.)/ &
2335  lamc(k)
2336 
2337 ! IMMERSION FREEZING (BIGG 1953)
2338 
2339 ! MNUCCC(K) = MNUCCC(K)+CONS39* &
2340 ! EXP(LOG(CDIST1(K))+LOG(GAMMA(7.+PGAM(K)))-6.*LOG(LAMC(K)))* &
2341 ! EXP(AIMM*(273.15-T3D(K)))
2342 
2343 ! NNUCCC(K) = NNUCCC(K)+ &
2344 ! CONS40*EXP(LOG(CDIST1(K))+LOG(GAMMA(PGAM(K)+4.))-3.*LOG(LAMC(K))) &
2345 ! *EXP(AIMM*(273.15-T3D(K)))
2346 
2347 ! hm 7/15/13 fix for consistency w/ original formula
2348  mnuccc(k) = mnuccc(k)+cons39* &
2349  exp(log(cdist1(k))+log(gamma(7.+pgam(k)))-6.*log(lamc(k)))* &
2350  (exp(aimm*(273.15-t3d(k)))-1.)
2351 
2352  nnuccc(k) = nnuccc(k)+ &
2353  cons40*exp(log(cdist1(k))+log(gamma(pgam(k)+4.))-3.*log(lamc(k))) &
2354  *(exp(aimm*(273.15-t3d(k)))-1.)
2355 
2356 ! PUT IN A CATCH HERE TO PREVENT DIVERGENCE BETWEEN NUMBER CONC. AND
2357 ! MIXING RATIO, SINCE STRICT CONSERVATION NOT CHECKED FOR NUMBER CONC
2358 
2359  nnuccc(k) = min(nnuccc(k),nc3d(k)/dt)
2360 
2361  END IF
2362 
2363 !.................................................................
2364 !.......................................................................
2365 ! AUTOCONVERSION OF CLOUD LIQUID WATER TO RAIN
2366 ! FORMULA FROM BEHENG (1994)
2367 ! USING NUMERICAL SIMULATION OF STOCHASTIC COLLECTION EQUATION
2368 ! AND INITIAL CLOUD DROPLET SIZE DISTRIBUTION SPECIFIED
2369 ! AS A GAMMA DISTRIBUTION
2370 
2371 ! USE MINIMUM VALUE OF 1.E-6 TO PREVENT FLOATING POINT ERROR
2372 
2373  IF (qc3d(k).GE.1.e-6) THEN
2374 
2375 ! HM ADD 12/13/06, REPLACE WITH NEWER FORMULA
2376 ! FROM KHAIROUTDINOV AND KOGAN 2000, MWR
2377 
2378  prc(k)=1350.*qc3d(k)**2.47* &
2379  (nc3d(k)/1.e6*rho(k))**(-1.79)
2380 
2381 ! note: nprc1 is change in Nr,
2382 ! nprc is change in Nc
2383 
2384  nprc1(k) = prc(k)/cons29
2385  nprc(k) = prc(k)/(qc3d(k)/nc3d(k))
2386 
2387 ! hm bug fix 3/20/12
2388  nprc(k) = min(nprc(k),nc3d(k)/dt)
2389  nprc1(k) = min(nprc1(k),nprc(k))
2390 
2391  END IF
2392 !.......................................................................
2393 ! SELF-COLLECTION OF DROPLET NOT INCLUDED IN KK2000 SCHEME
2394 
2395 ! SNOW AGGREGATION FROM PASSARELLI, 1978, USED BY REISNER, 1998
2396 ! THIS IS HARD-WIRED FOR BS = 0.4 FOR NOW
2397 
2398  IF (qni3d(k).GE.1.e-8) THEN
2399  nsagg(k) = cons15*asn(k)*rho(k)** &
2400  ((2.+bs)/3.)*qni3d(k)**((2.+bs)/3.)* &
2401  (ns3d(k)*rho(k))**((4.-bs)/3.)/ &
2402  (rho(k))
2403  END IF
2404 
2405 !.......................................................................
2406 ! ACCRETION OF CLOUD DROPLETS ONTO SNOW/GRAUPEL
2407 ! HERE USE CONTINUOUS COLLECTION EQUATION WITH
2408 ! SIMPLE GRAVITATIONAL COLLECTION KERNEL IGNORING
2409 
2410 ! SNOW
2411 
2412  IF (qni3d(k).GE.1.e-8 .AND. qc3d(k).GE.qsmall) THEN
2413 
2414  psacws(k) = cons13*asn(k)*qc3d(k)*rho(k)* &
2415  n0s(k)/ &
2416  lams(k)**(bs+3.)
2417  npsacws(k) = cons13*asn(k)*nc3d(k)*rho(k)* &
2418  n0s(k)/ &
2419  lams(k)**(bs+3.)
2420 
2421  END IF
2422 
2423 !............................................................................
2424 ! COLLECTION OF CLOUD WATER BY GRAUPEL
2425 
2426  IF (qg3d(k).GE.1.e-8 .AND. qc3d(k).GE.qsmall) THEN
2427 
2428  psacwg(k) = cons14*agn(k)*qc3d(k)*rho(k)* &
2429  n0g(k)/ &
2430  lamg(k)**(bg+3.)
2431  npsacwg(k) = cons14*agn(k)*nc3d(k)*rho(k)* &
2432  n0g(k)/ &
2433  lamg(k)**(bg+3.)
2434  END IF
2435 !.......................................................................
2436 ! HM, ADD 12/13/06
2437 ! CLOUD ICE COLLECTING DROPLETS, ASSUME THAT CLOUD ICE MEAN DIAM > 100 MICRON
2438 ! BEFORE RIMING CAN OCCUR
2439 ! ASSUME THAT RIME COLLECTED ON CLOUD ICE DOES NOT LEAD
2440 ! TO HALLET-MOSSOP SPLINTERING
2441 
2442  IF (qi3d(k).GE.1.e-8 .AND. qc3d(k).GE.qsmall) THEN
2443 
2444 ! PUT IN SIZE DEPENDENT COLLECTION EFFICIENCY BASED ON STOKES LAW
2445 ! FROM THOMPSON ET AL. 2004, MWR
2446 
2447  IF (1./lami(k).GE.100.e-6) THEN
2448 
2449  psacwi(k) = cons16*ain(k)*qc3d(k)*rho(k)* &
2450  n0i(k)/ &
2451  lami(k)**(bi+3.)
2452  npsacwi(k) = cons16*ain(k)*nc3d(k)*rho(k)* &
2453  n0i(k)/ &
2454  lami(k)**(bi+3.)
2455  END IF
2456  END IF
2457 
2458 !.......................................................................
2459 ! ACCRETION OF RAIN WATER BY SNOW
2460 ! FORMULA FROM IKAWA AND SAITO, 1991, USED BY REISNER ET AL, 1998
2461 
2462  IF (qr3d(k).GE.1.e-8.AND.qni3d(k).GE.1.e-8) THEN
2463 
2464  ums = asn(k)*cons3/(lams(k)**bs)
2465  umr = arn(k)*cons4/(lamr(k)**br)
2466  uns = asn(k)*cons5/lams(k)**bs
2467  unr = arn(k)*cons6/lamr(k)**br
2468 
2469 ! SET REASLISTIC LIMITS ON FALLSPEEDS
2470 
2471 ! bug fix, 10/08/09
2472  dum=(rhosu/rho(k))**0.54
2473  ums=min(ums,1.2*dum)
2474  uns=min(uns,1.2*dum)
2475  umr=min(umr,9.1*dum)
2476  unr=min(unr,9.1*dum)
2477 
2478  pracs(k) = cons41*(((1.2*umr-0.95*ums)**2+ &
2479  0.08*ums*umr)**0.5*rho(k)* &
2480  n0rr(k)*n0s(k)/lamr(k)**3* &
2481  (5./(lamr(k)**3*lams(k))+ &
2482  2./(lamr(k)**2*lams(k)**2)+ &
2483  0.5/(lamr(k)*lams(k)**3)))
2484 
2485  npracs(k) = cons32*rho(k)*(1.7*(unr-uns)**2+ &
2486  0.3*unr*uns)**0.5*n0rr(k)*n0s(k)* &
2487  (1./(lamr(k)**3*lams(k))+ &
2488  1./(lamr(k)**2*lams(k)**2)+ &
2489  1./(lamr(k)*lams(k)**3))
2490 
2491 ! MAKE SURE PRACS DOESN'T EXCEED TOTAL RAIN MIXING RATIO
2492 ! AS THIS MAY OTHERWISE RESULT IN TOO MUCH TRANSFER OF WATER DURING
2493 ! RIME-SPLINTERING
2494 
2495  pracs(k) = min(pracs(k),qr3d(k)/dt)
2496 
2497 ! COLLECTION OF SNOW BY RAIN - NEEDED FOR GRAUPEL CONVERSION CALCULATIONS
2498 ! ONLY CALCULATE IF SNOW AND RAIN MIXING RATIOS EXCEED 0.1 G/KG
2499 
2500 ! HM MODIFY FOR WRFV3.1
2501 ! IF (IHAIL.EQ.0) THEN
2502  IF (qni3d(k).GE.0.1e-3.AND.qr3d(k).GE.0.1e-3) THEN
2503  psacr(k) = cons31*(((1.2*umr-0.95*ums)**2+ &
2504  0.08*ums*umr)**0.5*rho(k)* &
2505  n0rr(k)*n0s(k)/lams(k)**3* &
2506  (5./(lams(k)**3*lamr(k))+ &
2507  2./(lams(k)**2*lamr(k)**2)+ &
2508  0.5/(lams(k)*lamr(k)**3)))
2509  END IF
2510 ! END IF
2511 
2512  END IF
2513 
2514 !.......................................................................
2515 
2516 ! COLLECTION OF RAINWATER BY GRAUPEL, FROM IKAWA AND SAITO 1990,
2517 ! USED BY REISNER ET AL 1998
2518  IF (qr3d(k).GE.1.e-8.AND.qg3d(k).GE.1.e-8) THEN
2519 
2520  umg = agn(k)*cons7/(lamg(k)**bg)
2521  umr = arn(k)*cons4/(lamr(k)**br)
2522  ung = agn(k)*cons8/lamg(k)**bg
2523  unr = arn(k)*cons6/lamr(k)**br
2524 
2525 ! SET REASLISTIC LIMITS ON FALLSPEEDS
2526 ! bug fix, 10/08/09
2527  dum=(rhosu/rho(k))**0.54
2528  umg=min(umg,20.*dum)
2529  ung=min(ung,20.*dum)
2530  umr=min(umr,9.1*dum)
2531  unr=min(unr,9.1*dum)
2532 
2533  pracg(k) = cons41*(((1.2*umr-0.95*umg)**2+ &
2534  0.08*umg*umr)**0.5*rho(k)* &
2535  n0rr(k)*n0g(k)/lamr(k)**3* &
2536  (5./(lamr(k)**3*lamg(k))+ &
2537  2./(lamr(k)**2*lamg(k)**2)+ &
2538  0.5/(lamr(k)*lamg(k)**3)))
2539 
2540  npracg(k) = cons32*rho(k)*(1.7*(unr-ung)**2+ &
2541  0.3*unr*ung)**0.5*n0rr(k)*n0g(k)* &
2542  (1./(lamr(k)**3*lamg(k))+ &
2543  1./(lamr(k)**2*lamg(k)**2)+ &
2544  1./(lamr(k)*lamg(k)**3))
2545 
2546 ! MAKE SURE PRACG DOESN'T EXCEED TOTAL RAIN MIXING RATIO
2547 ! AS THIS MAY OTHERWISE RESULT IN TOO MUCH TRANSFER OF WATER DURING
2548 ! RIME-SPLINTERING
2549 
2550  pracg(k) = min(pracg(k),qr3d(k)/dt)
2551 
2552  END IF
2553 
2554 !.......................................................................
2555 ! RIME-SPLINTERING - SNOW
2556 ! HALLET-MOSSOP (1974)
2557 ! NUMBER OF SPLINTERS FORMED IS BASED ON MASS OF RIMED WATER
2558 
2559 ! DUM1 = MASS OF INDIVIDUAL SPLINTERS
2560 
2561 ! HM ADD THRESHOLD SNOW AND DROPLET MIXING RATIO FOR RIME-SPLINTERING
2562 ! TO LIMIT RIME-SPLINTERING IN STRATIFORM CLOUDS
2563 ! THESE THRESHOLDS CORRESPOND WITH GRAUPEL THRESHOLDS IN RH 1984
2564 
2565 !v1.4
2566  IF (qni3d(k).GE.0.1e-3) THEN
2567  IF (qc3d(k).GE.0.5e-3.OR.qr3d(k).GE.0.1e-3) THEN
2568  IF (psacws(k).GT.0..OR.pracs(k).GT.0.) THEN
2569  IF (t3d(k).LT.270.16 .AND. t3d(k).GT.265.16) THEN
2570 
2571  IF (t3d(k).GT.270.16) THEN
2572  fmult = 0.
2573  ELSE IF (t3d(k).LE.270.16.AND.t3d(k).GT.268.16) THEN
2574  fmult = (270.16-t3d(k))/2.
2575  ELSE IF (t3d(k).GE.265.16.AND.t3d(k).LE.268.16) THEN
2576  fmult = (t3d(k)-265.16)/3.
2577  ELSE IF (t3d(k).LT.265.16) THEN
2578  fmult = 0.
2579  END IF
2580 
2581 ! 1000 IS TO CONVERT FROM KG TO G
2582 
2583 ! SPLINTERING FROM DROPLETS ACCRETED ONTO SNOW
2584 
2585  IF (psacws(k).GT.0.) THEN
2586  nmults(k) = 35.e4*psacws(k)*fmult*1000.
2587  qmults(k) = nmults(k)*mmult
2588 
2589 ! CONSTRAIN SO THAT TRANSFER OF MASS FROM SNOW TO ICE CANNOT BE MORE MASS
2590 ! THAN WAS RIMED ONTO SNOW
2591 
2592  qmults(k) = min(qmults(k),psacws(k))
2593  psacws(k) = psacws(k)-qmults(k)
2594 
2595  END IF
2596 
2597 ! RIMING AND SPLINTERING FROM ACCRETED RAINDROPS
2598 
2599  IF (pracs(k).GT.0.) THEN
2600  nmultr(k) = 35.e4*pracs(k)*fmult*1000.
2601  qmultr(k) = nmultr(k)*mmult
2602 
2603 ! CONSTRAIN SO THAT TRANSFER OF MASS FROM SNOW TO ICE CANNOT BE MORE MASS
2604 ! THAN WAS RIMED ONTO SNOW
2605 
2606  qmultr(k) = min(qmultr(k),pracs(k))
2607 
2608  pracs(k) = pracs(k)-qmultr(k)
2609 
2610  END IF
2611 
2612  END IF
2613  END IF
2614  END IF
2615  END IF
2616 
2617 !.......................................................................
2618 ! RIME-SPLINTERING - GRAUPEL
2619 ! HALLET-MOSSOP (1974)
2620 ! NUMBER OF SPLINTERS FORMED IS BASED ON MASS OF RIMED WATER
2621 
2622 ! DUM1 = MASS OF INDIVIDUAL SPLINTERS
2623 
2624 ! HM ADD THRESHOLD SNOW MIXING RATIO FOR RIME-SPLINTERING
2625 ! TO LIMIT RIME-SPLINTERING IN STRATIFORM CLOUDS
2626 
2627 ! IF (IHAIL.EQ.0) THEN
2628 ! v1.4
2629  IF (qg3d(k).GE.0.1e-3) THEN
2630  IF (qc3d(k).GE.0.5e-3.OR.qr3d(k).GE.0.1e-3) THEN
2631  IF (psacwg(k).GT.0..OR.pracg(k).GT.0.) THEN
2632  IF (t3d(k).LT.270.16 .AND. t3d(k).GT.265.16) THEN
2633 
2634  IF (t3d(k).GT.270.16) THEN
2635  fmult = 0.
2636  ELSE IF (t3d(k).LE.270.16.AND.t3d(k).GT.268.16) THEN
2637  fmult = (270.16-t3d(k))/2.
2638  ELSE IF (t3d(k).GE.265.16.AND.t3d(k).LE.268.16) THEN
2639  fmult = (t3d(k)-265.16)/3.
2640  ELSE IF (t3d(k).LT.265.16) THEN
2641  fmult = 0.
2642  END IF
2643 
2644 ! 1000 IS TO CONVERT FROM KG TO G
2645 
2646 ! SPLINTERING FROM DROPLETS ACCRETED ONTO GRAUPEL
2647 
2648  IF (psacwg(k).GT.0.) THEN
2649  nmultg(k) = 35.e4*psacwg(k)*fmult*1000.
2650  qmultg(k) = nmultg(k)*mmult
2651 
2652 ! CONSTRAIN SO THAT TRANSFER OF MASS FROM GRAUPEL TO ICE CANNOT BE MORE MASS
2653 ! THAN WAS RIMED ONTO GRAUPEL
2654 
2655  qmultg(k) = min(qmultg(k),psacwg(k))
2656  psacwg(k) = psacwg(k)-qmultg(k)
2657 
2658  END IF
2659 
2660 ! RIMING AND SPLINTERING FROM ACCRETED RAINDROPS
2661 
2662  IF (pracg(k).GT.0.) THEN
2663  nmultrg(k) = 35.e4*pracg(k)*fmult*1000.
2664  qmultrg(k) = nmultrg(k)*mmult
2665 
2666 ! CONSTRAIN SO THAT TRANSFER OF MASS FROM GRAUPEL TO ICE CANNOT BE MORE MASS
2667 ! THAN WAS RIMED ONTO GRAUPEL
2668 
2669  qmultrg(k) = min(qmultrg(k),pracg(k))
2670  pracg(k) = pracg(k)-qmultrg(k)
2671 
2672  END IF
2673  END IF
2674  END IF
2675  END IF
2676  END IF
2677 ! END IF
2678 
2679 !........................................................................
2680 ! CONVERSION OF RIMED CLOUD WATER ONTO SNOW TO GRAUPEL/HAIL
2681 
2682 ! IF (IHAIL.EQ.0) THEN
2683  IF (psacws(k).GT.0.) THEN
2684 ! ONLY ALLOW CONVERSION IF QNI > 0.1 AND QC > 0.5 G/KG FOLLOWING RUTLEDGE AND HOBBS (1984)
2685  IF (qni3d(k).GE.0.1e-3.AND.qc3d(k).GE.0.5e-3) THEN
2686 
2687 ! PORTION OF RIMING CONVERTED TO GRAUPEL (REISNER ET AL. 1998, ORIGINALLY IS1991)
2688  pgsacw(k) = min(psacws(k),cons17*dt*n0s(k)*qc3d(k)*qc3d(k)* &
2689  asn(k)*asn(k)/ &
2690  (rho(k)*lams(k)**(2.*bs+2.)))
2691 
2692 ! MIX RAT CONVERTED INTO GRAUPEL AS EMBRYO (REISNER ET AL. 1998, ORIG M1990)
2693  dum = max(rhosn/(rhog-rhosn)*pgsacw(k),0.)
2694 
2695 ! NUMBER CONCENTRAITON OF EMBRYO GRAUPEL FROM RIMING OF SNOW
2696  nscng(k) = dum/mg0*rho(k)
2697 ! LIMIT MAX NUMBER CONVERTED TO SNOW NUMBER
2698  nscng(k) = min(nscng(k),ns3d(k)/dt)
2699 
2700 ! PORTION OF RIMING LEFT FOR SNOW
2701  psacws(k) = psacws(k) - pgsacw(k)
2702  END IF
2703  END IF
2704 
2705 ! CONVERSION OF RIMED RAINWATER ONTO SNOW CONVERTED TO GRAUPEL
2706 
2707  IF (pracs(k).GT.0.) THEN
2708 ! ONLY ALLOW CONVERSION IF QNI > 0.1 AND QR > 0.1 G/KG FOLLOWING RUTLEDGE AND HOBBS (1984)
2709  IF (qni3d(k).GE.0.1e-3.AND.qr3d(k).GE.0.1e-3) THEN
2710 ! PORTION OF COLLECTED RAINWATER CONVERTED TO GRAUPEL (REISNER ET AL. 1998)
2711  dum = cons18*(4./lams(k))**3*(4./lams(k))**3 &
2712  /(cons18*(4./lams(k))**3*(4./lams(k))**3+ &
2713  cons19*(4./lamr(k))**3*(4./lamr(k))**3)
2714  dum=min(dum,1.)
2715  dum=max(dum,0.)
2716  pgracs(k) = (1.-dum)*pracs(k)
2717  ngracs(k) = (1.-dum)*npracs(k)
2718 ! LIMIT MAX NUMBER CONVERTED TO MIN OF EITHER RAIN OR SNOW NUMBER CONCENTRATION
2719  ngracs(k) = min(ngracs(k),nr3d(k)/dt)
2720  ngracs(k) = min(ngracs(k),ns3d(k)/dt)
2721 
2722 ! AMOUNT LEFT FOR SNOW PRODUCTION
2723  pracs(k) = pracs(k) - pgracs(k)
2724  npracs(k) = npracs(k) - ngracs(k)
2725 ! CONVERSION TO GRAUPEL DUE TO COLLECTION OF SNOW BY RAIN
2726  psacr(k)=psacr(k)*(1.-dum)
2727  END IF
2728  END IF
2729 ! END IF
2730 
2731 !.......................................................................
2732 ! FREEZING OF RAIN DROPS
2733 ! FREEZING ALLOWED BELOW -4 C
2734 
2735  IF (t3d(k).LT.269.15.AND.qr3d(k).GE.qsmall) THEN
2736 
2737 ! IMMERSION FREEZING (BIGG 1953)
2738 ! MNUCCR(K) = CONS20*NR3D(K)*EXP(AIMM*(273.15-T3D(K)))/LAMR(K)**3 &
2739 ! /LAMR(K)**3
2740 
2741 ! NNUCCR(K) = PI*NR3D(K)*BIMM*EXP(AIMM*(273.15-T3D(K)))/LAMR(K)**3
2742 
2743 ! hm fix 7/15/13 for consistency w/ original formula
2744  mnuccr(k) = cons20*nr3d(k)*(exp(aimm*(273.15-t3d(k)))-1.)/lamr(k)**3 &
2745  /lamr(k)**3
2746 
2747  nnuccr(k) = pi*nr3d(k)*bimm*(exp(aimm*(273.15-t3d(k)))-1.)/lamr(k)**3
2748 
2749 ! PREVENT DIVERGENCE BETWEEN MIXING RATIO AND NUMBER CONC
2750  nnuccr(k) = min(nnuccr(k),nr3d(k)/dt)
2751 
2752  END IF
2753 
2754 !.......................................................................
2755 ! ACCRETION OF CLOUD LIQUID WATER BY RAIN
2756 ! CONTINUOUS COLLECTION EQUATION WITH
2757 ! GRAVITATIONAL COLLECTION KERNEL, DROPLET FALL SPEED NEGLECTED
2758 
2759  IF (qr3d(k).GE.1.e-8 .AND. qc3d(k).GE.1.e-8) THEN
2760 
2761 ! 12/13/06 HM ADD, REPLACE WITH NEWER FORMULA FROM
2762 ! KHAIROUTDINOV AND KOGAN 2000, MWR
2763 
2764  dum=(qc3d(k)*qr3d(k))
2765  pra(k) = 67.*(dum)**1.15
2766  npra(k) = pra(k)/(qc3d(k)/nc3d(k))
2767 
2768  END IF
2769 !.......................................................................
2770 ! SELF-COLLECTION OF RAIN DROPS
2771 ! FROM BEHENG(1994)
2772 ! FROM NUMERICAL SIMULATION OF THE STOCHASTIC COLLECTION EQUATION
2773 ! AS DESCRINED ABOVE FOR AUTOCONVERSION
2774 
2775  IF (qr3d(k).GE.1.e-8) THEN
2776 ! include breakup add 10/09/09
2777  dum1=300.e-6
2778  if (1./lamr(k).lt.dum1) then
2779  dum=1.
2780  else if (1./lamr(k).ge.dum1) then
2781  dum=2.-exp(2300.*(1./lamr(k)-dum1))
2782  end if
2783 ! NRAGG(K) = -8.*NR3D(K)*QR3D(K)*RHO(K)
2784  nragg(k) = -5.78*dum*nr3d(k)*qr3d(k)*rho(k)
2785  END IF
2786 
2787 !.......................................................................
2788 ! AUTOCONVERSION OF CLOUD ICE TO SNOW
2789 ! FOLLOWING HARRINGTON ET AL. (1995) WITH MODIFICATION
2790 ! HERE IT IS ASSUMED THAT AUTOCONVERSION CAN ONLY OCCUR WHEN THE
2791 ! ICE IS GROWING, I.E. IN CONDITIONS OF ICE SUPERSATURATION
2792 
2793  IF (qi3d(k).GE.1.e-8 .AND.qvqvsi(k).GE.1.) THEN
2794 
2795 ! COFFI = 2./LAMI(K)
2796 ! IF (COFFI.GE.DCS) THEN
2797  nprci(k) = cons21*(qv3d(k)-qvi(k))*rho(k) &
2798  *n0i(k)*exp(-lami(k)*dcs)*dv(k)/abi(k)
2799  prci(k) = cons22*nprci(k)
2800  nprci(k) = min(nprci(k),ni3d(k)/dt)
2801 
2802 ! END IF
2803  END IF
2804 
2805 !.......................................................................
2806 ! ACCRETION OF CLOUD ICE BY SNOW
2807 ! FOR THIS CALCULATION, IT IS ASSUMED THAT THE VS >> VI
2808 ! AND DS >> DI FOR CONTINUOUS COLLECTION
2809 
2810  IF (qni3d(k).GE.1.e-8 .AND. qi3d(k).GE.qsmall) THEN
2811  prai(k) = cons23*asn(k)*qi3d(k)*rho(k)*n0s(k)/ &
2812  lams(k)**(bs+3.)
2813  nprai(k) = cons23*asn(k)*ni3d(k)* &
2814  rho(k)*n0s(k)/ &
2815  lams(k)**(bs+3.)
2816  nprai(k)=min(nprai(k),ni3d(k)/dt)
2817  END IF
2818 
2819 !.......................................................................
2820 ! HM, ADD 12/13/06, COLLISION OF RAIN AND ICE TO PRODUCE SNOW OR GRAUPEL
2821 ! FOLLOWS REISNER ET AL. 1998
2822 ! ASSUMED FALLSPEED AND SIZE OF ICE CRYSTAL << THAN FOR RAIN
2823 
2824  IF (qr3d(k).GE.1.e-8.AND.qi3d(k).GE.1.e-8.AND.t3d(k).LE.273.15) THEN
2825 
2826 ! ALLOW GRAUPEL FORMATION FROM RAIN-ICE COLLISIONS ONLY IF RAIN MIXING RATIO > 0.1 G/KG,
2827 ! OTHERWISE ADD TO SNOW
2828 
2829  IF (qr3d(k).GE.0.1e-3) THEN
2830  niacr(k)=cons24*ni3d(k)*n0rr(k)*arn(k) &
2831  /lamr(k)**(br+3.)*rho(k)
2832  piacr(k)=cons25*ni3d(k)*n0rr(k)*arn(k) &
2833  /lamr(k)**(br+3.)/lamr(k)**3*rho(k)
2834  praci(k)=cons24*qi3d(k)*n0rr(k)*arn(k)/ &
2835  lamr(k)**(br+3.)*rho(k)
2836  niacr(k)=min(niacr(k),nr3d(k)/dt)
2837  niacr(k)=min(niacr(k),ni3d(k)/dt)
2838  ELSE
2839  niacrs(k)=cons24*ni3d(k)*n0rr(k)*arn(k) &
2840  /lamr(k)**(br+3.)*rho(k)
2841  piacrs(k)=cons25*ni3d(k)*n0rr(k)*arn(k) &
2842  /lamr(k)**(br+3.)/lamr(k)**3*rho(k)
2843  pracis(k)=cons24*qi3d(k)*n0rr(k)*arn(k)/ &
2844  lamr(k)**(br+3.)*rho(k)
2845  niacrs(k)=min(niacrs(k),nr3d(k)/dt)
2846  niacrs(k)=min(niacrs(k),ni3d(k)/dt)
2847  END IF
2848  END IF
2849 
2850 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
2851 ! NUCLEATION OF CLOUD ICE FROM HOMOGENEOUS AND HETEROGENEOUS FREEZING ON AEROSOL
2852 
2853  IF (inuc.EQ.0) THEN
2854 
2855 ! add threshold according to Greg Thomspon
2856 
2857  if ((qvqvs(k).GE.0.999.and.t3d(k).le.265.15).or. &
2858  qvqvsi(k).ge.1.08) then
2859 
2860 ! hm, modify dec. 5, 2006, replace with cooper curve
2861  kc2 = 0.005*exp(0.304*(273.15-t3d(k)))*1000. ! convert from L-1 to m-3
2862 ! limit to 500 L-1
2863  kc2 = min(kc2,500.e3)
2864  kc2=max(kc2/rho(k),0.) ! convert to kg-1
2865 
2866  IF (kc2.GT.ni3d(k)+ns3d(k)+ng3d(k)) THEN
2867  nnuccd(k) = (kc2-ni3d(k)-ns3d(k)-ng3d(k))/dt
2868  mnuccd(k) = nnuccd(k)*mi0
2869  END IF
2870 
2871  END IF
2872 
2873  ELSE IF (inuc.EQ.1) THEN
2874 
2875  IF (t3d(k).LT.273.15.AND.qvqvsi(k).GT.1.) THEN
2876 
2877  kc2 = 0.16*1000./rho(k) ! CONVERT FROM L-1 TO KG-1
2878  IF (kc2.GT.ni3d(k)+ns3d(k)+ng3d(k)) THEN
2879  nnuccd(k) = (kc2-ni3d(k)-ns3d(k)-ng3d(k))/dt
2880  mnuccd(k) = nnuccd(k)*mi0
2881  END IF
2882  END IF
2883 
2884  END IF
2885 
2886 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
2887 
2888  101 CONTINUE
2889 
2890 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
2891 ! CALCULATE EVAP/SUB/DEP TERMS FOR QI,QNI,QR
2892 
2893 ! NO VENTILATION FOR CLOUD ICE
2894 
2895  IF (qi3d(k).GE.qsmall) THEN
2896 
2897  epsi = 2.*pi*n0i(k)*rho(k)*dv(k)/(lami(k)*lami(k))
2898 
2899  ELSE
2900  epsi = 0.
2901  END IF
2902 
2903  IF (qni3d(k).GE.qsmall) THEN
2904  epss = 2.*pi*n0s(k)*rho(k)*dv(k)* &
2905  (f1s/(lams(k)*lams(k))+ &
2906  f2s*(asn(k)*rho(k)/mu(k))**0.5* &
2907  sc(k)**(1./3.)*cons10/ &
2908  (lams(k)**cons35))
2909  ELSE
2910  epss = 0.
2911  END IF
2912 
2913  IF (qg3d(k).GE.qsmall) THEN
2914  epsg = 2.*pi*n0g(k)*rho(k)*dv(k)* &
2915  (f1s/(lamg(k)*lamg(k))+ &
2916  f2s*(agn(k)*rho(k)/mu(k))**0.5* &
2917  sc(k)**(1./3.)*cons11/ &
2918  (lamg(k)**cons36))
2919 
2920 
2921  ELSE
2922  epsg = 0.
2923  END IF
2924 
2925  IF (qr3d(k).GE.qsmall) THEN
2926  epsr = 2.*pi*n0rr(k)*rho(k)*dv(k)* &
2927  (f1r/(lamr(k)*lamr(k))+ &
2928  f2r*(arn(k)*rho(k)/mu(k))**0.5* &
2929  sc(k)**(1./3.)*cons9/ &
2930  (lamr(k)**cons34))
2931  ELSE
2932  epsr = 0.
2933  END IF
2934 
2935 ! ONLY INCLUDE REGION OF ICE SIZE DIST < DCS
2936 ! DUM IS FRACTION OF D*N(D) < DCS
2937 
2938 ! LOGIC BELOW FOLLOWS THAT OF HARRINGTON ET AL. 1995 (JAS)
2939  IF (qi3d(k).GE.qsmall) THEN
2940  dum=(1.-exp(-lami(k)*dcs)*(1.+lami(k)*dcs))
2941  prd(k) = epsi*(qv3d(k)-qvi(k))/abi(k)*dum
2942  ELSE
2943  dum=0.
2944  END IF
2945 ! ADD DEPOSITION IN TAIL OF ICE SIZE DIST TO SNOW IF SNOW IS PRESENT
2946  IF (qni3d(k).GE.qsmall) THEN
2947  prds(k) = epss*(qv3d(k)-qvi(k))/abi(k)+ &
2948  epsi*(qv3d(k)-qvi(k))/abi(k)*(1.-dum)
2949 ! OTHERWISE ADD TO CLOUD ICE
2950  ELSE
2951  prd(k) = prd(k)+epsi*(qv3d(k)-qvi(k))/abi(k)*(1.-dum)
2952  END IF
2953 ! VAPOR DPEOSITION ON GRAUPEL
2954  prdg(k) = epsg*(qv3d(k)-qvi(k))/abi(k)
2955 
2956 ! NO CONDENSATION ONTO RAIN, ONLY EVAP
2957 
2958  IF (qv3d(k).LT.qvs(k)) THEN
2959  pre(k) = epsr*(qv3d(k)-qvs(k))/ab(k)
2960  pre(k) = min(pre(k),0.)
2961  ELSE
2962  pre(k) = 0.
2963  END IF
2964 
2965 ! MAKE SURE NOT PUSHED INTO ICE SUPERSAT/SUBSAT
2966 ! FORMULA FROM REISNER 2 SCHEME
2967 
2968  dum = (qv3d(k)-qvi(k))/dt
2969 
2970  fudgef = 0.9999
2971  sum_dep = prd(k)+prds(k)+mnuccd(k)+prdg(k)
2972 
2973  IF( (dum.GT.0. .AND. sum_dep.GT.dum*fudgef) .OR. &
2974  (dum.LT.0. .AND. sum_dep.LT.dum*fudgef) ) THEN
2975  mnuccd(k) = fudgef*mnuccd(k)*dum/sum_dep
2976  prd(k) = fudgef*prd(k)*dum/sum_dep
2977  prds(k) = fudgef*prds(k)*dum/sum_dep
2978  prdg(k) = fudgef*prdg(k)*dum/sum_dep
2979  ENDIF
2980 
2981 ! IF CLOUD ICE/SNOW/GRAUPEL VAP DEPOSITION IS NEG, THEN ASSIGN TO SUBLIMATION PROCESSES
2982 
2983  IF (prd(k).LT.0.) THEN
2984  eprd(k)=prd(k)
2985  prd(k)=0.
2986  END IF
2987  IF (prds(k).LT.0.) THEN
2988  eprds(k)=prds(k)
2989  prds(k)=0.
2990  END IF
2991  IF (prdg(k).LT.0.) THEN
2992  eprdg(k)=prdg(k)
2993  prdg(k)=0.
2994  END IF
2995 !.......................................................................
2996 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
2997 
2998 ! CONSERVATION OF WATER
2999 ! THIS IS ADOPTED LOOSELY FROM MM5 RESINER CODE. HOWEVER, HERE WE
3000 ! ONLY ADJUST PROCESSES THAT ARE NEGATIVE, RATHER THAN ALL PROCESSES.
3001 
3002 ! IF MIXING RATIOS LESS THAN QSMALL, THEN NO DEPLETION OF WATER
3003 ! THROUGH MICROPHYSICAL PROCESSES, SKIP CONSERVATION
3004 
3005 ! NOTE: CONSERVATION CHECK NOT APPLIED TO NUMBER CONCENTRATION SPECIES. ADDITIONAL CATCH
3006 ! BELOW WILL PREVENT NEGATIVE NUMBER CONCENTRATION
3007 ! FOR EACH MICROPHYSICAL PROCESS WHICH PROVIDES A SOURCE FOR NUMBER, THERE IS A CHECK
3008 ! TO MAKE SURE THAT CAN'T EXCEED TOTAL NUMBER OF DEPLETED SPECIES WITH THE TIME
3009 ! STEP
3010 
3011 !****SENSITIVITY - NO ICE
3012 
3013  IF (iliq.EQ.1) THEN
3014  mnuccc(k)=0.
3015  nnuccc(k)=0.
3016  mnuccr(k)=0.
3017  nnuccr(k)=0.
3018  mnuccd(k)=0.
3019  nnuccd(k)=0.
3020  END IF
3021 
3022 ! ****SENSITIVITY - NO GRAUPEL
3023  IF (igraup.EQ.1) THEN
3024  pracg(k) = 0.
3025  psacr(k) = 0.
3026  psacwg(k) = 0.
3027  prdg(k) = 0.
3028  eprdg(k) = 0.
3029  evpmg(k) = 0.
3030  pgmlt(k) = 0.
3031  npracg(k) = 0.
3032  npsacwg(k) = 0.
3033  nscng(k) = 0.
3034  ngracs(k) = 0.
3035  nsubg(k) = 0.
3036  ngmltg(k) = 0.
3037  ngmltr(k) = 0.
3038 ! fix 053011
3039  piacrs(k)=piacrs(k)+piacr(k)
3040  piacr(k) = 0.
3041 ! fix 070713
3042  pracis(k)=pracis(k)+praci(k)
3043  praci(k) = 0.
3044  psacws(k)=psacws(k)+pgsacw(k)
3045  pgsacw(k) = 0.
3046  pracs(k)=pracs(k)+pgracs(k)
3047  pgracs(k) = 0.
3048  END IF
3049 
3050 ! CONSERVATION OF QC
3051 
3052  dum = (prc(k)+pra(k)+mnuccc(k)+psacws(k)+psacwi(k)+qmults(k)+psacwg(k)+pgsacw(k)+qmultg(k))*dt
3053 
3054  IF (dum.GT.qc3d(k).AND.qc3d(k).GE.qsmall) THEN
3055  ratio = qc3d(k)/dum
3056 
3057  prc(k) = prc(k)*ratio
3058  pra(k) = pra(k)*ratio
3059  mnuccc(k) = mnuccc(k)*ratio
3060  psacws(k) = psacws(k)*ratio
3061  psacwi(k) = psacwi(k)*ratio
3062  qmults(k) = qmults(k)*ratio
3063  qmultg(k) = qmultg(k)*ratio
3064  psacwg(k) = psacwg(k)*ratio
3065  pgsacw(k) = pgsacw(k)*ratio
3066  END IF
3067 
3068 ! CONSERVATION OF QI
3069 
3070  dum = (-prd(k)-mnuccc(k)+prci(k)+prai(k)-qmults(k)-qmultg(k)-qmultr(k)-qmultrg(k) &
3071  -mnuccd(k)+praci(k)+pracis(k)-eprd(k)-psacwi(k))*dt
3072 
3073  IF (dum.GT.qi3d(k).AND.qi3d(k).GE.qsmall) THEN
3074 
3075  ratio = (qi3d(k)/dt+prd(k)+mnuccc(k)+qmults(k)+qmultg(k)+qmultr(k)+qmultrg(k)+ &
3076  mnuccd(k)+psacwi(k))/ &
3077  (prci(k)+prai(k)+praci(k)+pracis(k)-eprd(k))
3078 
3079  prci(k) = prci(k)*ratio
3080  prai(k) = prai(k)*ratio
3081  praci(k) = praci(k)*ratio
3082  pracis(k) = pracis(k)*ratio
3083  eprd(k) = eprd(k)*ratio
3084 
3085  END IF
3086 
3087 ! CONSERVATION OF QR
3088 
3089  dum=((pracs(k)-pre(k))+(qmultr(k)+qmultrg(k)-prc(k))+(mnuccr(k)-pra(k))+ &
3090  piacr(k)+piacrs(k)+pgracs(k)+pracg(k))*dt
3091 
3092  IF (dum.GT.qr3d(k).AND.qr3d(k).GE.qsmall) THEN
3093 
3094  ratio = (qr3d(k)/dt+prc(k)+pra(k))/ &
3095  (-pre(k)+qmultr(k)+qmultrg(k)+pracs(k)+mnuccr(k)+piacr(k)+piacrs(k)+pgracs(k)+pracg(k))
3096 
3097  pre(k) = pre(k)*ratio
3098  pracs(k) = pracs(k)*ratio
3099  qmultr(k) = qmultr(k)*ratio
3100  qmultrg(k) = qmultrg(k)*ratio
3101  mnuccr(k) = mnuccr(k)*ratio
3102  piacr(k) = piacr(k)*ratio
3103  piacrs(k) = piacrs(k)*ratio
3104  pgracs(k) = pgracs(k)*ratio
3105  pracg(k) = pracg(k)*ratio
3106 
3107  END IF
3108 
3109 ! CONSERVATION OF QNI
3110 ! CONSERVATION FOR GRAUPEL SCHEME
3111 
3112  IF (igraup.EQ.0) THEN
3113 
3114  dum = (-prds(k)-psacws(k)-prai(k)-prci(k)-pracs(k)-eprds(k)+psacr(k)-piacrs(k)-pracis(k))*dt
3115 
3116  IF (dum.GT.qni3d(k).AND.qni3d(k).GE.qsmall) THEN
3117 
3118  ratio = (qni3d(k)/dt+prds(k)+psacws(k)+prai(k)+prci(k)+pracs(k)+piacrs(k)+pracis(k))/(-eprds(k)+psacr(k))
3119 
3120  eprds(k) = eprds(k)*ratio
3121  psacr(k) = psacr(k)*ratio
3122 
3123  END IF
3124 
3125 ! FOR NO GRAUPEL, NEED TO INCLUDE FREEZING OF RAIN FOR SNOW
3126  ELSE IF (igraup.EQ.1) THEN
3127 
3128  dum = (-prds(k)-psacws(k)-prai(k)-prci(k)-pracs(k)-eprds(k)+psacr(k)-piacrs(k)-pracis(k)-mnuccr(k))*dt
3129 
3130  IF (dum.GT.qni3d(k).AND.qni3d(k).GE.qsmall) THEN
3131 
3132  ratio = (qni3d(k)/dt+prds(k)+psacws(k)+prai(k)+prci(k)+pracs(k)+piacrs(k)+pracis(k)+mnuccr(k))/(-eprds(k)+psacr(k))
3133 
3134  eprds(k) = eprds(k)*ratio
3135  psacr(k) = psacr(k)*ratio
3136 
3137  END IF
3138 
3139  END IF
3140 
3141 ! CONSERVATION OF QG
3142 
3143  dum = (-psacwg(k)-pracg(k)-pgsacw(k)-pgracs(k)-prdg(k)-mnuccr(k)-eprdg(k)-piacr(k)-praci(k)-psacr(k))*dt
3144 
3145  IF (dum.GT.qg3d(k).AND.qg3d(k).GE.qsmall) THEN
3146 
3147  ratio = (qg3d(k)/dt+psacwg(k)+pracg(k)+pgsacw(k)+pgracs(k)+prdg(k)+mnuccr(k)+psacr(k)+&
3148  piacr(k)+praci(k))/(-eprdg(k))
3149 
3150  eprdg(k) = eprdg(k)*ratio
3151 
3152  END IF
3153 
3154 ! TENDENCIES
3155 
3156  qv3dten(k) = qv3dten(k)+(-pre(k)-prd(k)-prds(k)-mnuccd(k)-eprd(k)-eprds(k)-prdg(k)-eprdg(k))
3157 
3158 ! BUG FIX HM, 3/1/11, INCLUDE PIACR AND PIACRS
3159  t3dten(k) = t3dten(k)+(pre(k) &
3160  *xxlv(k)+(prd(k)+prds(k)+ &
3161  mnuccd(k)+eprd(k)+eprds(k)+prdg(k)+eprdg(k))*xxls(k)+ &
3162  (psacws(k)+psacwi(k)+mnuccc(k)+mnuccr(k)+ &
3163  qmults(k)+qmultg(k)+qmultr(k)+qmultrg(k)+pracs(k) &
3164  +psacwg(k)+pracg(k)+pgsacw(k)+pgracs(k)+piacr(k)+piacrs(k))*xlf(k))/cpm(k)
3165 
3166  qc3dten(k) = qc3dten(k)+ &
3167  (-pra(k)-prc(k)-mnuccc(k)+pcc(k)- &
3168  psacws(k)-psacwi(k)-qmults(k)-qmultg(k)-psacwg(k)-pgsacw(k))
3169  qi3dten(k) = qi3dten(k)+ &
3170  (prd(k)+eprd(k)+psacwi(k)+mnuccc(k)-prci(k)- &
3171  prai(k)+qmults(k)+qmultg(k)+qmultr(k)+qmultrg(k)+mnuccd(k)-praci(k)-pracis(k))
3172  qr3dten(k) = qr3dten(k)+ &
3173  (pre(k)+pra(k)+prc(k)-pracs(k)-mnuccr(k)-qmultr(k)-qmultrg(k) &
3174  -piacr(k)-piacrs(k)-pracg(k)-pgracs(k))
3175  IF (igraup.EQ.0) THEN
3176 
3177  qni3dten(k) = qni3dten(k)+ &
3178  (prai(k)+psacws(k)+prds(k)+pracs(k)+prci(k)+eprds(k)-psacr(k)+piacrs(k)+pracis(k))
3179  ns3dten(k) = ns3dten(k)+(nsagg(k)+nprci(k)-nscng(k)-ngracs(k)+niacrs(k))
3180  qg3dten(k) = qg3dten(k)+(pracg(k)+psacwg(k)+pgsacw(k)+pgracs(k)+ &
3181  prdg(k)+eprdg(k)+mnuccr(k)+piacr(k)+praci(k)+psacr(k))
3182  ng3dten(k) = ng3dten(k)+(nscng(k)+ngracs(k)+nnuccr(k)+niacr(k))
3183 
3184 ! FOR NO GRAUPEL, NEED TO INCLUDE FREEZING OF RAIN FOR SNOW
3185  ELSE IF (igraup.EQ.1) THEN
3186 
3187  qni3dten(k) = qni3dten(k)+ &
3188  (prai(k)+psacws(k)+prds(k)+pracs(k)+prci(k)+eprds(k)-psacr(k)+piacrs(k)+pracis(k)+mnuccr(k))
3189  ns3dten(k) = ns3dten(k)+(nsagg(k)+nprci(k)-nscng(k)-ngracs(k)+niacrs(k)+nnuccr(k))
3190 
3191  END IF
3192 
3193  nc3dten(k) = nc3dten(k)+(-nnuccc(k)-npsacws(k) &
3194  -npra(k)-nprc(k)-npsacwi(k)-npsacwg(k))
3195 
3196  ni3dten(k) = ni3dten(k)+ &
3197  (nnuccc(k)-nprci(k)-nprai(k)+nmults(k)+nmultg(k)+nmultr(k)+nmultrg(k)+ &
3198  nnuccd(k)-niacr(k)-niacrs(k))
3199 
3200  nr3dten(k) = nr3dten(k)+(nprc1(k)-npracs(k)-nnuccr(k) &
3201  +nragg(k)-niacr(k)-niacrs(k)-npracg(k)-ngracs(k))
3202 
3203 ! HM ADD, WRF-CHEM, ADD TENDENCIES FOR C2PREC
3204 
3205  c2prec(k) = pra(k)+prc(k)+psacws(k)+qmults(k)+qmultg(k)+psacwg(k)+ &
3206  pgsacw(k)+mnuccc(k)+psacwi(k)
3207 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
3208 ! NOW CALCULATE SATURATION ADJUSTMENT TO CONDENSE EXTRA VAPOR ABOVE
3209 ! WATER SATURATION
3210 
3211  dumt = t3d(k)+dt*t3dten(k)
3212  dumqv = qv3d(k) + dt * qv3dten(k)
3213 
3214  ! hm, add fix for low pressure, 5/12/10
3215  dum=min(0.99*pres(k),polysvp(dumt,0))
3216  dumqss = ep_2*dum/(pres(k)-dum)
3217 
3218  dumqc = qc3d(k) + dt * qc3dten(k)
3219 
3220  dumqc = max(dumqc,0.)
3221 
3222 ! SATURATION ADJUSTMENT FOR LIQUID
3223 
3224  dums = dumqv-dumqss
3225 
3226  pcc(k) = dums/(1.+xxlv(k)**2*dumqss/(cpm(k)*rv*dumt**2))/dt
3227 
3228  IF (pcc(k)*dt+dumqc.LT.0.) THEN
3229  pcc(k) = -dumqc/dt
3230  END IF
3231 
3232  qv3dten(k) = qv3dten(k)-pcc(k)
3233  t3dten(k) = t3dten(k)+pcc(k)*xxlv(k)/cpm(k)
3234  qc3dten(k) = qc3dten(k)+pcc(k)
3235 ! Fortran version
3236 !.......................................................................
3237 ! ACTIVATION OF CLOUD DROPLETS
3238 ! ACTIVATION OF DROPLET CURRENTLY NOT CALCULATED
3239 ! DROPLET CONCENTRATION IS SPECIFIED !!!!!
3240 
3241 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
3242 ! SUBLIMATE, MELT, OR EVAPORATE NUMBER CONCENTRATION
3243 ! THIS FORMULATION ASSUMES 1:1 RATIO BETWEEN MASS LOSS AND
3244 ! LOSS OF NUMBER CONCENTRATION
3245 
3246 ! IF (PCC(K).LT.0.) THEN
3247 ! DUM = PCC(K)*DT/QC3D(K)
3248 ! DUM = MAX(-1.,DUM)
3249 ! NSUBC(K) = DUM*NC3D(K)/DT
3250 ! END IF
3251 
3252  IF (eprd(k).LT.0.) THEN
3253  dum = eprd(k)*dt/qi3d(k)
3254  dum = max(-1.,dum)
3255  nsubi(k) = dum*ni3d(k)/dt
3256  END IF
3257  IF (eprds(k).LT.0.) THEN
3258  dum = eprds(k)*dt/qni3d(k)
3259  dum = max(-1.,dum)
3260  nsubs(k) = dum*ns3d(k)/dt
3261  END IF
3262  IF (pre(k).LT.0.) THEN
3263  dum = pre(k)*dt/qr3d(k)
3264  dum = max(-1.,dum)
3265  nsubr(k) = dum*nr3d(k)/dt
3266  END IF
3267  IF (eprdg(k).LT.0.) THEN
3268  dum = eprdg(k)*dt/qg3d(k)
3269  dum = max(-1.,dum)
3270  nsubg(k) = dum*ng3d(k)/dt
3271  END IF
3272 
3273 ! nsubr(k)=0.
3274 ! nsubs(k)=0.
3275 ! nsubg(k)=0.
3276 
3277 ! UPDATE TENDENCIES
3278 
3279 ! NC3DTEN(K) = NC3DTEN(K)+NSUBC(K)
3280  ni3dten(k) = ni3dten(k)+nsubi(k)
3281  ns3dten(k) = ns3dten(k)+nsubs(k)
3282  ng3dten(k) = ng3dten(k)+nsubg(k)
3283  nr3dten(k) = nr3dten(k)+nsubr(k)
3284 
3285  END IF !!!!!! TEMPERATURE
3286 
3287 ! SWITCH LTRUE TO 1, SINCE HYDROMETEORS ARE PRESENT
3288  ltrue = 1
3289 
3290  200 CONTINUE
3291 
3292  END DO
3293 
3294 ! INITIALIZE PRECIP AND SNOW RATES
3295  precrt = 0.
3296  snowrt = 0.
3297 ! hm added 7/13/13
3298  snowprt = 0.
3299  grplprt = 0.
3300 
3301 ! IF THERE ARE NO HYDROMETEORS, THEN SKIP TO END OF SUBROUTINE
3302  DO k = kts,kte
3303 ! Fortran version
3304  END DO
3305  IF (ltrue.EQ.0) GOTO 400
3306 
3307 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
3308 !.......................................................................
3309 ! CALCULATE SEDIMENATION
3310 ! THE NUMERICS HERE FOLLOW FROM REISNER ET AL. (1998)
3311 ! FALLOUT TERMS ARE CALCULATED ON SPLIT TIME STEPS TO ENSURE NUMERICAL
3312 ! STABILITY, I.E. COURANT# < 1
3313 
3314 !.......................................................................
3315 
3316  nstep = 1
3317 
3318  DO k = kte,kts,-1
3319 
3320  dumi(k) = qi3d(k)+qi3dten(k)*dt
3321  dumqs(k) = qni3d(k)+qni3dten(k)*dt
3322  dumr(k) = qr3d(k)+qr3dten(k)*dt
3323  dumfni(k) = ni3d(k)+ni3dten(k)*dt
3324  dumfns(k) = ns3d(k)+ns3dten(k)*dt
3325  dumfnr(k) = nr3d(k)+nr3dten(k)*dt
3326  dumc(k) = qc3d(k)+qc3dten(k)*dt
3327  dumfnc(k) = nc3d(k)+nc3dten(k)*dt
3328  dumg(k) = qg3d(k)+qg3dten(k)*dt
3329  dumfng(k) = ng3d(k)+ng3dten(k)*dt
3330 
3331 ! SWITCH FOR CONSTANT DROPLET NUMBER
3332  IF (iinum.EQ.1) THEN
3333  dumfnc(k) = nc3d(k)
3334  END IF
3335 
3336 ! GET DUMMY LAMDA FOR SEDIMENTATION CALCULATIONS
3337 
3338 ! MAKE SURE NUMBER CONCENTRATIONS ARE POSITIVE
3339  dumfni(k) = max(0.,dumfni(k))
3340  dumfns(k) = max(0.,dumfns(k))
3341  dumfnc(k) = max(0.,dumfnc(k))
3342  dumfnr(k) = max(0.,dumfnr(k))
3343  dumfng(k) = max(0.,dumfng(k))
3344 
3345 !......................................................................
3346 ! CLOUD ICE
3347 
3348  IF (dumi(k).GE.qsmall) THEN
3349  dlami = (cons12*dumfni(k)/dumi(k))**(1./di)
3350  dlami=max(dlami,lammini)
3351  dlami=min(dlami,lammaxi)
3352  END IF
3353 !......................................................................
3354 ! RAIN
3355 
3356  IF (dumr(k).GE.qsmall) THEN
3357  dlamr = (pi*rhow*dumfnr(k)/dumr(k))**(1./3.)
3358  dlamr=max(dlamr,lamminr)
3359  dlamr=min(dlamr,lammaxr)
3360  END IF
3361 !......................................................................
3362 ! CLOUD DROPLETS
3363 
3364  IF (dumc(k).GE.qsmall) THEN
3365  dum = pres(k)/(287.15*t3d(k))
3366  pgam(k)=0.0005714*(nc3d(k)/1.e6*dum)+0.2714
3367  pgam(k)=1./(pgam(k)**2)-1.
3368  pgam(k)=max(pgam(k),2.)
3369  pgam(k)=min(pgam(k),10.)
3370 
3371  dlamc = (cons26*dumfnc(k)*gamma(pgam(k)+4.)/(dumc(k)*gamma(pgam(k)+1.)))**(1./3.)
3372  lammin = (pgam(k)+1.)/60.e-6
3373  lammax = (pgam(k)+1.)/1.e-6
3374  dlamc=max(dlamc,lammin)
3375  dlamc=min(dlamc,lammax)
3376  END IF
3377 !......................................................................
3378 ! SNOW
3379 
3380  IF (dumqs(k).GE.qsmall) THEN
3381  dlams = (cons1*dumfns(k)/ dumqs(k))**(1./ds)
3382  dlams=max(dlams,lammins)
3383  dlams=min(dlams,lammaxs)
3384  END IF
3385 !......................................................................
3386 ! GRAUPEL
3387 
3388  IF (dumg(k).GE.qsmall) THEN
3389  dlamg = (cons2*dumfng(k)/ dumg(k))**(1./dg)
3390  dlamg=max(dlamg,lamming)
3391  dlamg=min(dlamg,lammaxg)
3392  END IF
3393 !......................................................................
3394 ! CALCULATE NUMBER-WEIGHTED AND MASS-WEIGHTED TERMINAL FALL SPEEDS
3395 
3396 ! CLOUD WATER
3397 
3398  IF (dumc(k).GE.qsmall) THEN
3399  unc = acn(k)*gamma(1.+bc+pgam(k))/ (dlamc**bc*gamma(pgam(k)+1.))
3400  umc = acn(k)*gamma(4.+bc+pgam(k))/ (dlamc**bc*gamma(pgam(k)+4.))
3401  ELSE
3402  umc = 0.
3403  unc = 0.
3404  END IF
3405 
3406  IF (dumi(k).GE.qsmall) THEN
3407  uni = ain(k)*cons27/dlami**bi
3408  umi = ain(k)*cons28/(dlami**bi)
3409  ELSE
3410  umi = 0.
3411  uni = 0.
3412  END IF
3413 
3414  IF (dumr(k).GE.qsmall) THEN
3415  unr = arn(k)*cons6/dlamr**br
3416  umr = arn(k)*cons4/(dlamr**br)
3417  ELSE
3418  umr = 0.
3419  unr = 0.
3420  END IF
3421 
3422  IF (dumqs(k).GE.qsmall) THEN
3423  ums = asn(k)*cons3/(dlams**bs)
3424  uns = asn(k)*cons5/dlams**bs
3425  ELSE
3426  ums = 0.
3427  uns = 0.
3428  END IF
3429 
3430  IF (dumg(k).GE.qsmall) THEN
3431  umg = agn(k)*cons7/(dlamg**bg)
3432  ung = agn(k)*cons8/dlamg**bg
3433  ELSE
3434  umg = 0.
3435  ung = 0.
3436  END IF
3437 
3438 ! SET REALISTIC LIMITS ON FALLSPEED
3439 
3440 ! bug fix, 10/08/09
3441  dum=(rhosu/rho(k))**0.54
3442  ums=min(ums,1.2*dum)
3443  uns=min(uns,1.2*dum)
3444 ! fix 053011
3445 ! fix for correction by AA 4/6/11
3446  umi=min(umi,1.2*(rhosu/rho(k))**0.35)
3447  uni=min(uni,1.2*(rhosu/rho(k))**0.35)
3448  umr=min(umr,9.1*dum)
3449  unr=min(unr,9.1*dum)
3450  umg=min(umg,20.*dum)
3451  ung=min(ung,20.*dum)
3452 
3453  fr(k) = umr
3454  fi(k) = umi
3455  fni(k) = uni
3456  fs(k) = ums
3457  fns(k) = uns
3458  fnr(k) = unr
3459  fc(k) = umc
3460  fnc(k) = unc
3461  fg(k) = umg
3462  fng(k) = ung
3463 
3464 ! V3.3 MODIFY FALLSPEED BELOW LEVEL OF PRECIP
3465 
3466  IF (k.LE.kte-1) THEN
3467  IF (fr(k).LT.1.e-10) THEN
3468  fr(k)=fr(k+1)
3469  END IF
3470  IF (fi(k).LT.1.e-10) THEN
3471  fi(k)=fi(k+1)
3472  END IF
3473  IF (fni(k).LT.1.e-10) THEN
3474  fni(k)=fni(k+1)
3475  END IF
3476  IF (fs(k).LT.1.e-10) THEN
3477  fs(k)=fs(k+1)
3478  END IF
3479  IF (fns(k).LT.1.e-10) THEN
3480  fns(k)=fns(k+1)
3481  END IF
3482  IF (fnr(k).LT.1.e-10) THEN
3483  fnr(k)=fnr(k+1)
3484  END IF
3485  IF (fc(k).LT.1.e-10) THEN
3486  fc(k)=fc(k+1)
3487  END IF
3488  IF (fnc(k).LT.1.e-10) THEN
3489  fnc(k)=fnc(k+1)
3490  END IF
3491  IF (fg(k).LT.1.e-10) THEN
3492  fg(k)=fg(k+1)
3493  END IF
3494  IF (fng(k).LT.1.e-10) THEN
3495  fng(k)=fng(k+1)
3496  END IF
3497  END IF ! K LE KTE-1
3498 
3499 ! CALCULATE NUMBER OF SPLIT TIME STEPS
3500 
3501  rgvm = max(fr(k),fi(k),fs(k),fc(k),fni(k),fnr(k),fns(k),fnc(k),fg(k),fng(k))
3502 ! VVT CHANGED IFIX -> INT (GENERIC FUNCTION)
3503  nstep = max(int(rgvm*dt/dzq(k)+1.),nstep)
3504 ! Fortran version
3505 ! MULTIPLY VARIABLES BY RHO
3506  dumr(k) = dumr(k)*rho(k)
3507  dumi(k) = dumi(k)*rho(k)
3508  dumfni(k) = dumfni(k)*rho(k)
3509  dumqs(k) = dumqs(k)*rho(k)
3510  dumfns(k) = dumfns(k)*rho(k)
3511  dumfnr(k) = dumfnr(k)*rho(k)
3512  dumc(k) = dumc(k)*rho(k)
3513  dumfnc(k) = dumfnc(k)*rho(k)
3514  dumg(k) = dumg(k)*rho(k)
3515  dumfng(k) = dumfng(k)*rho(k)
3516 
3517  END DO
3518 
3519  DO n = 1,nstep
3520 
3521  DO k = kts,kte
3522  faloutr(k) = fr(k)*dumr(k)
3523  falouti(k) = fi(k)*dumi(k)
3524  faloutni(k) = fni(k)*dumfni(k)
3525  falouts(k) = fs(k)*dumqs(k)
3526  faloutns(k) = fns(k)*dumfns(k)
3527  faloutnr(k) = fnr(k)*dumfnr(k)
3528  faloutc(k) = fc(k)*dumc(k)
3529  faloutnc(k) = fnc(k)*dumfnc(k)
3530  faloutg(k) = fg(k)*dumg(k)
3531  faloutng(k) = fng(k)*dumfng(k)
3532  END DO
3533 
3534 ! TOP OF MODEL
3535 
3536  k = kte
3537  faltndr = faloutr(k)/dzq(k)
3538  faltndi = falouti(k)/dzq(k)
3539  faltndni = faloutni(k)/dzq(k)
3540  faltnds = falouts(k)/dzq(k)
3541  faltndns = faloutns(k)/dzq(k)
3542  faltndnr = faloutnr(k)/dzq(k)
3543  faltndc = faloutc(k)/dzq(k)
3544  faltndnc = faloutnc(k)/dzq(k)
3545  faltndg = faloutg(k)/dzq(k)
3546  faltndng = faloutng(k)/dzq(k)
3547 ! ADD FALLOUT TERMS TO EULERIAN TENDENCIES
3548 
3549  qrsten(k) = qrsten(k)-faltndr/nstep/rho(k)
3550  qisten(k) = qisten(k)-faltndi/nstep/rho(k)
3551  ni3dten(k) = ni3dten(k)-faltndni/nstep/rho(k)
3552  qnisten(k) = qnisten(k)-faltnds/nstep/rho(k)
3553  ns3dten(k) = ns3dten(k)-faltndns/nstep/rho(k)
3554  nr3dten(k) = nr3dten(k)-faltndnr/nstep/rho(k)
3555  qcsten(k) = qcsten(k)-faltndc/nstep/rho(k)
3556  nc3dten(k) = nc3dten(k)-faltndnc/nstep/rho(k)
3557  qgsten(k) = qgsten(k)-faltndg/nstep/rho(k)
3558  ng3dten(k) = ng3dten(k)-faltndng/nstep/rho(k)
3559 
3560  dumr(k) = dumr(k)-faltndr*dt/nstep
3561  dumi(k) = dumi(k)-faltndi*dt/nstep
3562  dumfni(k) = dumfni(k)-faltndni*dt/nstep
3563  dumqs(k) = dumqs(k)-faltnds*dt/nstep
3564  dumfns(k) = dumfns(k)-faltndns*dt/nstep
3565  dumfnr(k) = dumfnr(k)-faltndnr*dt/nstep
3566  dumc(k) = dumc(k)-faltndc*dt/nstep
3567  dumfnc(k) = dumfnc(k)-faltndnc*dt/nstep
3568  dumg(k) = dumg(k)-faltndg*dt/nstep
3569  dumfng(k) = dumfng(k)-faltndng*dt/nstep
3570 
3571  DO k = kte-1,kts,-1
3572  faltndr = (faloutr(k+1)-faloutr(k))/dzq(k)
3573  faltndi = (falouti(k+1)-falouti(k))/dzq(k)
3574  faltndni = (faloutni(k+1)-faloutni(k))/dzq(k)
3575  faltnds = (falouts(k+1)-falouts(k))/dzq(k)
3576  faltndns = (faloutns(k+1)-faloutns(k))/dzq(k)
3577  faltndnr = (faloutnr(k+1)-faloutnr(k))/dzq(k)
3578  faltndc = (faloutc(k+1)-faloutc(k))/dzq(k)
3579  faltndnc = (faloutnc(k+1)-faloutnc(k))/dzq(k)
3580  faltndg = (faloutg(k+1)-faloutg(k))/dzq(k)
3581  faltndng = (faloutng(k+1)-faloutng(k))/dzq(k)
3582 
3583 ! ADD FALLOUT TERMS TO EULERIAN TENDENCIES
3584 
3585  qrsten(k) = qrsten(k)+faltndr/nstep/rho(k)
3586  qisten(k) = qisten(k)+faltndi/nstep/rho(k)
3587  ni3dten(k) = ni3dten(k)+faltndni/nstep/rho(k)
3588  qnisten(k) = qnisten(k)+faltnds/nstep/rho(k)
3589  ns3dten(k) = ns3dten(k)+faltndns/nstep/rho(k)
3590  nr3dten(k) = nr3dten(k)+faltndnr/nstep/rho(k)
3591  qcsten(k) = qcsten(k)+faltndc/nstep/rho(k)
3592  nc3dten(k) = nc3dten(k)+faltndnc/nstep/rho(k)
3593  qgsten(k) = qgsten(k)+faltndg/nstep/rho(k)
3594  ng3dten(k) = ng3dten(k)+faltndng/nstep/rho(k)
3595 ! Fortran version
3596  dumr(k) = dumr(k)+faltndr*dt/nstep
3597  dumi(k) = dumi(k)+faltndi*dt/nstep
3598  dumfni(k) = dumfni(k)+faltndni*dt/nstep
3599  dumqs(k) = dumqs(k)+faltnds*dt/nstep
3600  dumfns(k) = dumfns(k)+faltndns*dt/nstep
3601  dumfnr(k) = dumfnr(k)+faltndnr*dt/nstep
3602  dumc(k) = dumc(k)+faltndc*dt/nstep
3603  dumfnc(k) = dumfnc(k)+faltndnc*dt/nstep
3604  dumg(k) = dumg(k)+faltndg*dt/nstep
3605  dumfng(k) = dumfng(k)+faltndng*dt/nstep
3606 
3607 ! FOR WRF-CHEM, NEED PRECIP RATES (UNITS OF KG/M^2/S)
3608  csed(k)=csed(k)+faloutc(k)/nstep
3609  ised(k)=ised(k)+falouti(k)/nstep
3610  ssed(k)=ssed(k)+falouts(k)/nstep
3611  gsed(k)=gsed(k)+faloutg(k)/nstep
3612  rsed(k)=rsed(k)+faloutr(k)/nstep
3613  END DO
3614 
3615 ! GET PRECIPITATION AND SNOWFALL ACCUMULATION DURING THE TIME STEP
3616 ! FACTOR OF 1000 CONVERTS FROM M TO MM, BUT DIVISION BY DENSITY
3617 ! OF LIQUID WATER CANCELS THIS FACTOR OF 1000
3618 
3619  precrt = precrt+(faloutr(kts)+faloutc(kts)+falouts(kts)+falouti(kts)+faloutg(kts)) &
3620  *dt/nstep
3621  snowrt = snowrt+(falouts(kts)+falouti(kts)+faloutg(kts))*dt/nstep
3622 ! hm added 7/13/13
3623  snowprt = snowprt+(falouti(kts)+falouts(kts))*dt/nstep
3624  grplprt = grplprt+(faloutg(kts))*dt/nstep
3625  END DO
3626 
3627  DO k=kts,kte
3628 ! Fortran version
3629 ! ADD ON SEDIMENTATION TENDENCIES FOR MIXING RATIO TO REST OF TENDENCIES
3630 
3631  qr3dten(k)=qr3dten(k)+qrsten(k)
3632  qi3dten(k)=qi3dten(k)+qisten(k)
3633  qc3dten(k)=qc3dten(k)+qcsten(k)
3634  qg3dten(k)=qg3dten(k)+qgsten(k)
3635  qni3dten(k)=qni3dten(k)+qnisten(k)
3636 ! Fortran version
3637 ! PUT ALL CLOUD ICE IN SNOW CATEGORY IF MEAN DIAMETER EXCEEDS 2 * dcs
3638 
3639 !hm 4/7/09 bug fix
3640 ! IF (QI3D(K).GE.QSMALL.AND.T3D(K).LT.273.15) THEN
3641  IF (qi3d(k).GE.qsmall.AND.t3d(k).LT.273.15.AND.lami(k).GE.1.e-10) THEN
3642  IF (1./lami(k).GE.2.*dcs) THEN
3643  qni3dten(k) = qni3dten(k)+qi3d(k)/dt+ qi3dten(k)
3644  ns3dten(k) = ns3dten(k)+ni3d(k)/dt+ ni3dten(k)
3645  qi3dten(k) = -qi3d(k)/dt
3646  ni3dten(k) = -ni3d(k)/dt
3647  END IF
3648  END IF
3649 
3650 ! hm add tendencies here, then call sizeparameter
3651 ! to ensure consisitency between mixing ratio and number concentration
3652 
3653  qc3d(k) = qc3d(k)+qc3dten(k)*dt
3654  qi3d(k) = qi3d(k)+qi3dten(k)*dt
3655  qni3d(k) = qni3d(k)+qni3dten(k)*dt
3656  qr3d(k) = qr3d(k)+qr3dten(k)*dt
3657  nc3d(k) = nc3d(k)+nc3dten(k)*dt
3658  ni3d(k) = ni3d(k)+ni3dten(k)*dt
3659  ns3d(k) = ns3d(k)+ns3dten(k)*dt
3660  nr3d(k) = nr3d(k)+nr3dten(k)*dt
3661 ! Fortran version
3662  IF (igraup.EQ.0) THEN
3663  qg3d(k) = qg3d(k)+qg3dten(k)*dt
3664  ng3d(k) = ng3d(k)+ng3dten(k)*dt
3665  END IF
3666 
3667 ! ADD TEMPERATURE AND WATER VAPOR TENDENCIES FROM MICROPHYSICS
3668  t3d(k) = t3d(k)+t3dten(k)*dt
3669  qv3d(k) = qv3d(k)+qv3dten(k)*dt
3670 ! Fortran version
3671 ! SATURATION VAPOR PRESSURE AND MIXING RATIO
3672 
3673 ! hm, add fix for low pressure, 5/12/10
3674  evs(k) = min(0.99*pres(k),polysvp(t3d(k),0)) ! PA
3675  eis(k) = min(0.99*pres(k),polysvp(t3d(k),1)) ! PA
3676 
3677 ! MAKE SURE ICE SATURATION DOESN'T EXCEED WATER SAT. NEAR FREEZING
3678 
3679  IF (eis(k).GT.evs(k)) eis(k) = evs(k)
3680 
3681  qvs(k) = ep_2*evs(k)/(pres(k)-evs(k))
3682  qvi(k) = ep_2*eis(k)/(pres(k)-eis(k))
3683 
3684  qvqvs(k) = qv3d(k)/qvs(k)
3685  qvqvsi(k) = qv3d(k)/qvi(k)
3686 ! Fortran version
3687 ! AT SUBSATURATION, REMOVE SMALL AMOUNTS OF CLOUD/PRECIP WATER
3688 ! hm 7/9/09 change limit to 1.e-8
3689 
3690  IF (qvqvs(k).LT.0.9) THEN
3691  IF (qr3d(k).LT.1.e-8) THEN
3692  qv3d(k)=qv3d(k)+qr3d(k)
3693  t3d(k)=t3d(k)-qr3d(k)*xxlv(k)/cpm(k)
3694  qr3d(k)=0.
3695  END IF
3696  IF (qc3d(k).LT.1.e-8) THEN
3697  qv3d(k)=qv3d(k)+qc3d(k)
3698  t3d(k)=t3d(k)-qc3d(k)*xxlv(k)/cpm(k)
3699  qc3d(k)=0.
3700  END IF
3701  END IF
3702  IF (qvqvsi(k).LT.0.9) THEN
3703  IF (qi3d(k).LT.1.e-8) THEN
3704  qv3d(k)=qv3d(k)+qi3d(k)
3705  t3d(k)=t3d(k)-qi3d(k)*xxls(k)/cpm(k)
3706  qi3d(k)=0.
3707  END IF
3708  IF (qni3d(k).LT.1.e-8) THEN
3709  qv3d(k)=qv3d(k)+qni3d(k)
3710  t3d(k)=t3d(k)-qni3d(k)*xxls(k)/cpm(k)
3711  qni3d(k)=0.
3712  END IF
3713  IF (qg3d(k).LT.1.e-8) THEN
3714  qv3d(k)=qv3d(k)+qg3d(k)
3715  t3d(k)=t3d(k)-qg3d(k)*xxls(k)/cpm(k)
3716  qg3d(k)=0.
3717  END IF
3718  END IF
3719 ! Fortran version
3720 !..................................................................
3721 ! IF MIXING RATIO < QSMALL SET MIXING RATIO AND NUMBER CONC TO ZERO
3722 
3723  IF (qc3d(k).LT.qsmall) THEN
3724  qc3d(k) = 0.
3725  nc3d(k) = 0.
3726  effc(k) = 0.
3727  END IF
3728  IF (qr3d(k).LT.qsmall) THEN
3729  qr3d(k) = 0.
3730  nr3d(k) = 0.
3731  effr(k) = 0.
3732  END IF
3733  IF (qi3d(k).LT.qsmall) THEN
3734  qi3d(k) = 0.
3735  ni3d(k) = 0.
3736  effi(k) = 0.
3737  END IF
3738  IF (qni3d(k).LT.qsmall) THEN
3739  qni3d(k) = 0.
3740  ns3d(k) = 0.
3741  effs(k) = 0.
3742  END IF
3743  IF (qg3d(k).LT.qsmall) THEN
3744  qg3d(k) = 0.
3745  ng3d(k) = 0.
3746  effg(k) = 0.
3747  END IF
3748 ! Fortran version
3749 !..................................
3750 ! IF THERE IS NO CLOUD/PRECIP WATER, THEN SKIP CALCULATIONS
3751 
3752  IF (qc3d(k).LT.qsmall.AND.qi3d(k).LT.qsmall.AND.qni3d(k).LT.qsmall &
3753  .AND.qr3d(k).LT.qsmall.AND.qg3d(k).LT.qsmall) GOTO 500
3754 
3755 !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
3756 ! CALCULATE INSTANTANEOUS PROCESSES
3757 
3758 ! ADD MELTING OF CLOUD ICE TO FORM RAIN
3759 
3760  IF (qi3d(k).GE.qsmall.AND.t3d(k).GE.273.15) THEN
3761  qr3d(k) = qr3d(k)+qi3d(k)
3762  t3d(k) = t3d(k)-qi3d(k)*xlf(k)/cpm(k)
3763  qi3d(k) = 0.
3764  nr3d(k) = nr3d(k)+ni3d(k)
3765  ni3d(k) = 0.
3766  END IF
3767 ! Fortran version
3768 ! ****SENSITIVITY - NO ICE
3769  IF (iliq.EQ.1) GOTO 778
3770 
3771 ! HOMOGENEOUS FREEZING OF CLOUD WATER
3772 
3773  IF (t3d(k).LE.233.15.AND.qc3d(k).GE.qsmall) THEN
3774  qi3d(k)=qi3d(k)+qc3d(k)
3775  t3d(k)=t3d(k)+qc3d(k)*xlf(k)/cpm(k)
3776  qc3d(k)=0.
3777  ni3d(k)=ni3d(k)+nc3d(k)
3778  nc3d(k)=0.
3779  END IF
3780 ! Fortran version
3781 ! HOMOGENEOUS FREEZING OF RAIN
3782 
3783  IF (igraup.EQ.0) THEN
3784 
3785  IF (t3d(k).LE.233.15.AND.qr3d(k).GE.qsmall) THEN
3786  qg3d(k) = qg3d(k)+qr3d(k)
3787  t3d(k) = t3d(k)+qr3d(k)*xlf(k)/cpm(k)
3788  qr3d(k) = 0.
3789  ng3d(k) = ng3d(k)+ nr3d(k)
3790  nr3d(k) = 0.
3791  END IF
3792 
3793  ELSE IF (igraup.EQ.1) THEN
3794 
3795  IF (t3d(k).LE.233.15.AND.qr3d(k).GE.qsmall) THEN
3796  qni3d(k) = qni3d(k)+qr3d(k)
3797  t3d(k) = t3d(k)+qr3d(k)*xlf(k)/cpm(k)
3798  qr3d(k) = 0.
3799  ns3d(k) = ns3d(k)+nr3d(k)
3800  nr3d(k) = 0.
3801  END IF
3802 
3803  END IF
3804 ! Fortran version
3805  778 CONTINUE
3806 
3807 ! MAKE SURE NUMBER CONCENTRATIONS AREN'T NEGATIVE
3808 
3809  ni3d(k) = max(0.,ni3d(k))
3810  ns3d(k) = max(0.,ns3d(k))
3811  nc3d(k) = max(0.,nc3d(k))
3812  nr3d(k) = max(0.,nr3d(k))
3813  ng3d(k) = max(0.,ng3d(k))
3814 
3815 !......................................................................
3816 ! CLOUD ICE
3817 
3818  IF (qi3d(k).GE.qsmall) THEN
3819  lami(k) = (cons12* &
3820  ni3d(k)/qi3d(k))**(1./di)
3821 
3822 ! CHECK FOR SLOPE
3823 
3824 ! ADJUST VARS
3825 
3826  IF (lami(k).LT.lammini) THEN
3827 
3828  lami(k) = lammini
3829 
3830  n0i(k) = lami(k)**4*qi3d(k)/cons12
3831 
3832  ni3d(k) = n0i(k)/lami(k)
3833  ELSE IF (lami(k).GT.lammaxi) THEN
3834  lami(k) = lammaxi
3835  n0i(k) = lami(k)**4*qi3d(k)/cons12
3836 
3837  ni3d(k) = n0i(k)/lami(k)
3838  END IF
3839  END IF
3840 
3841 !......................................................................
3842 ! RAIN
3843 
3844  IF (qr3d(k).GE.qsmall) THEN
3845  lamr(k) = (pi*rhow*nr3d(k)/qr3d(k))**(1./3.)
3846 
3847 ! CHECK FOR SLOPE
3848 
3849 ! ADJUST VARS
3850 
3851  IF (lamr(k).LT.lamminr) THEN
3852 
3853  lamr(k) = lamminr
3854 
3855  n0rr(k) = lamr(k)**4*qr3d(k)/(pi*rhow)
3856 
3857  nr3d(k) = n0rr(k)/lamr(k)
3858  ELSE IF (lamr(k).GT.lammaxr) THEN
3859  lamr(k) = lammaxr
3860  n0rr(k) = lamr(k)**4*qr3d(k)/(pi*rhow)
3861 
3862  nr3d(k) = n0rr(k)/lamr(k)
3863  END IF
3864 
3865  END IF
3866 
3867 !......................................................................
3868 ! CLOUD DROPLETS
3869 
3870 ! MARTIN ET AL. (1994) FORMULA FOR PGAM
3871 
3872  IF (qc3d(k).GE.qsmall) THEN
3873 
3874  dum = pres(k)/(287.15*t3d(k))
3875  pgam(k)=0.0005714*(nc3d(k)/1.e6*dum)+0.2714
3876  pgam(k)=1./(pgam(k)**2)-1.
3877  pgam(k)=max(pgam(k),2.)
3878  pgam(k)=min(pgam(k),10.)
3879 
3880 ! CALCULATE LAMC
3881 
3882  lamc(k) = (cons26*nc3d(k)*gamma(pgam(k)+4.)/ &
3883  (qc3d(k)*gamma(pgam(k)+1.)))**(1./3.)
3884 
3885 ! LAMMIN, 60 MICRON DIAMETER
3886 ! LAMMAX, 1 MICRON
3887 
3888  lammin = (pgam(k)+1.)/60.e-6
3889  lammax = (pgam(k)+1.)/1.e-6
3890 
3891  IF (lamc(k).LT.lammin) THEN
3892  lamc(k) = lammin
3893  nc3d(k) = exp(3.*log(lamc(k))+log(qc3d(k))+ &
3894  log(gamma(pgam(k)+1.))-log(gamma(pgam(k)+4.)))/cons26
3895 
3896  ELSE IF (lamc(k).GT.lammax) THEN
3897  lamc(k) = lammax
3898  nc3d(k) = exp(3.*log(lamc(k))+log(qc3d(k))+ &
3899  log(gamma(pgam(k)+1.))-log(gamma(pgam(k)+4.)))/cons26
3900 
3901  END IF
3902 
3903  END IF
3904 
3905 !......................................................................
3906 ! SNOW
3907 
3908  IF (qni3d(k).GE.qsmall) THEN
3909  lams(k) = (cons1*ns3d(k)/qni3d(k))**(1./ds)
3910 
3911 ! CHECK FOR SLOPE
3912 
3913 ! ADJUST VARS
3914 
3915  IF (lams(k).LT.lammins) THEN
3916  lams(k) = lammins
3917  n0s(k) = lams(k)**4*qni3d(k)/cons1
3918 
3919  ns3d(k) = n0s(k)/lams(k)
3920 
3921  ELSE IF (lams(k).GT.lammaxs) THEN
3922 
3923  lams(k) = lammaxs
3924  n0s(k) = lams(k)**4*qni3d(k)/cons1
3925  ns3d(k) = n0s(k)/lams(k)
3926  END IF
3927 
3928  END IF
3929 
3930 !......................................................................
3931 ! GRAUPEL
3932 
3933  IF (qg3d(k).GE.qsmall) THEN
3934  lamg(k) = (cons2*ng3d(k)/qg3d(k))**(1./dg)
3935 
3936 ! CHECK FOR SLOPE
3937 
3938 ! ADJUST VARS
3939 
3940  IF (lamg(k).LT.lamming) THEN
3941  lamg(k) = lamming
3942  n0g(k) = lamg(k)**4*qg3d(k)/cons2
3943 
3944  ng3d(k) = n0g(k)/lamg(k)
3945 
3946  ELSE IF (lamg(k).GT.lammaxg) THEN
3947 
3948  lamg(k) = lammaxg
3949  n0g(k) = lamg(k)**4*qg3d(k)/cons2
3950 
3951  ng3d(k) = n0g(k)/lamg(k)
3952  END IF
3953 
3954  END IF
3955 
3956  500 CONTINUE
3957 
3958 ! CALCULATE EFFECTIVE RADIUS
3959 
3960  IF (qi3d(k).GE.qsmall) THEN
3961  effi(k) = 3./lami(k)/2.*1.e6
3962  ELSE
3963  effi(k) = 25.
3964  END IF
3965 
3966  IF (qni3d(k).GE.qsmall) THEN
3967  effs(k) = 3./lams(k)/2.*1.e6
3968  ELSE
3969  effs(k) = 25.
3970  END IF
3971 
3972  IF (qr3d(k).GE.qsmall) THEN
3973  effr(k) = 3./lamr(k)/2.*1.e6
3974  ELSE
3975  effr(k) = 25.
3976  END IF
3977 
3978  IF (qc3d(k).GE.qsmall) THEN
3979  effc(k) = gamma(pgam(k)+4.)/ &
3980  gamma(pgam(k)+3.)/lamc(k)/2.*1.e6
3981  ELSE
3982  effc(k) = 25.
3983  END IF
3984 
3985  IF (qg3d(k).GE.qsmall) THEN
3986  effg(k) = 3./lamg(k)/2.*1.e6
3987  ELSE
3988  effg(k) = 25.
3989  END IF
3990 
3991 ! HM ADD 1/10/06, ADD UPPER BOUND ON ICE NUMBER, THIS IS NEEDED
3992 ! TO PREVENT VERY LARGE ICE NUMBER DUE TO HOMOGENEOUS FREEZING
3993 ! OF DROPLETS, ESPECIALLY WHEN INUM = 1, SET MAX AT 10 CM-3
3994 ! NI3D(K) = MIN(NI3D(K),10.E6/RHO(K))
3995 ! HM, 12/28/12, LOWER MAXIMUM ICE CONCENTRATION TO ADDRESS PROBLEM
3996 ! OF EXCESSIVE AND PERSISTENT ANVIL
3997 ! NOTE: THIS MAY CHANGE/REDUCE SENSITIVITY TO AEROSOL/CCN CONCENTRATION
3998  ni3d(k) = min(ni3d(k),0.3e6/rho(k))
3999 
4000 ! ADD BOUND ON DROPLET NUMBER - CANNOT EXCEED AEROSOL CONCENTRATION
4001  IF (iinum.EQ.0.AND.iact.EQ.2) THEN
4002  nc3d(k) = min(nc3d(k),(nanew1+nanew2)/rho(k))
4003  END IF
4004 ! SWITCH FOR CONSTANT DROPLET NUMBER
4005  IF (iinum.EQ.1) THEN
4006 ! CHANGE NDCNST FROM CM-3 TO KG-1
4007  nc3d(k) = ndcnst*1.e6/rho(k)
4008  END IF
4009 
4010  END DO !!! K LOOP
4011 
4012  400 CONTINUE
4013 
4014 ! ALL DONE !!!!!!!!!!!
4015  RETURN

Referenced by mp_morr_two_moment().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ mp_morr_two_moment()

subroutine, public module_mp_morr_two_moment::mp_morr_two_moment ( integer, intent(in)  ITIMESTEP,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  TH,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  QV,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  QC,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  QR,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  QI,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  QS,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  QG,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  NI,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  NS,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  NR,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  NG,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(in)  RHO,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(in)  PII,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(in)  P,
real(c_double), intent(in)  DT_IN,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(in)  DZ,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(in)  W,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  RAINNC,
real(c_double), dimension(ims:ime, jms:jme), intent(inout)  RAINNCV,
real(c_double), dimension(ims:ime, jms:jme), intent(inout)  SR,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  SNOWNC,
real(c_double), dimension(ims:ime, jms:jme), intent(inout)  SNOWNCV,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout)  GRAUPELNC,
real(c_double), dimension(ims:ime, jms:jme), intent(inout)  GRAUPELNCV,
real(c_double), dimension(ims:ime, jms:kme, kms:jme), intent(inout)  refl_10cm,
logical, intent(in), optional  diagflag,
integer, intent(in), optional  do_radar_ref,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(in)  qrcuten,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(in)  qscuten,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(in)  qicuten,
logical, intent(in), optional  F_QNDROP,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout), optional  qndrop,
integer, intent(in)  IDS,
integer, intent(in)  IDE,
integer, intent(in)  JDS,
integer, intent(in)  JDE,
integer, intent(in)  KDS,
integer, intent(in)  KDE,
integer, intent(in)  IMS,
integer, intent(in)  IME,
integer, intent(in)  JMS,
integer, intent(in)  JME,
integer, intent(in)  KMS,
integer, intent(in)  KME,
integer, intent(in)  ITS,
integer, intent(in)  ITE,
integer, intent(in)  JTS,
integer, intent(in)  JTE,
integer, intent(in)  KTS,
integer, intent(in)  KTE,
logical, intent(in), optional  wetscav_on,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout), optional  rainprod,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout), optional  evapprod,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout), optional  QLSINK,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout), optional  PRECR,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout), optional  PRECI,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout), optional  PRECS,
real(c_double), dimension(ims:ime, jms:jme, kms:kme), intent(inout), optional  PRECG 
)
614 
615 ! QV - water vapor mixing ratio (kg/kg)
616 ! QC - cloud water mixing ratio (kg/kg)
617 ! QR - rain water mixing ratio (kg/kg)
618 ! QI - cloud ice mixing ratio (kg/kg)
619 ! QS - snow mixing ratio (kg/kg)
620 ! QG - graupel mixing ratio (KG/KG)
621 ! NI - cloud ice number concentration (1/kg)
622 ! NS - Snow Number concentration (1/kg)
623 ! NR - Rain Number concentration (1/kg)
624 ! NG - Graupel number concentration (1/kg)
625 ! NOTE: RHO AND HT NOT USED BY THIS SCHEME AND DO NOT NEED TO BE PASSED INTO SCHEME!!!!
626 ! P - AIR PRESSURE (PA)
627 ! W - VERTICAL AIR VELOCITY (M/S)
628 ! TH - POTENTIAL TEMPERATURE (K)
629 ! PII - exner function - used to convert potential temp to temp
630 ! DZ - difference in height over interface (m)
631 ! DT_IN - model time step (sec)
632 ! ITIMESTEP - time step counter
633 ! RAINNC - accumulated grid-scale precipitation (mm)
634 ! RAINNCV - one time step grid scale precipitation (mm/time step)
635 ! SNOWNC - accumulated grid-scale snow plus cloud ice (mm)
636 ! SNOWNCV - one time step grid scale snow plus cloud ice (mm/time step)
637 ! GRAUPELNC - accumulated grid-scale graupel (mm)
638 ! GRAUPELNCV - one time step grid scale graupel (mm/time step)
639 ! SR - one time step mass ratio of snow to total precip
640 ! qrcuten, rain tendency from parameterized cumulus convection
641 ! qscuten, snow tendency from parameterized cumulus convection
642 ! qicuten, cloud ice tendency from parameterized cumulus convection
643 
644 ! variables below currently not in use, not coupled to PBL or radiation codes
645 ! TKE - turbulence kinetic energy (m^2 s-2), NEEDED FOR DROPLET ACTIVATION (SEE CODE BELOW)
646 ! NCTEND - droplet concentration tendency from pbl (kg-1 s-1)
647 ! NCTEND - CLOUD ICE concentration tendency from pbl (kg-1 s-1)
648 ! KZH - heat eddy diffusion coefficient from YSU scheme (M^2 S-1), NEEDED FOR DROPLET ACTIVATION (SEE CODE BELOW)
649 ! EFFCS - CLOUD DROPLET EFFECTIVE RADIUS OUTPUT TO RADIATION CODE (micron)
650 ! EFFIS - CLOUD DROPLET EFFECTIVE RADIUS OUTPUT TO RADIATION CODE (micron)
651 ! HM, ADDED FOR WRF-CHEM COUPLING
652 ! QLSINK - TENDENCY OF CLOUD WATER TO RAIN, SNOW, GRAUPEL (KG/KG/S)
653 ! CSED,ISED,SSED,GSED,RSED - SEDIMENTATION FLUXES (KG/M^2/S) FOR CLOUD WATER, ICE, SNOW, GRAUPEL, RAIN
654 ! PRECI,PRECS,PRECG,PRECR - SEDIMENTATION FLUXES (KG/M^2/S) FOR ICE, SNOW, GRAUPEL, RAIN
655 
656 ! rainprod - total tendency of conversion of cloud water/ice and graupel to rain (kg kg-1 s-1)
657 ! evapprod - tendency of evaporation of rain (kg kg-1 s-1)
658 
659 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
660 
661 ! reflectivity currently not included!!!!
662 ! REFL_10CM - CALCULATED RADAR REFLECTIVITY AT 10 CM (DBZ)
663 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
664 
665 ! EFFC - DROPLET EFFECTIVE RADIUS (MICRON)
666 ! EFFR - RAIN EFFECTIVE RADIUS (MICRON)
667 ! EFFS - SNOW EFFECTIVE RADIUS (MICRON)
668 ! EFFI - CLOUD ICE EFFECTIVE RADIUS (MICRON)
669 
670 ! ADDITIONAL OUTPUT FROM MICRO - SEDIMENTATION TENDENCIES, NEEDED FOR LIQUID-ICE STATIC ENERGY
671 
672 ! QGSTEN - GRAUPEL SEDIMENTATION TEND (KG/KG/S)
673 ! QRSTEN - RAIN SEDIMENTATION TEND (KG/KG/S)
674 ! QISTEN - CLOUD ICE SEDIMENTATION TEND (KG/KG/S)
675 ! QNISTEN - SNOW SEDIMENTATION TEND (KG/KG/S)
676 ! QCSTEN - CLOUD WATER SEDIMENTATION TEND (KG/KG/S)
677 
678  IMPLICIT NONE
679 
680  INTEGER, INTENT(IN ) :: ids, ide, jds, jde, kds, kde , &
681  ims, ime, jms, jme, kms, kme , &
682  its, ite, jts, jte, kts, kte
683 
684  REAL(C_DOUBLE), DIMENSION(ims:ime, jms:jme, kms:kme), INTENT(INOUT) :: qv, qc, qr, qi, qs, qg, ni, ns, nr, TH, NG
685  REAL(C_DOUBLE), DIMENSION(ims:ime, jms:jme, kms:kme), optional,INTENT(INOUT) :: qndrop
686  REAL(C_DOUBLE), DIMENSION(ims:ime, jms:jme, kms:kme), optional,INTENT(INOUT) :: QLSINK, rainprod, evapprod, PRECI,PRECS,PRECG, &
687  & PRECR
688  REAL(C_DOUBLE), DIMENSION(ims:ime, jms:jme, kms:kme), INTENT(IN ) :: pii, p, dz, rho, w
689 
690  REAL(C_DOUBLE), INTENT(IN):: dt_in
691  INTEGER, INTENT(IN):: ITIMESTEP
692  REAL(C_DOUBLE), INTENT(INOUT), DIMENSION(ims:ime, jms:jme, kms:kme) :: rainnc, snownc, graupelnc
693  REAL(C_DOUBLE), INTENT(INOUT), DIMENSION(ims:ime, jms:jme) :: rainncv, sr, snowncv, graupelncv
694 ! REAL(C_DOUBLE), DIMENSION(ims:ime, jms:jme ), INTENT(INOUT) :: RAINNC, RAINNCV, SR, SNOWNC,SNOWNCV,GRAUPELNC,GRAUPELNCV
695  REAL(C_DOUBLE), DIMENSION(ims:ime, jms:kme, kms:jme), INTENT(INOUT) :: refl_10cm
696 ! REAL(C_DOUBLE), DIMENSION(ims:ime ,jms:jme ), INTENT(IN ) :: ht
697 
698  LOGICAL, optional, INTENT(IN) :: wetscav_on
699 
700  ! LOCAL VARIABLES
701 
702  REAL(C_DOUBLE), DIMENSION(its:ite, jts:jte, kts:kte) :: effi, effs, effr, EFFG
703  REAL(C_DOUBLE), DIMENSION(its:ite, jts:jte, kts:kte) :: T, EFFC
704 
705  REAL(C_DOUBLE), DIMENSION(kts:kte) :: &
706  QC_TEND1D, QI_TEND1D, QNI_TEND1D, QR_TEND1D, &
707  NI_TEND1D, NS_TEND1D, NR_TEND1D, &
708  QC1D, QI1D, QR1D,NI1D, NS1D, NR1D, QS1D, &
709  T_TEND1D,QV_TEND1D, T1D, QV1D, P1D, W1D, &
710  EFFC1D, EFFI1D, EFFS1D, EFFR1D,DZ1D, &
711  QG_TEND1D, NG_TEND1D, QG1D, NG1D, EFFG1D, &
712 ! ADD SEDIMENTATION TENDENCIES (UNITS OF KG/KG/S)
713  qgsten,qrsten, qisten, qnisten, qcsten, &
714 ! ADD CUMULUS TENDENCIES
715  qrcu1d, qscu1d, qicu1d
716 ! add cumulus tendencies
717 
718  REAL(C_DOUBLE), DIMENSION(ims:ime, jms:jme, kms:kme), INTENT(IN):: qrcuten, qscuten, qicuten
719 
720  LOGICAL, INTENT(IN), OPTIONAL :: F_QNDROP
721  LOGICAL :: flag_qndrop ! wrf-chem
722  integer :: iinum ! wrf-chem
723 
724 ! wrf-chem
725  REAL(C_DOUBLE), DIMENSION(kts:kte) :: nc1d, nc_tend1d,C2PREC,CSED,ISED,SSED,GSED,RSED
726  REAL(C_DOUBLE), DIMENSION(kts:kte) :: rainprod1d, evapprod1d
727 ! HM add reflectivity
728  REAL(C_DOUBLE), DIMENSION(kts:kte) :: dBZ
729 
730  REAL(C_DOUBLE) PRECPRT1D, SNOWRT1D, SNOWPRT1D, GRPLPRT1D ! hm added 7/13/13
731 
732  INTEGER I,J,K
733 
734  REAL(C_DOUBLE) DT
735 
736  LOGICAL, OPTIONAL, INTENT(IN) :: diagflag
737  INTEGER, OPTIONAL, INTENT(IN) :: do_radar_ref
738 
739  LOGICAL :: has_wetscav
740 
741  ! return
742  print *,'IN MORR_TWO_MOMENT MAIN ROUTINE '
743 
744 ! below for wrf-chem
745  flag_qndrop = .false.
746  IF ( PRESENT ( f_qndrop ) ) flag_qndrop = f_qndrop
747 !!!!!!!!!!!!!!!!!!!!!!
748 
749  IF( PRESENT( wetscav_on ) ) THEN
750  has_wetscav = wetscav_on
751  ELSE
752  has_wetscav = .false.
753  ENDIF
754 
755  ! Initialize tendencies (all set to 0) and transfer
756  ! array to local variables
757  dt = dt_in
758 
759  ! Currently mixing of number concentrations also is neglected (not coupled with PBL schemes)
760 
761  ! print *,'IMS IME ',ims,ime
762  ! print *,'JMS JME ',jms,jme
763  ! print *,'KMS KME ',kms,kme
764 
765  ! print *,'ITS ITE ',its,ite
766  ! print *,'JTS JTE ',jts,jte
767  ! print *,'KTS KTE ',kts,kte
768 
769  DO i=its,ite
770  DO j=jts,jte
771  DO k=kts,kte
772  t(i,j,k) = th(i,j,k) * pii(i,j,k)
773  END DO
774  END DO
775  END DO
776 
777  do i=its,ite ! i loop (east-west)
778  do j=jts,jte ! j loop (north-south)
779  !
780  ! Transfer 3D arrays into 1D for microphysical calculations
781  !
782 
783 ! hm , initialize 1d tendency arrays to zero
784 
785  do k=kts,kte ! k loop (vertical)
786 
787  qc_tend1d(k) = 0.
788  qi_tend1d(k) = 0.
789  qni_tend1d(k) = 0.
790  qr_tend1d(k) = 0.
791  ni_tend1d(k) = 0.
792  ns_tend1d(k) = 0.
793  nr_tend1d(k) = 0.
794  t_tend1d(k) = 0.
795  qv_tend1d(k) = 0.
796  nc_tend1d(k) = 0. ! wrf-chem
797 
798  qc1d(k) = qc(i,j,k)
799  qi1d(k) = qi(i,j,k)
800  qs1d(k) = qs(i,j,k)
801  qr1d(k) = qr(i,j,k)
802 
803  ni1d(k) = ni(i,j,k)
804 
805  ns1d(k) = ns(i,j,k)
806  nr1d(k) = nr(i,j,k)
807 ! HM ADD GRAUPEL
808  qg1d(k) = qg(i,j,k)
809  ng1d(k) = ng(i,j,k)
810  qg_tend1d(k) = 0.
811  ng_tend1d(k) = 0.
812 
813  t1d(k) = t(i,j,k)
814  qv1d(k) = qv(i,j,k)
815  p1d(k) = p(i,j,k)
816  dz1d(k) = dz(i,j,k)
817  w1d(k) = w(i,j,k)
818 ! add cumulus tendencies, already decoupled
819  qrcu1d(k) = qrcuten(i,j,k)
820  qscu1d(k) = qscuten(i,j,k)
821  qicu1d(k) = qicuten(i,j,k)
822  end do !jdf added this
823 
824 ! below for wrf-chem
825  IF (flag_qndrop .AND. PRESENT( qndrop )) THEN
826  iact = 3
827  DO k = kts, kte
828  nc1d(k)=qndrop(i,j,k)
829  iinum=0
830  ENDDO
831  ELSE
832  DO k = kts, kte
833  nc1d(k)=0. ! temporary placeholder, set to constant in microphysics subroutine
834  iinum=1
835  ENDDO
836  ENDIF
837 
838 !jdf end do
839 
840  call morr_two_moment_micro(i,j,itimestep,kts,kte,qc_tend1d, qi_tend1d, qni_tend1d, qr_tend1d, &
841  ni_tend1d, ns_tend1d, nr_tend1d, &
842  qc1d, qi1d, qs1d, qr1d,ni1d, ns1d, nr1d, &
843  t_tend1d,qv_tend1d, t1d, qv1d, p1d, dz1d, w1d, &
844  precprt1d,snowrt1d, &
845  snowprt1d,grplprt1d, & ! hm added 7/13/13
846  effc1d,effi1d,effs1d,effr1d,dt, &
847  qg_tend1d,ng_tend1d,qg1d,ng1d,effg1d, &
848  qrcu1d, qscu1d, qicu1d, &
849 ! ADD SEDIMENTATION TENDENCIES
850  qgsten,qrsten,qisten,qnisten,qcsten, &
851  nc1d, nc_tend1d, iinum, c2prec,csed,ised,ssed,gsed,rsed)
852 
853  !
854  ! Transfer 1D arrays back into 3D arrays
855  !
856  do k=kts,kte
857 
858  ! hm, add tendencies to update global variables
859  ! HM, TENDENCIES FOR Q AND N NOW ADDED IN M2005MICRO, SO WE
860  ! ONLY NEED TO TRANSFER 1D VARIABLES BACK TO 3D
861 
862  qc(i,j,k) = qc1d(k)
863  qi(i,j,k) = qi1d(k)
864  qs(i,j,k) = qs1d(k)
865  qr(i,j,k) = qr1d(k)
866  ni(i,j,k) = ni1d(k)
867  ns(i,j,k) = ns1d(k)
868  nr(i,j,k) = nr1d(k)
869  qg(i,j,k) = qg1d(k)
870  ng(i,j,k) = ng1d(k)
871 
872  t(i,j,k) = t1d(k)
873  th(i,j,k) = t(i,j,k)/pii(i,j,k) ! CONVERT TEMP BACK TO POTENTIAL TEMP
874  qv(i,j,k) = qv1d(k)
875 
876  effc(i,j,k) = effc1d(k)
877  effi(i,j,k) = effi1d(k)
878  effs(i,j,k) = effs1d(k)
879  effr(i,j,k) = effr1d(k)
880  effg(i,j,k) = effg1d(k)
881 
882 ! wrf-chem
883  IF (flag_qndrop .AND. PRESENT( qndrop )) THEN
884  qndrop(i,j,k) = nc1d(k)
885 !jdf CSED3D(i,j,k) = CSED(k)
886  END IF
887  IF ( PRESENT( qlsink ) ) THEN
888  if(qc(i,j,k)>1.e-10) then
889  qlsink(i,j,k) = c2prec(k)/qc(i,j,k)
890  else
891  qlsink(i,j,k) = 0.0
892  endif
893  END IF
894  IF ( PRESENT( precr ) ) precr(i,j,k) = rsed(k)
895  IF ( PRESENT( preci ) ) preci(i,j,k) = ised(k)
896  IF ( PRESENT( precs ) ) precs(i,j,k) = ssed(k)
897  IF ( PRESENT( precg ) ) precg(i,j,k) = gsed(k)
898 ! EFFECTIVE RADIUS FOR RADIATION CODE (currently not coupled)
899 ! HM, ADD LIMIT TO PREVENT BLOWING UP OPTICAL PROPERTIES, 8/18/07
900 ! EFFCS(I,J,K) = MIN(EFFC(I,J,K),50.)
901 ! EFFCS(I,J,K) = MAX(EFFCS(I,J,K),1.)
902 ! EFFIS(I,J,K) = MIN(EFFI(I,J,K),130.)
903 ! EFFIS(I,J,K) = MAX(EFFIS(I,J,K),13.)
904 
905  end do
906 
907 ! hm modified so that m2005 precip variables correctly match wrf precip variables
908  rainnc(i,j,kts) = rainnc(i,j,kts)+precprt1d
909  rainncv(i,j) = precprt1d
910  snownc(i,j,kts) = snownc(i,j,kts)+snowprt1d
911  snowncv(i,j) = snowprt1d
912  graupelnc(i,j,kts) = graupelnc(i,j,kts)+grplprt1d
913  graupelncv(i,j) = grplprt1d
914  sr(i,j) = snowrt1d/(precprt1d+1.e-12)
915 
916  end do
917  end do
918 

Referenced by mp_morr_two_moment_isohelper::mp_morr_two_moment_c().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ polysvp()

real(c_double) function, public module_mp_morr_two_moment::polysvp ( real(c_double)  T,
integer  TYPE 
)
4021 
4022 !-------------------------------------------
4023 
4024 ! COMPUTE SATURATION VAPOR PRESSURE
4025 
4026 ! POLYSVP RETURNED IN UNITS OF PA.
4027 ! T IS INPUT IN UNITS OF K.
4028 ! TYPE REFERS TO SATURATION WITH RESPECT TO LIQUID (0) OR ICE (1)
4029 
4030 ! REPLACE GOFF-GRATCH WITH FASTER FORMULATION FROM FLATAU ET AL. 1992, TABLE 4 (RIGHT-HAND COLUMN)
4031 
4032  IMPLICIT NONE
4033 
4034  REAL(C_DOUBLE) DUM
4035  REAL(C_DOUBLE) T
4036  INTEGER TYPE
4037 ! ice
4038  real a0i,a1i,a2i,a3i,a4i,a5i,a6i,a7i,a8i
4039  data a0i,a1i,a2i,a3i,a4i,a5i,a6i,a7i,a8i /&
4040  6.11147274, 0.503160820, 0.188439774e-1, &
4041  0.420895665e-3, 0.615021634e-5,0.602588177e-7, &
4042  0.385852041e-9, 0.146898966e-11, 0.252751365e-14/
4043 
4044 ! liquid
4045  real a0,a1,a2,a3,a4,a5,a6,a7,a8
4046 
4047 ! V1.7
4048  data a0,a1,a2,a3,a4,a5,a6,a7,a8 /&
4049  6.11239921, 0.443987641, 0.142986287e-1, &
4050  0.264847430e-3, 0.302950461e-5, 0.206739458e-7, &
4051  0.640689451e-10,-0.952447341e-13,-0.976195544e-15/
4052  real dt
4053 
4054 ! ICE
4055 
4056  IF (type.EQ.1) THEN
4057 
4058 ! POLYSVP = 10.**(-9.09718*(273.16/T-1.)-3.56654* &
4059 ! LOG10(273.16/T)+0.876793*(1.-T/273.16)+ &
4060 ! LOG10(6.1071))*100.
4061 ! hm 11/16/20, use Goff-Gratch for T < 195.8 K and Flatau et al. equal or above 195.8 K
4062  if (t.ge.195.8) then
4063  dt=t-273.15
4064  polysvp = a0i + dt*(a1i+dt*(a2i+dt*(a3i+dt*(a4i+dt*(a5i+dt*(a6i+dt*(a7i+a8i*dt)))))))
4065  polysvp = polysvp*100.
4066  else
4067 
4068  polysvp = 10.**(-9.09718*(273.16/t-1.)-3.56654* &
4069  dlog10(273.16d0/t)+0.876793*(1.-t/273.16)+ &
4070  dlog10(6.1071d0))*100.
4071 
4072  end if
4073 
4074  END IF
4075 
4076 ! LIQUID
4077 
4078  IF (type.EQ.0) THEN
4079 
4080 ! POLYSVP = 10.**(-7.90298*(373.16/T-1.)+ &
4081 ! 5.02808*LOG10(373.16/T)- &
4082 ! 1.3816E-7*(10**(11.344*(1.-T/373.16))-1.)+ &
4083 ! 8.1328E-3*(10**(-3.49149*(373.16/T-1.))-1.)+ &
4084 ! LOG10(1013.246))*100.
4085 ! hm 11/16/20, use Goff-Gratch for T < 202.0 K and Flatau et al. equal or above 202.0 K
4086  if (t.ge.202.0) then
4087  dt = t-273.15
4088  polysvp = a0 + dt*(a1+dt*(a2+dt*(a3+dt*(a4+dt*(a5+dt*(a6+dt*(a7+a8*dt)))))))
4089  polysvp = polysvp*100.
4090  else
4091 
4092 ! note: uncertain below -70 C, but produces physical values (non-negative) unlike flatau
4093  polysvp = 10.**(-7.90298*(373.16/t-1.)+ &
4094  5.02808*dlog10(373.16d0/t)- &
4095  1.3816e-7*(10**(11.344*(1.-t/373.16))-1.)+ &
4096  8.1328e-3*(10**(-3.49149*(373.16/t-1.))-1.)+ &
4097  dlog10(1013.246d0))*100.
4098  end if
4099 
4100  END IF
4101 
4102 

Referenced by calc_saturation_vapor_pressure(), and morr_two_moment_micro().

Here is the caller graph for this function:

◆ set_morrison_ndcnst()

subroutine, public module_mp_morr_two_moment::set_morrison_ndcnst ( real(c_double), intent(in), value  ndcnst_in)
574 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
575 ! THIS SUBROUTINE SETS THE CONSTANT DROPLET CONCENTRATION (NDCNST)
576 ! ALLOWS RUNTIME CONFIGURATION FROM INPUT FILE
577 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
578  IMPLICIT NONE
579  REAL(C_DOUBLE), INTENT(IN), VALUE :: ndcnst_in
580 
581  ndcnst = ndcnst_in
582 

Variable Documentation

◆ ac

real(c_double), private module_mp_morr_two_moment::ac
private

Referenced by erf_slow_rhs_post(), and morr_two_moment_init().

◆ ag

real(c_double), private module_mp_morr_two_moment::ag
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ ai

real(c_double), private module_mp_morr_two_moment::ai
private
181  REAL(C_DOUBLE), PRIVATE :: AI,AC,AS,AR,AG ! 'A' PARAMETER IN FALLSPEED-DIAM RELATIONSHIP

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ aimm

real(c_double), private module_mp_morr_two_moment::aimm
private
191  REAL(C_DOUBLE), PRIVATE :: AIMM ! PARAMETER IN BIGG IMMERSION FREEZING

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ ar

real(c_double), private module_mp_morr_two_moment::ar
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ as

real(c_double), private module_mp_morr_two_moment::as
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ bact

real(c_double), private module_mp_morr_two_moment::bact
private
225  REAL(C_DOUBLE), PRIVATE :: BACT ! ACTIVATION PARAMETER

Referenced by morr_two_moment_init().

◆ bc

real(c_double), private module_mp_morr_two_moment::bc
private

Referenced by morr_two_moment_init(), morr_two_moment_micro(), and redistribute_term().

◆ bg

real(c_double), private module_mp_morr_two_moment::bg
private

Referenced by AddPertToBckgnd(), morr_two_moment_init(), and morr_two_moment_micro().

◆ bi

real(c_double), private module_mp_morr_two_moment::bi
private
182  REAL(C_DOUBLE), PRIVATE :: BI,BC,BS,BR,BG ! 'B' PARAMETER IN FALLSPEED-DIAM RELATIONSHIP

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ bimm

real(c_double), private module_mp_morr_two_moment::bimm
private
192  REAL(C_DOUBLE), PRIVATE :: BIMM ! PARAMETER IN BIGG IMMERSION FREEZING

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ br

real(c_double), private module_mp_morr_two_moment::br
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ bs

real(c_double), private module_mp_morr_two_moment::bs
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ c1

real(c_double), private module_mp_morr_two_moment::c1
private
212  REAL(C_DOUBLE), PRIVATE :: C1 ! 'C' IN NCCN = CS^K (CM-3)

Referenced by SatMethods::Bolton_svp_water(), TerminalVelocity< RT >::CloudRainShima(), ComputeStrain_S(), ComputeStrain_T(), DiffusionSrcForState_S(), DiffusionSrcForState_T(), interp_trilinear(), and morr_two_moment_init().

◆ cg

real(c_double), private module_mp_morr_two_moment::cg
private

Referenced by morr_two_moment_init().

◆ ci

real(c_double), private module_mp_morr_two_moment::ci
private
203  REAL(C_DOUBLE), PRIVATE :: CI,DI,CS,DS,CG,DG ! SIZE DISTRIBUTION PARAMETERS FOR CLOUD ICE, SNOW, GRAUPEL

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons1

real(c_double), private module_mp_morr_two_moment::cons1
private
241  REAL(C_DOUBLE), PRIVATE :: CONS1,CONS2,CONS3,CONS4,CONS5,CONS6,CONS7,CONS8,CONS9,CONS10

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons10

real(c_double), private module_mp_morr_two_moment::cons10
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons11

real(c_double), private module_mp_morr_two_moment::cons11
private
242  REAL(C_DOUBLE), PRIVATE :: CONS11,CONS12,CONS13,CONS14,CONS15,CONS16,CONS17,CONS18,CONS19,CONS20

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons12

real(c_double), private module_mp_morr_two_moment::cons12
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons13

real(c_double), private module_mp_morr_two_moment::cons13
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons14

real(c_double), private module_mp_morr_two_moment::cons14
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons15

real(c_double), private module_mp_morr_two_moment::cons15
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons16

real(c_double), private module_mp_morr_two_moment::cons16
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons17

real(c_double), private module_mp_morr_two_moment::cons17
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons18

real(c_double), private module_mp_morr_two_moment::cons18
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons19

real(c_double), private module_mp_morr_two_moment::cons19
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons2

real(c_double), private module_mp_morr_two_moment::cons2
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons20

real(c_double), private module_mp_morr_two_moment::cons20
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons21

real(c_double), private module_mp_morr_two_moment::cons21
private
243  REAL(C_DOUBLE), PRIVATE :: CONS21,CONS22,CONS23,CONS24,CONS25,CONS26,CONS27,CONS28,CONS29,CONS30

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons22

real(c_double), private module_mp_morr_two_moment::cons22
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons23

real(c_double), private module_mp_morr_two_moment::cons23
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons24

real(c_double), private module_mp_morr_two_moment::cons24
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons25

real(c_double), private module_mp_morr_two_moment::cons25
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons26

real(c_double), private module_mp_morr_two_moment::cons26
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons27

real(c_double), private module_mp_morr_two_moment::cons27
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons28

real(c_double), private module_mp_morr_two_moment::cons28
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons29

real(c_double), private module_mp_morr_two_moment::cons29
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons3

real(c_double), private module_mp_morr_two_moment::cons3
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons30

real(c_double), private module_mp_morr_two_moment::cons30
private

Referenced by morr_two_moment_init().

◆ cons31

real(c_double), private module_mp_morr_two_moment::cons31
private
244  REAL(C_DOUBLE), PRIVATE :: CONS31,CONS32,CONS33,CONS34,CONS35,CONS36,CONS37,CONS38,CONS39,CONS40

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons32

real(c_double), private module_mp_morr_two_moment::cons32
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons33

real(c_double), private module_mp_morr_two_moment::cons33
private

Referenced by morr_two_moment_init().

◆ cons34

real(c_double), private module_mp_morr_two_moment::cons34
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons35

real(c_double), private module_mp_morr_two_moment::cons35
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons36

real(c_double), private module_mp_morr_two_moment::cons36
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons37

real(c_double), private module_mp_morr_two_moment::cons37
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons38

real(c_double), private module_mp_morr_two_moment::cons38
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons39

real(c_double), private module_mp_morr_two_moment::cons39
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons4

real(c_double), private module_mp_morr_two_moment::cons4
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons40

real(c_double), private module_mp_morr_two_moment::cons40
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons41

real(c_double), private module_mp_morr_two_moment::cons41
private
245  REAL(C_DOUBLE), PRIVATE :: CONS41

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons5

real(c_double), private module_mp_morr_two_moment::cons5
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons6

real(c_double), private module_mp_morr_two_moment::cons6
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons7

real(c_double), private module_mp_morr_two_moment::cons7
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons8

real(c_double), private module_mp_morr_two_moment::cons8
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cons9

real(c_double), private module_mp_morr_two_moment::cons9
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cpw

real(c_double), private module_mp_morr_two_moment::cpw
private
208  REAL(C_DOUBLE), PRIVATE :: CPW ! SPECIFIC HEAT OF LIQUID WATER

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ cs

real(c_double), private module_mp_morr_two_moment::cs
private

Referenced by morr_two_moment_init(), Radiation::run_impl(), and ERF::WriteSubvolume().

◆ dcs

real(c_double), private module_mp_morr_two_moment::dcs
private
194  REAL(C_DOUBLE), PRIVATE :: DCS ! THRESHOLD SIZE FOR CLOUD ICE AUTOCONVERSION

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ dg

real(c_double), private module_mp_morr_two_moment::dg
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ di

real(c_double), private module_mp_morr_two_moment::di
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ ds

real(c_double), private module_mp_morr_two_moment::ds
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ eci

real(c_double), private module_mp_morr_two_moment::eci
private
205  REAL(C_DOUBLE), PRIVATE :: ECI ! COLLECTION EFFICIENCY, ICE-DROPLET COLLISIONS

Referenced by morr_two_moment_init().

◆ ecr

real(c_double), private module_mp_morr_two_moment::ecr
private
193  REAL(C_DOUBLE), PRIVATE :: ECR ! COLLECTION EFFICIENCY BETWEEN DROPLETS/RAIN AND SNOW/RAIN

Referenced by morr_two_moment_init().

◆ eii

real(c_double), private module_mp_morr_two_moment::eii
private
204  REAL(C_DOUBLE), PRIVATE :: EII ! COLLECTION EFFICIENCY, ICE-ICE COLLISIONS

Referenced by morr_two_moment_init().

◆ epsm

real(c_double), private module_mp_morr_two_moment::epsm
private
220  REAL(C_DOUBLE), PRIVATE :: EPSM ! AEROSOL SOLUBLE FRACTION

Referenced by morr_two_moment_init().

◆ f11

real(c_double), private module_mp_morr_two_moment::f11
private
232  REAL(C_DOUBLE), PRIVATE :: F11 ! CORRECTION FACTOR FOR ACTIVATION, MODE 1

Referenced by morr_two_moment_init().

◆ f12

real(c_double), private module_mp_morr_two_moment::f12
private
233  REAL(C_DOUBLE), PRIVATE :: F12 ! CORRECTION FACTOR FOR ACTIVATION, MODE 1

Referenced by morr_two_moment_init().

◆ f1r

real(c_double), private module_mp_morr_two_moment::f1r
private
199  REAL(C_DOUBLE), PRIVATE :: F1R ! VENTILATION PARAMETER FOR RAIN

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ f1s

real(c_double), private module_mp_morr_two_moment::f1s
private
197  REAL(C_DOUBLE), PRIVATE :: F1S ! VENTILATION PARAMETER FOR SNOW

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ f21

real(c_double), private module_mp_morr_two_moment::f21
private
234  REAL(C_DOUBLE), PRIVATE :: F21 ! CORRECTION FACTOR FOR ACTIVATION, MODE 2

Referenced by morr_two_moment_init().

◆ f22

real(c_double), private module_mp_morr_two_moment::f22
private
235  REAL(C_DOUBLE), PRIVATE :: F22 ! CORRECTION FACTOR FOR ACTIVATION, MODE 2

Referenced by morr_two_moment_init().

◆ f2r

real(c_double), private module_mp_morr_two_moment::f2r
private
200  REAL(C_DOUBLE), PRIVATE :: F2R ! VENTILATION PARAMETER FOR RAIN

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ f2s

real(c_double), private module_mp_morr_two_moment::f2s
private
198  REAL(C_DOUBLE), PRIVATE :: F2S ! VENTILATION PARAMETER FOR SNOW

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ iact

integer, private module_mp_morr_two_moment::iact
private
118  INTEGER, PRIVATE :: IACT

Referenced by morr_two_moment_init(), morr_two_moment_micro(), and mp_morr_two_moment().

◆ ibase

integer, private module_mp_morr_two_moment::ibase
private
156  INTEGER, PRIVATE :: IBASE

Referenced by morr_two_moment_init().

◆ igraup

integer, private module_mp_morr_two_moment::igraup
private
170  INTEGER, PRIVATE :: IGRAUP

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ ihail

integer, private module_mp_morr_two_moment::ihail
private
177  INTEGER, PRIVATE :: IHAIL

Referenced by morr_two_moment_init().

◆ iliq

integer, private module_mp_morr_two_moment::iliq
private
134  INTEGER, PRIVATE :: ILIQ

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ inuc

integer, private module_mp_morr_two_moment::inuc
private
140  INTEGER, PRIVATE :: INUC

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ inum

integer, private module_mp_morr_two_moment::inum
private
125  INTEGER, PRIVATE :: INUM

Referenced by morr_two_moment_init().

◆ isub

integer, private module_mp_morr_two_moment::isub
private
164  INTEGER, PRIVATE :: ISUB

Referenced by ERF::ErrorEst(), TurbulentPerturbation::init_tpi(), morr_two_moment_init(), ERF::project_momenta(), and ERF::WriteSubvolume().

◆ k1

real(c_double), private module_mp_morr_two_moment::k1
private
213  REAL(C_DOUBLE), PRIVATE :: K1 ! 'K' IN NCCN = CS^K

Referenced by interp_trilinear(), morr_two_moment_init(), TerrainIF::operator()(), TerminalVelocity< RT >::RogersYau(), and writeNCPlotFile().

◆ lammaxg

real(c_double), private module_mp_morr_two_moment::lammaxg
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ lammaxi

real(c_double), private module_mp_morr_two_moment::lammaxi
private
237  REAL(C_DOUBLE), PRIVATE :: LAMMAXI,LAMMINI,LAMMAXR,LAMMINR,LAMMAXS,LAMMINS,LAMMAXG,LAMMING

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ lammaxr

real(c_double), private module_mp_morr_two_moment::lammaxr
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ lammaxs

real(c_double), private module_mp_morr_two_moment::lammaxs
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ lamming

real(c_double), private module_mp_morr_two_moment::lamming
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ lammini

real(c_double), private module_mp_morr_two_moment::lammini
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ lamminr

real(c_double), private module_mp_morr_two_moment::lamminr
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ lammins

real(c_double), private module_mp_morr_two_moment::lammins
private

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ ma

real(c_double), private module_mp_morr_two_moment::ma
private
223  REAL(C_DOUBLE), PRIVATE :: MA ! MOLECULAR WEIGHT OF 'AIR' (KG/MOL)

Referenced by morr_two_moment_init(), and reduce_to_max_per_height().

◆ map

real(c_double), private module_mp_morr_two_moment::map
private
222  REAL(C_DOUBLE), PRIVATE :: MAP ! MOLECULAR WEIGHT AEROSOL (KG/MOL)

Referenced by morr_two_moment_init().

◆ mg0

real(c_double), private module_mp_morr_two_moment::mg0
private
196  REAL(C_DOUBLE), PRIVATE :: MG0 ! MASS OF EMBRYO GRAUPEL

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ mi0

real(c_double), private module_mp_morr_two_moment::mi0
private
195  REAL(C_DOUBLE), PRIVATE :: MI0 ! INITIAL SIZE OF NUCLEATED CRYSTAL

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ mmult

real(c_double), private module_mp_morr_two_moment::mmult
private
236  REAL(C_DOUBLE), PRIVATE :: MMULT ! MASS OF SPLINTERED ICE PARTICLE

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ mw

real(c_double), private module_mp_morr_two_moment::mw
private
217  REAL(C_DOUBLE), PRIVATE :: MW ! MOLECULAR WEIGHT WATER (KG/MOL)

Referenced by morr_two_moment_init().

◆ nanew1

real(c_double), private module_mp_morr_two_moment::nanew1
private
228  REAL(C_DOUBLE), PRIVATE :: NANEW1 ! TOTAL AEROSOL CONCENTRATION, MODE 1 (M^-3)

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ nanew2

real(c_double), private module_mp_morr_two_moment::nanew2
private
229  REAL(C_DOUBLE), PRIVATE :: NANEW2 ! TOTAL AEROSOL CONCENTRATION, MODE 2 (M^-3)

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ ndcnst

real(c_double), private module_mp_morr_two_moment::ndcnst
private
128  REAL(C_DOUBLE), PRIVATE :: NDCNST

Referenced by morr_two_moment_init(), morr_two_moment_micro(), and set_morrison_ndcnst().

◆ osm

real(c_double), private module_mp_morr_two_moment::osm
private
218  REAL(C_DOUBLE), PRIVATE :: OSM ! OSMOTIC COEFFICIENT

Referenced by morr_two_moment_init().

◆ pi

real(c_double), parameter, private module_mp_morr_two_moment::pi = 3.1415926535897932384626434D0
private
100  REAL(C_DOUBLE), PARAMETER :: PI = 3.1415926535897932384626434d0

Referenced by WSM6::Advance(), compute_max_reflectivity_dbz(), ERF::erf_enforce_hse(), gamma(), morr_two_moment_init(), morr_two_moment_micro(), wrf_gamma(), wsm6_nislfv_rain_plm6_scratch(), and wsm6_nislfv_rain_plm_scratch().

◆ qsmall

real(c_double), private module_mp_morr_two_moment::qsmall
private
202  REAL(C_DOUBLE), PRIVATE :: QSMALL ! SMALLEST ALLOWED HYDROMETEOR MIXING RATIO

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ rhoa

real(c_double), private module_mp_morr_two_moment::rhoa
private
221  REAL(C_DOUBLE), PRIVATE :: RHOA ! AEROSOL BULK DENSITY (KG/M3)

Referenced by morr_two_moment_init().

◆ rhog

real(c_double), private module_mp_morr_two_moment::rhog
private
190  REAL(C_DOUBLE), PRIVATE :: RHOG ! BULK DENSITY OF GRAUPEL

◆ rhoi

real(c_double), private module_mp_morr_two_moment::rhoi
private
188  REAL(C_DOUBLE), PRIVATE :: RHOI ! BULK DENSITY OF CLOUD ICE

Referenced by make_sources(), and morr_two_moment_init().

◆ rhosn

real(c_double), private module_mp_morr_two_moment::rhosn
private
189  REAL(C_DOUBLE), PRIVATE :: RHOSN ! BULK DENSITY OF SNOW

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ rhosu

real(c_double), private module_mp_morr_two_moment::rhosu
private
186  REAL(C_DOUBLE), PRIVATE :: RHOSU ! STANDARD AIR DENSITY AT 850 MB

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ rhow

real(c_double), private module_mp_morr_two_moment::rhow
private
187  REAL(C_DOUBLE), PRIVATE :: RHOW ! DENSITY OF LIQUID WATER

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ rin

real(c_double), private module_mp_morr_two_moment::rin
private
206  REAL(C_DOUBLE), PRIVATE :: RIN ! RADIUS OF CONTACT NUCLEI (M)

Referenced by morr_two_moment_init(), and morr_two_moment_micro().

◆ rm1

real(c_double), private module_mp_morr_two_moment::rm1
private
226  REAL(C_DOUBLE), PRIVATE :: RM1 ! GEOMETRIC MEAN RADIUS, MODE 1 (M)

Referenced by morr_two_moment_init().

◆ rm2

real(c_double), private module_mp_morr_two_moment::rm2
private
227  REAL(C_DOUBLE), PRIVATE :: RM2 ! GEOMETRIC MEAN RADIUS, MODE 2 (M)

Referenced by morr_two_moment_init().

◆ rr

real(c_double), private module_mp_morr_two_moment::rr
private
224  REAL(C_DOUBLE), PRIVATE :: RR ! UNIVERSAL GAS CONSTANT

Referenced by init_zlevels(), morr_two_moment_init(), ERF::refinement_criteria_setup(), ERF::Write3DPlotFile(), and ERF::WriteMultiLevelPlotfileWithTerrain().

◆ sig1

real(c_double), private module_mp_morr_two_moment::sig1
private
230  REAL(C_DOUBLE), PRIVATE :: SIG1 ! STANDARD DEVIATION OF AEROSOL S.D., MODE 1

Referenced by morr_two_moment_init().

◆ sig2

real(c_double), private module_mp_morr_two_moment::sig2
private
231  REAL(C_DOUBLE), PRIVATE :: SIG2 ! STANDARD DEVIATION OF AEROSOL S.D., MODE 2

Referenced by morr_two_moment_init().

◆ vi

real(c_double), private module_mp_morr_two_moment::vi
private
219  REAL(C_DOUBLE), PRIVATE :: VI ! NUMBER OF ION DISSOCIATED IN SOLUTION

Referenced by polygon_::define(), and morr_two_moment_init().

◆ xxx

real(c_double), parameter, private module_mp_morr_two_moment::xxx = 0.9189385332046727417803297D0
private
101  REAL(C_DOUBLE), PARAMETER :: xxx = 0.9189385332046727417803297d0

Referenced by gamma().