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

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

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

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::average_down_mass_flux(), 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 if_average_down(), IBFaceSet::map_buildings_to_level0(), morr_two_moment_init(), morr_two_moment_micro(), erf_sbm::ProlongCarrierRelativeSpectrum(), redistribute_term(), ERF::turbPert_update_automatic(), and ERF::WarnSteepTerrainAtOpenFaces().

◆ 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 locate_latlon_on_grid(), ERF::make_subdomains(), 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(), ShocPDF::diagnose_pdf(), interp_trilinear(), morr_two_moment_init(), plotfile2d::sampled_target_is_bracketed(), TerrainScalarDiffusionFluxAtKFace(), and WallGradientCoeffs().

◆ cg

real(c_double), private module_mp_morr_two_moment::cg
private

Referenced by make_terrain_fitted_coords(), and 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 erf_sbm::SBMTransport::advance_stage(), ERF::cf_flux_history_targets(), ERF::cf_zface_layout(), ibseb::RayMap::column_top(), erf_sbm::flatten_constraint_groups(), ibseb::RayMap::locate(), morr_two_moment_init(), morr_two_moment_micro(), ibseb::ray_blocked(), ibseb::ray_hit(), and SurfaceLayer::set_pblh().

◆ 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 ComputeVerticalShear2_EB(), fill_fire_wind_from_interpolation(), ERF::ibseb_check_refined_levels(), ERFFillPatcher::InterpFace(), invert_latlon_offset(), locate_latlon_on_grid(), morr_two_moment_init(), morr_two_moment_micro(), ScalarDiffusionFaceCoefficient(), and erf_topography::shadowed().

◆ ds

real(c_double), private module_mp_morr_two_moment::ds
private

Referenced by morr_two_moment_init(), morr_two_moment_micro(), and erf_topography::shadowed().

◆ 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(), ERF::turbPert_update_automatic(), 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 ERF::advance_dycore(), interp_trilinear(), make_mom_sources(), morr_two_moment_init(), TerrainIF::operator()(), TerminalVelocity< RT >::RogersYau(), ERF::sample_stations(), 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().

◆ 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 erf_auxiliary::BuildMappedDryAirCarrierFluxRate(), SurfaceModel::is_field_mapped(), morr_two_moment_init(), SurfaceModel::register_radiation_input(), and SurfaceModel::register_radiation_output().

◆ 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 SLM::AdvanceSLM(), ERF::cf_flux_history_targets(), SLM::fluxes_canopy(), SLM::Init(), SLM::init_from_inputs(), SLM::initialize_surface_outputs(), morr_two_moment_init(), SLM::ReadCheckpoint(), SLM::WriteCheckpoint(), and SLM::writeSLM_Data().

◆ 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 IBFaceSet::add_heat_flux_to_source(), WSM6::AdvanceColumns(), IBFaceSet::compute_longwave(), IBFaceSet::compute_longwave_two_stream(), compute_max_reflectivity_dbz(), IBFaceSet::compute_shortwave_two_stream(), erf_cloud_chamber::deterministic_perturbation(), ERF::erf_enforce_hse(), ERF::fill_from_bndryregs(), gamma(), get_inverse_exner(), ERF::initHSE(), erf_cloud_chamber::initial_u_velocity_perturbation(), erf_cloud_chamber::initial_v_velocity_perturbation(), morr_two_moment_init(), morr_two_moment_micro(), IBFaceSet::scatter_diagnostics(), IBFaceSet::scatter_field(), erf_cloud_chamber::staggered_velocity_v_scale(), 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 SLM::AdvanceSLM(), SLM::fluxes_canopy(), morr_two_moment_init(), morr_two_moment_micro(), and SLM::solve_ground_skin_temperature().

◆ 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 ERF::average_down_mass_flux(), ERF::average_down_precip_accumulators(), ERF::cf_flux_history(), ERF::cf_flux_history_targets(), ERF::cf_zface_coarse_data(), ERF::cf_zface_layout(), ERF::fill_from_bndryregs(), if_cf_surface_mismatch(), init_zlevels(), SurfaceLayer::inject_pblh_from_coarser(), morr_two_moment_init(), ERF::poisson_wall_dist(), ERF::read_box_for_refinement(), 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().