Line data Source code
1 : module nucleate_ice
2 :
3 : !-------------------------------------------------------------------------------
4 : ! Purpose:
5 : ! A parameterization of ice nucleation.
6 : !
7 : ! *** This module is intended to be a "portable" code layer. Ideally it should
8 : ! *** not contain any use association of modules that belong to the model framework.
9 : !
10 : !
11 : ! Method:
12 : ! The current method is based on Liu & Penner (2005) & Liu et al. (2007)
13 : ! It related the ice nucleation with the aerosol number, temperature and the
14 : ! updraft velocity. It includes homogeneous freezing of sulfate & immersion
15 : ! freezing on mineral dust (soot disabled) in cirrus clouds, and
16 : ! Meyers et al. (1992) deposition nucleation in mixed-phase clouds
17 : !
18 : ! The effect of preexisting ice crystals on ice nucleation in cirrus clouds is included,
19 : ! and also consider the sub-grid variability of temperature in cirrus clouds,
20 : ! following X. Shi et al. ACP (2014).
21 : !
22 : ! Ice nucleation in mixed-phase clouds now uses classical nucleation theory (CNT),
23 : ! follows Y. Wang et al. ACP (2014), Hoose et al. (2010).
24 : !
25 : ! Authors:
26 : ! Xiaohong Liu, 01/2005, modifications by A. Gettelman 2009-2010
27 : ! Xiangjun Shi & Xiaohong Liu, 01/2014.
28 : !
29 : ! With help from C. C. Chen and B. Eaton (2014)
30 : !-------------------------------------------------------------------------------
31 :
32 : use wv_saturation, only: svp_water, svp_ice
33 :
34 : implicit none
35 : private
36 : save
37 :
38 : integer, parameter :: r8 = selected_real_kind(12)
39 :
40 : public :: nucleati_init, nucleati
41 :
42 : logical :: use_preexisting_ice
43 : logical :: use_hetfrz_classnuc
44 : logical :: use_incloud_nuc
45 : integer :: iulog
46 : real(r8) :: pi
47 : real(r8) :: mincld
48 :
49 : real(r8), parameter :: Shet = 1.3_r8 ! het freezing threshold
50 : real(r8), parameter :: rhoice = 0.5e3_r8 ! kg/m3, Wpice is not sensitive to rhoice
51 : real(r8), parameter :: minweff= 0.001_r8 ! m/s
52 : real(r8), parameter :: gamma1=1.0_r8
53 : real(r8), parameter :: gamma2=1.0_r8
54 : real(r8), parameter :: gamma3=2.0_r8
55 : real(r8), parameter :: gamma4=6.0_r8
56 :
57 : real(r8) :: ci
58 :
59 : !===============================================================================
60 : contains
61 : !===============================================================================
62 :
63 0 : subroutine nucleati_init( &
64 : use_preexisting_ice_in, use_hetfrz_classnuc_in, use_incloud_nuc_in, iulog_in, pi_in, &
65 : mincld_in)
66 :
67 : logical, intent(in) :: use_preexisting_ice_in
68 : logical, intent(in) :: use_hetfrz_classnuc_in
69 : logical, intent(in) :: use_incloud_nuc_in
70 : integer, intent(in) :: iulog_in
71 : real(r8), intent(in) :: pi_in
72 : real(r8), intent(in) :: mincld_in
73 :
74 0 : use_preexisting_ice = use_preexisting_ice_in
75 0 : use_hetfrz_classnuc = use_hetfrz_classnuc_in
76 0 : use_incloud_nuc = use_incloud_nuc_in
77 0 : iulog = iulog_in
78 0 : pi = pi_in
79 0 : mincld = mincld_in
80 :
81 0 : ci = rhoice*pi/6._r8
82 :
83 0 : end subroutine nucleati_init
84 :
85 : !===============================================================================
86 :
87 0 : subroutine nucleati( &
88 : wbar, tair, pmid, relhum, cldn, &
89 : qc, qi, ni_in, rhoair, &
90 : so4_num, dst_num, soot_num, subgrid, &
91 : nuci, onihf, oniimm, onidep, onimey, &
92 : wpice, weff, fhom, regm, &
93 : oso4_num, odst_num, osoot_num, &
94 : call_frm_zm_in, add_preexisting_ice_in)
95 :
96 : ! Input Arguments
97 : real(r8), intent(in) :: wbar ! grid cell mean vertical velocity (m/s)
98 : real(r8), intent(in) :: tair ! temperature (K)
99 : real(r8), intent(in) :: pmid ! pressure at layer midpoints (pa)
100 : real(r8), intent(in) :: relhum ! relative humidity with respective to liquid
101 : real(r8), intent(in) :: cldn ! new value of cloud fraction (fraction)
102 : real(r8), intent(in) :: qc ! liquid water mixing ratio (kg/kg)
103 : real(r8), intent(in) :: qi ! grid-mean preexisting cloud ice mass mixing ratio (kg/kg)
104 : real(r8), intent(in) :: ni_in ! grid-mean preexisting cloud ice number conc (#/kg)
105 : real(r8), intent(in) :: rhoair ! air density (kg/m3)
106 : real(r8), intent(in) :: so4_num ! so4 aerosol number (#/cm^3)
107 : real(r8), intent(in) :: dst_num ! total dust aerosol number (#/cm^3)
108 : real(r8), intent(in) :: soot_num ! soot (hydrophilic) aerosol number (#/cm^3)
109 : real(r8), intent(in) :: subgrid ! subgrid saturation scaling factor
110 :
111 : ! Output Arguments
112 : real(r8), intent(out) :: nuci ! ice number nucleated (#/kg)
113 : real(r8), intent(out) :: onihf ! nucleated number from homogeneous freezing of so4
114 : real(r8), intent(out) :: oniimm ! nucleated number from immersion freezing
115 : real(r8), intent(out) :: onidep ! nucleated number from deposition nucleation
116 : real(r8), intent(out) :: onimey ! nucleated number from deposition nucleation (meyers: mixed phase)
117 : real(r8), intent(out) :: wpice ! diagnosed Vertical velocity Reduction caused by preexisting ice (m/s), at Shom
118 : real(r8), intent(out) :: weff ! effective Vertical velocity for ice nucleation (m/s); weff=wbar-wpice
119 : real(r8), intent(out) :: fhom ! how much fraction of cloud can reach Shom
120 : real(r8), intent(out) :: regm ! nucleation regime indiator
121 : real(r8), intent(out) :: oso4_num ! so4 aerosol number (#/cm^3)
122 : real(r8), intent(out) :: odst_num ! total dust aerosol number (#/cm^3)
123 : real(r8), intent(out) :: osoot_num ! soot (hydrophilic) aerosol number (#/cm^3)
124 :
125 : ! Optional Arguments
126 : logical, intent(in), optional :: call_frm_zm_in ! true if called from ZM convection scheme
127 : logical, intent(in), optional :: add_preexisting_ice_in ! only false if called with pumas_v1.21+
128 :
129 : ! Local workspace
130 : real(r8) :: nihf ! nucleated number from homogeneous freezing of so4
131 : real(r8) :: niimm ! nucleated number from immersion freezing
132 : real(r8) :: nidep ! nucleated number from deposition nucleation
133 : real(r8) :: nimey ! nucleated number from deposition nucleation (meyers)
134 : real(r8) :: n1, ni ! nucleated number
135 : real(r8) :: tc, A, B ! work variable
136 : real(r8) :: esl, esi, deles ! work variable
137 : real(r8) :: wbar1, wbar2
138 :
139 : ! used in SUBROUTINE Vpreice
140 : real(r8) :: Ni_preice ! cloud ice number conc (1/m3)
141 : real(r8) :: lami,Ri_preice ! mean cloud ice radius (m)
142 : real(r8) :: Shom ! initial ice saturation ratio; if <1, use hom threshold Si
143 : real(r8) :: detaT,RHimean ! temperature standard deviation, mean cloudy RHi
144 : real(r8) :: wpicehet ! diagnosed Vertical velocity Reduction caused by preexisting ice (m/s), at shet
145 :
146 : real(r8) :: weffhet ! effective Vertical velocity for ice nucleation (m/s) weff=wbar-wpicehet
147 :
148 : logical :: call_frm_zm, add_preexisting_ice
149 : !-------------------------------------------------------------------------------
150 :
151 0 : nuci = 0._r8
152 :
153 0 : RHimean = relhum*svp_water(tair)/svp_ice(tair)*subgrid
154 :
155 : ! temp variables that depend on use_preexisting_ice
156 0 : wbar1 = wbar
157 0 : wbar2 = wbar
158 :
159 : ! If not using prexisting ice, the homogeneous freezing happens in the
160 : ! entire gridbox.
161 0 : fhom = 1._r8
162 :
163 0 : if (present(call_frm_zm_in)) then
164 0 : call_frm_zm = call_frm_zm_in
165 : else
166 : call_frm_zm = .false.
167 : end if
168 :
169 0 : if (present(add_preexisting_ice_in)) then
170 0 : add_preexisting_ice = add_preexisting_ice_in
171 : else
172 : add_preexisting_ice = .true.
173 : end if
174 :
175 0 : if (use_preexisting_ice .and. (.not. call_frm_zm)) then
176 :
177 0 : Ni_preice = ni_in*rhoair ! (convert from #/kg -> #/m3)
178 0 : Ni_preice = Ni_preice / max(mincld,cldn) ! in-cloud ice number density
179 :
180 0 : if (Ni_preice > 10.0_r8 .and. qi > 1.e-10_r8) then ! > 0.01/L = 10/m3
181 0 : Shom = -1.5_r8 ! if Shom<1 , Shom will be recalculated in SUBROUTINE Vpreice, according to Ren & McKenzie, 2005
182 0 : lami = (gamma4*ci*ni_in/qi)**(1._r8/3._r8)
183 0 : Ri_preice = 0.5_r8/lami ! radius
184 0 : Ri_preice = max(Ri_preice, 1e-8_r8) ! >0.01micron
185 0 : call Vpreice(pmid, tair, Ri_preice, Ni_preice, Shom, wpice)
186 0 : call Vpreice(pmid, tair, Ri_preice, Ni_preice, Shet, wpicehet)
187 : else
188 0 : wpice = 0.0_r8
189 0 : wpicehet = 0.0_r8
190 : endif
191 :
192 0 : weff = max(wbar-wpice, minweff)
193 0 : wpice = min(wpice, wbar)
194 0 : weffhet = max(wbar-wpicehet,minweff)
195 0 : wpicehet = min(wpicehet, wbar)
196 :
197 0 : wbar1 = weff
198 0 : wbar2 = weffhet
199 :
200 0 : detaT = wbar/0.23_r8
201 0 : if (use_incloud_nuc) then
202 0 : call frachom(tair, 1._r8, detaT, fhom)
203 : else
204 0 : call frachom(tair, RHimean, detaT, fhom)
205 : end if
206 : end if
207 :
208 0 : ni = 0._r8
209 0 : tc = tair - 273.15_r8
210 :
211 : ! initialize
212 0 : niimm = 0._r8
213 0 : nidep = 0._r8
214 0 : nihf = 0._r8
215 0 : deles = 0._r8
216 0 : esi = 0._r8
217 0 : regm = 0._r8
218 :
219 0 : oso4_num = 0._r8
220 0 : odst_num = 0._r8
221 0 : osoot_num = 0._r8
222 :
223 0 : if ((so4_num >= 1.0e-10_r8 .or. (soot_num+dst_num) >= 1.0e-10_r8) .and. cldn > 0._r8) then
224 :
225 0 : if (RHimean.ge.1.2_r8) then
226 :
227 0 : if ( ((tc.le.0.0_r8).and.(tc.ge.-37.0_r8).and.(qc.lt.1.e-12_r8)).or.(tc.le.-37.0_r8)) then
228 :
229 0 : A = -1.4938_r8 * log(soot_num+dst_num) + 12.884_r8
230 0 : B = -10.41_r8 * log(soot_num+dst_num) - 67.69_r8
231 0 : regm = A * log(wbar1) + B
232 :
233 : ! heterogeneous nucleation only
234 0 : if (tc .gt. regm .or. so4_num < 1.0e-10_r8) then
235 :
236 0 : if(tc.lt.-40._r8 .and. wbar1.gt.1._r8 .and. so4_num >= 1.0e-10_r8) then ! exclude T<-40 & W>1m/s from hetero. nucleation
237 :
238 0 : call hf(tc,wbar1,relhum*subgrid,so4_num,nihf)
239 0 : niimm=0._r8
240 0 : nidep=0._r8
241 :
242 : ! If some homogeneous nucleation happened, assume all of the that heterogeneous
243 : ! and coarse mode sulfate particles nucleated.
244 0 : if (nihf.gt.1e-3_r8) then ! hom occur, add preexisting ice
245 0 : niimm = dst_num + soot_num ! assuming dst_num freeze firstly
246 0 : odst_num = dst_num
247 0 : osoot_num = soot_num
248 :
249 0 : oso4_num = nihf
250 : endif
251 :
252 0 : nihf = nihf * fhom
253 0 : oso4_num = oso4_num * fhom
254 :
255 0 : n1 = nihf + niimm
256 : else
257 :
258 0 : call hetero(tc,wbar2,soot_num+dst_num,niimm,nidep)
259 :
260 0 : nihf = 0._r8
261 0 : n1 = niimm + nidep
262 :
263 0 : osoot_num = soot_num * (niimm + nidep) / (soot_num + dst_num)
264 0 : odst_num = dst_num * (niimm + nidep) / (soot_num + dst_num)
265 : endif
266 :
267 : ! homogeneous nucleation only
268 0 : else if (tc.lt.regm-5._r8 .or. (soot_num+dst_num) < 1.0e-10_r8) then
269 :
270 0 : call hf(tc,wbar1,relhum*subgrid,so4_num,nihf)
271 0 : niimm=0._r8
272 0 : nidep=0._r8
273 :
274 : ! If some homogeneous nucleation happened, assume all of the that
275 : ! heterogeneous and coarse mode sulfate particles nucleated.
276 0 : if (nihf.gt.1e-3_r8) then ! hom occur, add preexisting ice
277 0 : niimm = dst_num + soot_num ! assuming dst_num freeze firstly
278 0 : odst_num = dst_num
279 0 : osoot_num = soot_num
280 :
281 0 : oso4_num = nihf
282 : endif
283 :
284 0 : nihf = nihf * fhom
285 0 : oso4_num = oso4_num * fhom
286 :
287 0 : n1 = nihf + niimm
288 :
289 : ! transition between homogeneous and heterogeneous: interpolate in-between
290 : else
291 :
292 0 : if (tc.lt.-40._r8 .and. wbar1.gt.1._r8) then ! exclude T<-40 & W>1m/s from hetero. nucleation
293 :
294 0 : call hf(tc, wbar1, relhum*subgrid, so4_num, nihf)
295 0 : niimm = 0._r8
296 0 : nidep = 0._r8
297 :
298 : ! If some homogeneous nucleation happened, assume all of the
299 : ! that heterogeneous and coarse mode sulfate particles nucleated.
300 0 : if (nihf.gt.1e-3_r8) then ! hom occur, add preexisting ice
301 0 : niimm = dst_num + soot_num ! assuming dst_num freeze firstly
302 0 : odst_num = dst_num
303 0 : osoot_num = soot_num
304 :
305 0 : oso4_num = nihf
306 : endif
307 :
308 0 : nihf = nihf * fhom
309 0 : oso4_num = oso4_num * fhom
310 :
311 0 : n1 = nihf + niimm
312 :
313 : else
314 :
315 0 : call hf(regm-5._r8,wbar1,relhum*subgrid,so4_num,nihf)
316 0 : call hetero(regm,wbar2,soot_num+dst_num,niimm,nidep)
317 :
318 : ! If some homogeneous nucleation happened, assume all of the
319 : ! heterogeneous particles nucleated and add in a fraction of
320 : ! the homogeneous freezing.
321 0 : if (nihf.gt.1e-3_r8) then ! hom occur, add preexisting ice
322 0 : oso4_num = nihf
323 : endif
324 :
325 0 : osoot_num = soot_num * (niimm + nidep) / (soot_num + dst_num)
326 0 : odst_num = dst_num * (niimm + nidep) / (soot_num + dst_num)
327 :
328 0 : nihf = nihf * fhom * ((regm - tc) / 5._r8)**2
329 0 : oso4_num = oso4_num * fhom * ((regm - tc) / 5._r8)**2
330 :
331 0 : n1 = niimm + nidep + nihf
332 :
333 : end if
334 : end if
335 :
336 : ! Scale the rates for in-cloud number, since this is what
337 : ! MG is expecting to find.
338 0 : ni = n1
339 :
340 : ! If using prexsiting ice, and allowed to add, then add it to the total.
341 0 : if (use_preexisting_ice) then
342 0 : if (add_preexisting_ice .and. (.not. call_frm_zm)) then
343 0 : ni = ni + Ni_preice * 1e-6_r8
344 : end if
345 : end if
346 : end if
347 : end if
348 : end if
349 :
350 : ! deposition/condensation nucleation in mixed clouds (-37<T<0C) (Meyers, 1992)
351 0 : if(tc.lt.0._r8 .and. tc.gt.-37._r8 .and. qc.gt.1.e-12_r8) then
352 0 : esl = svp_water(tair) ! over water in mixed clouds
353 0 : esi = svp_ice(tair) ! over ice
354 0 : deles = (esl - esi)
355 0 : nimey=1.e-3_r8*exp(12.96_r8*deles/esi - 0.639_r8)
356 : else
357 : nimey=0._r8
358 : endif
359 :
360 0 : if (use_hetfrz_classnuc) nimey = 0._r8
361 :
362 0 : nuci=ni + nimey
363 :
364 0 : if(nuci.gt.9999._r8.or.nuci.lt.0._r8) then
365 0 : write(iulog, *) 'Warning: incorrect ice nucleation number (nuci reset =0)'
366 0 : write(iulog, *) ni, tair, relhum, wbar, nihf, niimm, nidep,deles,esi,dst_num,so4_num
367 0 : nuci=0._r8
368 : endif
369 :
370 0 : nuci = nuci*1.e+6_r8/rhoair ! change unit from #/cm3 to #/kg
371 0 : onimey = nimey*1.e+6_r8/rhoair
372 0 : onidep = nidep*1.e+6_r8/rhoair
373 0 : oniimm = niimm*1.e+6_r8/rhoair
374 0 : onihf = nihf*1.e+6_r8/rhoair
375 0 : end subroutine nucleati
376 :
377 : !===============================================================================
378 :
379 0 : subroutine hetero(T,ww,Ns,Nis,Nid)
380 :
381 : real(r8), intent(in) :: T, ww, Ns
382 : real(r8), intent(out) :: Nis, Nid
383 :
384 : real(r8) A11,A12,A21,A22,B11,B12,B21,B22
385 : real(r8) B,C
386 :
387 : !---------------------------------------------------------------------
388 : ! parameters
389 :
390 0 : A11 = 0.0263_r8
391 0 : A12 = -0.0185_r8
392 0 : A21 = 2.758_r8
393 0 : A22 = 1.3221_r8
394 0 : B11 = -0.008_r8
395 0 : B12 = -0.0468_r8
396 0 : B21 = -0.2667_r8
397 0 : B22 = -1.4588_r8
398 :
399 : ! ice from immersion nucleation (cm^-3)
400 :
401 0 : B = (A11+B11*log(Ns)) * log(ww) + (A12+B12*log(Ns))
402 0 : C = A21+B21*log(Ns)
403 :
404 0 : Nis = exp(A22) * Ns**B22 * exp(B*T) * ww**C
405 0 : Nis = min(Nis,Ns)
406 :
407 0 : Nid = 0.0_r8 ! don't include deposition nucleation for cirrus clouds when T<-37C
408 :
409 0 : end subroutine hetero
410 :
411 : !===============================================================================
412 :
413 0 : subroutine hf(T,ww,RH,Na,Ni)
414 :
415 : real(r8), intent(in) :: T, ww, RH, Na
416 : real(r8), intent(out) :: Ni
417 :
418 : real(r8) A1_fast,A21_fast,A22_fast,B1_fast,B21_fast,B22_fast
419 : real(r8) A2_fast,B2_fast
420 : real(r8) C1_fast,C2_fast,k1_fast,k2_fast
421 : real(r8) A1_slow,A2_slow,B1_slow,B2_slow,B3_slow
422 : real(r8) C1_slow,C2_slow,k1_slow,k2_slow
423 : real(r8) regm
424 : real(r8) A,B,C
425 : real(r8) RHw
426 :
427 : !---------------------------------------------------------------------
428 : ! parameters
429 :
430 0 : A1_fast =0.0231_r8
431 0 : A21_fast =-1.6387_r8 !(T>-64 deg)
432 0 : A22_fast =-6.045_r8 !(T<=-64 deg)
433 0 : B1_fast =-0.008_r8
434 0 : B21_fast =-0.042_r8 !(T>-64 deg)
435 0 : B22_fast =-0.112_r8 !(T<=-64 deg)
436 0 : C1_fast =0.0739_r8
437 0 : C2_fast =1.2372_r8
438 :
439 0 : A1_slow =-0.3949_r8
440 0 : A2_slow =1.282_r8
441 0 : B1_slow =-0.0156_r8
442 0 : B2_slow =0.0111_r8
443 0 : B3_slow =0.0217_r8
444 0 : C1_slow =0.120_r8
445 0 : C2_slow =2.312_r8
446 :
447 0 : Ni = 0.0_r8
448 :
449 : !----------------------------
450 : !RHw parameters
451 0 : A = 6.0e-4_r8*log(ww)+6.6e-3_r8
452 0 : B = 6.0e-2_r8*log(ww)+1.052_r8
453 0 : C = 1.68_r8 *log(ww)+129.35_r8
454 0 : RHw=(A*T*T+B*T+C)*0.01_r8
455 :
456 0 : if((T.le.-37.0_r8) .and. ((RH).ge.RHw)) then
457 :
458 0 : regm = 6.07_r8*log(ww)-55.0_r8
459 :
460 0 : if(T.ge.regm) then ! fast-growth regime
461 :
462 0 : if(T.gt.-64.0_r8) then
463 : A2_fast=A21_fast
464 : B2_fast=B21_fast
465 : else
466 0 : A2_fast=A22_fast
467 0 : B2_fast=B22_fast
468 : endif
469 :
470 0 : k1_fast = exp(A2_fast + B2_fast*T + C2_fast*log(ww))
471 0 : k2_fast = A1_fast+B1_fast*T+C1_fast*log(ww)
472 :
473 0 : Ni = k1_fast*Na**(k2_fast)
474 0 : Ni = min(Ni,Na)
475 :
476 : else ! slow-growth regime
477 :
478 0 : k1_slow = exp(A2_slow + (B2_slow+B3_slow*log(ww))*T + C2_slow*log(ww))
479 0 : k2_slow = A1_slow+B1_slow*T+C1_slow*log(ww)
480 :
481 0 : Ni = k1_slow*Na**(k2_slow)
482 0 : Ni = min(Ni,Na)
483 :
484 : endif
485 :
486 : end if
487 :
488 0 : end subroutine hf
489 :
490 : !===============================================================================
491 :
492 0 : SUBROUTINE Vpreice(P_in, T_in, R_in, C_in, S_in, V_out)
493 :
494 : ! based on Karcher et al. (2006)
495 : ! VERTICAL VELOCITY CALCULATED FROM DEPOSITIONAL LOSS TERM
496 :
497 : ! SUBROUTINE arguments
498 : REAL(r8), INTENT(in) :: P_in ! [Pa],INITIAL AIR pressure
499 : REAL(r8), INTENT(in) :: T_in ! [K] ,INITIAL AIR temperature
500 : REAL(r8), INTENT(in) :: R_in ! [m],INITIAL MEAN ICE CRYSTAL NUMBER RADIUS
501 : REAL(r8), INTENT(in) :: C_in ! [m-3],INITIAL TOTAL ICE CRYSTAL NUMBER DENSITY, [1/cm3]
502 : REAL(r8), INTENT(in) :: S_in ! [-],INITIAL ICE SATURATION RATIO;; if <1, use hom threshold Si
503 : REAL(r8), INTENT(out) :: V_out ! [m/s], VERTICAL VELOCITY REDUCTION (caused by preexisting ice)
504 :
505 : ! SUBROUTINE parameters
506 : REAL(r8), PARAMETER :: ALPHAc = 0.5_r8 ! density of ice (g/cm3), !!!V is not related to ALPHAc
507 : REAL(r8), PARAMETER :: FA1c = 0.601272523_r8
508 : REAL(r8), PARAMETER :: FA2c = 0.000342181855_r8
509 : REAL(r8), PARAMETER :: FA3c = 1.49236645E-12_r8
510 : REAL(r8), PARAMETER :: WVP1c = 3.6E+10_r8
511 : REAL(r8), PARAMETER :: WVP2c = 6145.0_r8
512 : REAL(r8), PARAMETER :: FVTHc = 11713803.0_r8
513 : REAL(r8), PARAMETER :: THOUBKc = 7.24637701E+18_r8
514 : REAL(r8), PARAMETER :: SVOLc = 3.23E-23_r8 ! SVOL=XMW/RHOICE
515 : REAL(r8), PARAMETER :: FDc = 249.239822_r8
516 : REAL(r8), PARAMETER :: FPIVOLc = 3.89051704E+23_r8
517 : REAL(r8) :: T,P,S,R,C
518 : REAL(r8) :: A1,A2,A3,B1,B2
519 : REAL(r8) :: T_1,PICE,FLUX,ALP4,CISAT,DLOSS,VICE
520 :
521 0 : T = T_in ! K , K
522 0 : P = P_in*1e-2_r8 ! Pa , hpa
523 :
524 0 : IF (S_in.LT.1.0_r8) THEN
525 0 : S = 2.349_r8 - (T/259.0_r8) ! homogeneous freezing threshold, according to Ren & McKenzie, 2005
526 : ELSE
527 : S = S_in ! INPUT ICE SATURATION RATIO, -, >1
528 : ENDIF
529 :
530 0 : R = R_in*1e2_r8 ! m => cm
531 0 : C = C_in*1e-6_r8 ! m-3 => cm-3
532 0 : T_1 = 1.0_r8/ T
533 0 : PICE = WVP1c * EXP(-(WVP2c*T_1))
534 0 : ALP4 = 0.25_r8 * ALPHAc
535 0 : FLUX = ALP4 * SQRT(FVTHc*T)
536 0 : CISAT = THOUBKc * PICE * T_1
537 0 : A1 = ( FA1c * T_1 - FA2c ) * T_1
538 0 : A2 = 1.0_r8/ CISAT
539 0 : A3 = FA3c * T_1 / P
540 0 : B1 = FLUX * SVOLc * CISAT * ( S-1.0_r8 )
541 0 : B2 = FLUX * FDc * P * T_1**1.94_r8
542 0 : DLOSS = FPIVOLc * C * B1 * R**2 / ( 1.0_r8+ B2 * R )
543 0 : VICE = ( A2 + A3 * S ) * DLOSS / ( A1 * S ) ! 2006,(19)
544 0 : V_out = VICE*1e-2_r8 ! cm/s => m/s
545 :
546 0 : END SUBROUTINE Vpreice
547 :
548 0 : subroutine frachom(Tmean,RHimean,detaT,fhom)
549 : ! How much fraction of cirrus might reach Shom
550 : ! base on "A cirrus cloud scheme for general circulation models",
551 : ! B. Karcher and U. Burkhardt 2008
552 :
553 : real(r8), intent(in) :: Tmean, RHimean, detaT
554 : real(r8), intent(out) :: fhom
555 :
556 : real(r8), parameter :: seta = 6132.9_r8 ! K
557 : integer, parameter :: Nbin=200 ! (Tmean - 3*detaT, Tmean + 3*detaT)
558 :
559 : real(r8) :: PDF_T(Nbin) ! temperature PDF; ! PDF_T=0 outside (Tmean-3*detaT, Tmean+3*detaT)
560 : real(r8) :: Sbin(Nbin) ! the fluctuations of Si that are driven by the T variations
561 : real(r8) :: Sihom, deta
562 : integer :: i
563 :
564 0 : Sihom = 2.349_r8-Tmean/259.0_r8 ! homogeneous freezing threshold, according to Ren & McKenzie, 2005
565 0 : fhom = 0.0_r8
566 :
567 0 : do i = Nbin, 1, -1
568 :
569 0 : deta = (i - 0.5_r8 - Nbin/2)*6.0_r8/Nbin ! PDF_T=0 outside (Tmean-3*detaT, Tmean+3*detaT)
570 0 : Sbin(i) = RHimean*exp(deta*detaT*seta/Tmean**2.0_r8)
571 0 : PDF_T(i) = exp(-deta**2.0_r8/2.0_r8)*6.0_r8/(sqrt(2.0_r8*Pi)*Nbin)
572 :
573 :
574 0 : if (Sbin(i).ge.Sihom) then
575 0 : fhom = fhom + PDF_T(i)
576 : else
577 : exit
578 : end if
579 : end do
580 :
581 0 : fhom = min(1.0_r8, fhom/0.997_r8) ! accounting for the finite limits (-3 , 3)
582 0 : end subroutine frachom
583 :
584 : end module nucleate_ice
585 :
|
