Line data Source code
1 : !---------------------------------------------------------------------------
2 : ! $Id$
3 : !===============================================================================
4 : module Nc_Ncn_eqns
5 :
6 : ! Description:
7 : ! Equations are provided to perform calculations back-and-forth between Nc and
8 : ! Ncn, where Nc is cloud droplet concentration and Ncn is simplified cloud
9 : ! nuclei concentration. The equation that relates the two is:
10 : !
11 : ! Nc = Ncn * H(chi);
12 : !
13 : ! where chi is extended liquid water mixing ratio, which is equal to cloud
14 : ! water mixing ratio, rc, when both are positive. However, chi is negative in
15 : ! subsaturated air.
16 : !
17 : ! Equation are provided relating mean cloud droplet concentration (overall),
18 : ! Ncm, and/or mean cloud droplet concentration (in-cloud), Nc_in_cloud, to
19 : ! mean simplified cloud nuclei concentration, Ncnm.
20 :
21 : ! Notes:
22 : !
23 : ! Meaning of Nc flag combinations:
24 : !
25 : ! l_const_Nc_in_cloud:
26 : ! When this flag is enabled, cloud droplet concentration (in-cloud) is
27 : ! constant (spatially) at a grid level (it is constant over the subgrid
28 : ! domain, but could vary over time depending on the value of l_predict_Nc).
29 : ! The value of in-cloud Nc does not vary at a grid level. This also means
30 : ! that Ncn is constant across the entire grid level. When this flag is turned
31 : ! off, both in-cloud Nc and Ncn vary at a grid level.
32 : !
33 : ! l_predict_Nc:
34 : ! When this flag is enabled, Nc_in_cloud (or alternatively Ncm) is predicted.
35 : ! It is advanced every time step by a predictive equation, and can change
36 : ! at every time step at a grid level. When this flag is turned off,
37 : ! Nc_in_cloud does not change at a grid level over the course of a model run.
38 : !
39 : ! 1) l_predict_Nc turned on and l_const_Nc_in_cloud turned on:
40 : ! The value of Nc_in_cloud (mean in-cloud Nc) is predicted and can change
41 : ! at every timestep at a grid level. However, the value of in-cloud Nc is
42 : ! constant (spatially) at a grid level (no subgrid variability).
43 : !
44 : ! 2) l_predict_Nc turned on and l_const_Nc_in_cloud turned off:
45 : ! The value of Nc_in_cloud (mean in-cloud Nc) is predicted and can change
46 : ! at every timestep at a grid level. The value of in-cloud Nc also varies
47 : ! (spatially) at a grid level (subgrid variability around mean
48 : ! in-cloud Nc).
49 : !
50 : ! 3) l_predict_Nc turned off and l_const_Nc_in_cloud turned on:
51 : ! The value of Nc_in_cloud (mean in-cloud Nc) is constant over time at a
52 : ! grid level. It retains its initial value. Additionally, the value of
53 : ! in-cloud Nc is constant (spatially) at a grid level (no subgrid
54 : ! variability). This configuration is used most often in idealized cases.
55 : !
56 : ! 4) l_predict_Nc turned off and l_const_Nc_in_cloud turned off:
57 : ! The value of Nc_in_cloud (mean in-cloud Nc) is constant over time at a
58 : ! grid level. It retains its initial value. However, the value of
59 : ! in-cloud Nc varies (spatially) at a grid level (subgrid variability
60 : ! around mean in-cloud Nc).
61 : !
62 : !
63 : !
64 : ! Nc_in_cloud/Nc - Ncn flow chart of CLUBB code:
65 : !
66 : ! (Please update when warranted).
67 : !
68 : !
69 : ! Ncm/Nc-in-cloud Ncnm/Ncn PDF params.
70 : ! --->
71 : ! | | Start of CLUBB main time step loop
72 : ! | |
73 : ! | | advance_clubb_core
74 : ! | |
75 : ! | |
76 : ! | |\
77 : ! | | \
78 : ! | | (intent in)-------setup_pdf_parameters-------->calc. Ncnm (local)
79 : ! | | |
80 : ! | | \ /
81 : ! | | mu_Ncn_i, sigma_Ncn_i,
82 : ! | | corr_xNcn_i
83 : ! | | |
84 : ! | | \ /
85 : ! | | PDF param. arrays:
86 : ! | | mu_x_i_n, sigma_x_i_n,
87 : ! | | corr_array_i_n
88 : ! | | (intent out)
89 : ! | | |
90 : ! | | |
91 : ! | | |
92 : ! | | |
93 : ! | | |
94 : ! | | |
95 : ! | |--(intent in)---calc_microphys_scheme_tendcies----(intent in)
96 : ! | | |
97 : ! | | |
98 : ! | | call a microphysics scheme
99 : ! | | |
100 : ! | | Local micro. scheme-----------Latin Hypercube-----------Upscaled KK
101 : ! | | | | |
102 : ! | | Ncm/Nc-in-cloud: Populate sample points Use PDF params.
103 : ! | | used to find micro. using PDF params (Ncn). of Ncn
104 : ! | | tendencies. At every sample point: (mu_Ncn_i, etc.)
105 : ! | | | Nc = Ncn * H(chi). to find micro.
106 : ! | | | Use sample-point Nc to tendencies.
107 : ! | | | find micro. tendencies |
108 : ! | | | when calling micro. scheme. |
109 : ! | | | | |
110 : ! | | hydromet_mc/-----------------hydromet_mc/-------------hydromet_mc
111 : ! | | Ncm_mc (intent out) | Ncm_mc (intent out) (intent out)
112 : ! | | |
113 : ! | | |
114 : ! | | |
115 : ! | | |
116 : ! | | |
117 : ! | | |
118 : ! | | |
119 : ! | | (intent in)
120 : ! | | |
121 : ! | |--(intent inout)----advance_microphys
122 : ! | |
123 : ! | |
124 : ! | | advance microphysics variables (hydromet, Nc_in_cloud/Ncm) one timestep
125 : ! | |
126 : ! | | l_predict_Nc = true:
127 : ! | | Nc_in_cloud/Ncm necessary for starting
128 : ! | | value of Nc_in_cloud/Ncm when advancing
129 : ! | | one timestep using predictive equation.
130 : ! | |
131 : ! | |
132 : ! | | End of CLUBB main time step loop
133 : ! <---
134 :
135 : ! References:
136 : !-------------------------------------------------------------------------
137 :
138 : implicit none
139 :
140 : private ! default scope
141 :
142 : public :: Ncnm_to_Nc_in_cloud, &
143 : Nc_in_cloud_to_Ncnm, &
144 : Ncnm_to_Ncm, &
145 : Ncm_to_Ncnm
146 :
147 : private :: bivar_NL_chi_Ncn_mean, &
148 : bivar_Ncnm_eqn_comp
149 :
150 : contains
151 :
152 : !=============================================================================
153 0 : elemental function Ncnm_to_Nc_in_cloud( mu_chi_1, mu_chi_2, mu_Ncn_1, &
154 : mu_Ncn_2, sigma_chi_1, sigma_chi_2, &
155 : sigma_Ncn_1, sigma_Ncn_2, &
156 : sigma_Ncn_1_n, sigma_Ncn_2_n, &
157 : corr_chi_Ncn_1_n, corr_chi_Ncn_2_n, &
158 : mixt_frac, cloud_frac_1, &
159 : cloud_frac_2 ) &
160 : result( Nc_in_cloud )
161 :
162 : ! Description:
163 : ! The in-cloud mean of cloud droplet concentration is calculated from the
164 : ! PDF parameters involving simplified cloud nuclei concentration, Ncn, and
165 : ! cloud fraction. At any point, cloud droplet concentration, Nc, is given
166 : ! by:
167 : !
168 : ! Nc = Ncn * H(chi);
169 : !
170 : ! where extended liquid water mixing ratio, chi, is equal to cloud water
171 : ! ratio, rc, when positive. When the atmosphere is saturated at this point,
172 : ! cloud water is found, and Nc = Ncn. Otherwise, only clear air is found,
173 : ! and Nc = 0.
174 : !
175 : ! The overall mean of cloud droplet concentration, <Nc>, is calculated from
176 : ! the PDF parameters involving Ncn. The in-cloud mean of cloud droplet
177 : ! concentration is calculated from <Nc> and cloud fraction.
178 :
179 : ! References:
180 : !-----------------------------------------------------------------------
181 :
182 : use constants_clubb, only: &
183 : one, & ! Constant(s)
184 : cloud_frac_min
185 :
186 : use clubb_precision, only: &
187 : core_rknd ! Variable(s)
188 :
189 : implicit none
190 :
191 : ! Input Variables
192 : real( kind = core_rknd ), intent(in) :: &
193 : mu_chi_1, & ! Mean of chi (old s) (1st PDF component) [kg/kg]
194 : mu_chi_2, & ! Mean of chi (old s) (2nd PDF component) [kg/kg]
195 : mu_Ncn_1, & ! Mean of Ncn (1st PDF component) [num/kg]
196 : mu_Ncn_2, & ! Mean of Ncn (2nd PDF component) [num/kg]
197 : sigma_chi_1, & ! Standard deviation of chi (1st PDF comp.) [kg/kg]
198 : sigma_chi_2, & ! Standard deviation of chi (2nd PDF comp.) [kg/kg]
199 : sigma_Ncn_1, & ! Standard deviation of Ncn (1st PDF comp.) [num/kg]
200 : sigma_Ncn_2, & ! Standard deviation of Ncn (2nd PDF comp.) [num/kg]
201 : sigma_Ncn_1_n, & ! Standard deviation of ln Ncn (1st PDF component) [-]
202 : sigma_Ncn_2_n, & ! Standard deviation of ln Ncn (2nd PDF component) [-]
203 : corr_chi_Ncn_1_n, & ! Correlation of chi and ln Ncn (1st PDF comp.) [-]
204 : corr_chi_Ncn_2_n, & ! Correlation of chi and ln Ncn (2nd PDF comp.) [-]
205 : mixt_frac, & ! Mixture fraction [-]
206 : cloud_frac_1, & ! Cloud fraction (1st PDF component) [-]
207 : cloud_frac_2 ! Cloud fraction (2nd PDF component) [-]
208 :
209 : ! Return Variable
210 : real( kind = core_rknd ) :: &
211 : Nc_in_cloud ! Mean cloud droplet concentration (in-cloud) [num/kg]
212 :
213 : ! Local Variable
214 : real( kind = core_rknd ) :: &
215 : Ncm, & ! Mean cloud droplet concentration (overall) [num/kg]
216 : cloud_frac ! Cloud fraction [-]
217 :
218 :
219 : ! Calculate overall cloud fraction as calculated by the PDF.
220 : ! The variable cloud_frac is not used here because it is altered by factors
221 : ! such as the trapezoidal rule calculation.
222 : ! Cloud fraction can be recalculated here from cloud_frac_1 and cloud_frac_2
223 : ! as long neither of these variables are altered by any factor. They can
224 : ! only be calculated from PDF.
225 0 : cloud_frac = mixt_frac * cloud_frac_1 + ( one - mixt_frac ) * cloud_frac_2
226 :
227 0 : if ( cloud_frac > cloud_frac_min ) then
228 :
229 : ! There is cloud found at this grid level. Calculate Nc_in_cloud.
230 : Ncm = Ncnm_to_Ncm( mu_chi_1, mu_chi_2, mu_Ncn_1, mu_Ncn_2, &
231 : sigma_chi_1, sigma_chi_2, sigma_Ncn_1, &
232 : sigma_Ncn_2, sigma_Ncn_1_n, sigma_Ncn_2_n, &
233 0 : corr_chi_Ncn_1_n, corr_chi_Ncn_2_n, mixt_frac )
234 :
235 0 : Nc_in_cloud = Ncm / cloud_frac
236 :
237 : else ! cloud_frac <= cloud_frac_min
238 :
239 : ! This level is entirely clear. Set Nc_in_cloud to <Ncn>.
240 : ! Since <Ncn> = mu_Ncn_1 = mu_Ncn_2, use mu_Ncn_1 here.
241 0 : Nc_in_cloud = mu_Ncn_1
242 :
243 : endif
244 :
245 :
246 : return
247 :
248 : end function Ncnm_to_Nc_in_cloud
249 :
250 : !=============================================================================
251 0 : elemental function Nc_in_cloud_to_Ncnm( mu_chi_1, mu_chi_2, sigma_chi_1, &
252 : sigma_chi_2, mixt_frac, Nc_in_cloud, &
253 : cloud_frac_1, cloud_frac_2, &
254 : const_Ncnp2_on_Ncnm2, &
255 : const_corr_chi_Ncn ) &
256 : result( Ncnm )
257 :
258 : ! Description:
259 : ! The overall mean of simplified cloud nuclei concentration, <Ncn>, is
260 : ! calculated from the in-cloud mean of cloud droplet concentration, <Nc>,
261 : ! cloud fraction, and some of the PDF parameters.
262 : !
263 : ! At any point, cloud droplet concentration, Nc, is given by:
264 : !
265 : ! Nc = Ncn * H(chi);
266 : !
267 : ! where extended liquid water mixing ratio, chi, is equal to cloud water
268 : ! ratio, rc, when positive. When the atmosphere is saturated at this point,
269 : ! cloud water is found, and Nc = Ncn. Otherwise, only clear air is found,
270 : ! and Nc = 0.
271 : !
272 : ! The overall mean of cloud droplet concentration, <Nc>, is calculated from
273 : ! Nc_in_cloud and cloud fraction. The value of <Ncn> is calculated from
274 : ! <Nc> and PDF parameters.
275 :
276 : ! References:
277 : !-----------------------------------------------------------------------
278 :
279 : use constants_clubb, only: &
280 : one, & ! Constant(s)
281 : cloud_frac_min, &
282 : eps
283 :
284 : use clubb_precision, only: &
285 : core_rknd ! Variable(s)
286 :
287 : implicit none
288 :
289 : ! Input Variables
290 : real( kind = core_rknd ), intent(in) :: &
291 : mu_chi_1, & ! Mean of chi (old s) (1st PDF component) [kg/kg]
292 : mu_chi_2, & ! Mean of chi (old s) (2nd PDF component) [kg/kg]
293 : sigma_chi_1, & ! Standard deviation of chi (1st PDF component) [kg/kg]
294 : sigma_chi_2, & ! Standard deviation of chi (2nd PDF component) [kg/kg]
295 : mixt_frac ! Mixture fraction [-]
296 :
297 : real( kind = core_rknd ), intent(in) :: &
298 : Nc_in_cloud, & ! Mean cloud droplet conc. (in-cloud) [num/kg]
299 : cloud_frac_1, & ! Cloud fraction (1st PDF component) [-]
300 : cloud_frac_2, & ! Cloud fraction (2nd PDF component) [-]
301 : const_Ncnp2_on_Ncnm2, & ! Prescribed ratio of <Ncn'^2> to <Ncn>^2 [-]
302 : const_corr_chi_Ncn ! Prescribed correlation of chi and Ncn [-]
303 :
304 : ! Return Variable
305 : real( kind = core_rknd ) :: &
306 : Ncnm ! Mean simplified cloud nuclei concentration (overall) [num/kg]
307 :
308 : ! Local Variable
309 : real( kind = core_rknd ) :: &
310 : Ncm, & ! Mean cloud droplet concentration (overall) [num/kg]
311 : cloud_frac ! Cloud fraction [-]
312 :
313 :
314 : ! Calculate overall cloud fraction as calculated by the PDF.
315 : ! The variable cloud_frac is not used here because it is altered by factors
316 : ! such as the trapezoidal rule calculation.
317 : ! Cloud fraction can be recalculated here from cloud_frac_1 and cloud_frac_2
318 : ! as long neither of these variables are altered by any factor. They can
319 : ! only be calculated from the PDF.
320 0 : cloud_frac = mixt_frac * cloud_frac_1 + ( one - mixt_frac ) * cloud_frac_2
321 :
322 : if ( cloud_frac > cloud_frac_min &
323 0 : .and. abs(const_corr_chi_Ncn * const_Ncnp2_on_Ncnm2) > eps ) then
324 :
325 : ! There is cloud found at this grid level. Additionally, Ncn varies.
326 : ! Calculate Nc_in_cloud.
327 0 : Ncm = Nc_in_cloud * cloud_frac
328 :
329 : Ncnm = Ncm_to_Ncnm( mu_chi_1, mu_chi_2, sigma_chi_1, sigma_chi_2, &
330 : mixt_frac, Ncm, const_Ncnp2_on_Ncnm2, &
331 0 : const_corr_chi_Ncn, Nc_in_cloud )
332 :
333 : else ! cloud_frac <= cloud_frac_min .or. const_Ncnp2_on_Ncnm2 = 0
334 :
335 : ! When Ncn is constant a a grid level, it is equal to Nc_in_cloud.
336 : ! Additionally, when a level is entirely clear, <Ncn>, which is based on
337 : ! Nc_in_cloud, here, must be set to something. Set <Ncn> to Nc_in_cloud.
338 0 : Ncnm = Nc_in_cloud
339 :
340 : endif
341 :
342 :
343 : return
344 :
345 : end function Nc_in_cloud_to_Ncnm
346 :
347 : !=============================================================================
348 0 : elemental function Ncnm_to_Ncm( mu_chi_1, mu_chi_2, mu_Ncn_1, mu_Ncn_2, &
349 : sigma_chi_1, sigma_chi_2, sigma_Ncn_1, &
350 : sigma_Ncn_2, sigma_Ncn_1_n, sigma_Ncn_2_n, &
351 : corr_chi_Ncn_1_n, corr_chi_Ncn_2_n, &
352 : mixt_frac ) &
353 : result( Ncm )
354 :
355 : ! Description:
356 : ! The overall mean of cloud droplet concentration, <Nc>, is calculated from
357 : ! the PDF parameters involving the simplified cloud nuclei concentration,
358 : ! Ncn. At any point, cloud droplet concentration, Nc, is given by:
359 : !
360 : ! Nc = Ncn * H(chi);
361 : !
362 : ! where extended liquid water mixing ratio, chi, is equal to cloud water
363 : ! ratio, rc, when positive. When the atmosphere is saturated at this point,
364 : ! cloud water is found, and Nc = Ncn. Otherwise, only clear air is found,
365 : ! and Nc = 0.
366 : !
367 : ! The overall mean of cloud droplet concentration, <Nc>, is found by
368 : ! integrating over the PDF of chi and Ncn, such that:
369 : !
370 : ! <Nc> = INT(-inf:inf) INT(0:inf) Ncn * H(chi) * P(chi,Ncn) dNcn dchi;
371 : !
372 : ! which can also be written as:
373 : !
374 : ! <Nc> = SUM(i=1,n) mixt_frac_i
375 : ! * INT(-inf:inf) INT(0:inf) Ncn * H(chi) * P_i(chi,Ncn) dNcn dchi;
376 : !
377 : ! where n is the number of multivariate joint PDF components, mixt_frac_i is
378 : ! the weight of the ith PDF component, and P_i is the functional form of the
379 : ! multivariate joint PDF in the ith PDF component.
380 : !
381 : ! This equation is rewritten as:
382 : !
383 : ! <Nc> = SUM(i=1,n) mixt_frac_i
384 : ! * INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi.
385 : !
386 : ! When both chi and Ncn vary in the ith PDF component, the integral is
387 : ! evaluated and the result is:
388 : !
389 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi
390 : ! = (1/2) * exp{ mu_Ncn_i_n + (1/2) * sigma_Ncn_i_n^2 }
391 : ! * erfc( - ( 1 / sqrt(2) ) * ( ( mu_chi_i / sigma_chi_i )
392 : ! + rho_chi_Ncn_i_n * sigma_Ncn_i_n ) );
393 : !
394 : ! which can be reduced to:
395 : !
396 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi
397 : ! = (1/2) * mu_Ncn_i
398 : ! * erfc( - ( 1 / sqrt(2) ) * ( ( mu_chi_i / sigma_chi_i )
399 : ! + rho_chi_Ncn_i_n * sigma_Ncn_i_n ) ).
400 : !
401 : ! When chi is constant, but Ncn varies, in the ith PDF component, the
402 : ! integral is evaluated and results in:
403 : !
404 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi = mu_Ncn_i;
405 : !
406 : ! when mu_chi_i > 0; and
407 : !
408 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi = 0;
409 : !
410 : ! when mu_chi_i <= 0.
411 : !
412 : ! When chi varies, but Ncn is constant, in the ith PDF component, the
413 : ! integral is evaluated and results in:
414 : !
415 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi
416 : ! = mu_Ncn_i * (1/2) * erfc( - ( mu_chi_i / ( sqrt(2) * sigma_chi_i ) ) ).
417 : !
418 : ! When both chi and Ncn are constant in the ith PDF component, the integral
419 : ! is evaluated and results in:
420 : !
421 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi = mu_Ncn_i;
422 : !
423 : ! when mu_chi_i > 0; and
424 : !
425 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi = 0;
426 : !
427 : ! when mu_chi_i <= 0.
428 :
429 : ! References:
430 : !-----------------------------------------------------------------------
431 :
432 : use constants_clubb, only: &
433 : one ! Constant(s)
434 :
435 : use clubb_precision, only: &
436 : core_rknd ! Variable(s)
437 :
438 : implicit none
439 :
440 : ! Input Variables
441 : real( kind = core_rknd ), intent(in) :: &
442 : mu_chi_1, & ! Mean of chi (old s) (1st PDF component) [kg/kg]
443 : mu_chi_2, & ! Mean of chi (old s) (2nd PDF component) [kg/kg]
444 : mu_Ncn_1, & ! Mean of Ncn (1st PDF component) [num/kg]
445 : mu_Ncn_2, & ! Mean of Ncn (2nd PDF component) [num/kg]
446 : sigma_chi_1, & ! Standard deviation of chi (1st PDF comp.) [kg/kg]
447 : sigma_chi_2, & ! Standard deviation of chi (2nd PDF comp.) [kg/kg]
448 : sigma_Ncn_1, & ! Standard deviation of Ncn (1st PDF comp.) [num/kg]
449 : sigma_Ncn_2, & ! Standard deviation of Ncn (2nd PDF comp.) [num/kg]
450 : sigma_Ncn_1_n, & ! Standard deviation of ln Ncn (1st PDF component) [-]
451 : sigma_Ncn_2_n, & ! Standard deviation of ln Ncn (2nd PDF component) [-]
452 : corr_chi_Ncn_1_n, & ! Correlation of chi and ln Ncn (1st PDF comp.) [-]
453 : corr_chi_Ncn_2_n, & ! Correlation of chi and ln Ncn (2nd PDF comp.) [-]
454 : mixt_frac ! Mixture fraction [-]
455 :
456 : ! Return Variable
457 : real( kind = core_rknd ) :: &
458 : Ncm ! Mean cloud droplet concentration (overall) [num/kg]
459 :
460 :
461 : ! Calculate mean cloud droplet concentration (overall), <Nc>.
462 : Ncm &
463 : = mixt_frac &
464 : * bivar_NL_chi_Ncn_mean( mu_chi_1, mu_Ncn_1, sigma_chi_1, &
465 : sigma_Ncn_1, sigma_Ncn_1_n, corr_chi_Ncn_1_n ) &
466 : + ( one - mixt_frac ) &
467 : * bivar_NL_chi_Ncn_mean( mu_chi_2, mu_Ncn_2, sigma_chi_2, &
468 0 : sigma_Ncn_2, sigma_Ncn_2_n, corr_chi_Ncn_2_n )
469 :
470 :
471 : return
472 :
473 : end function Ncnm_to_Ncm
474 :
475 : !=============================================================================
476 0 : elemental function Ncm_to_Ncnm( mu_chi_1, mu_chi_2, sigma_chi_1, &
477 : sigma_chi_2, mixt_frac, Ncm, &
478 : const_Ncnp2_on_Ncnm2, const_corr_chi_Ncn, &
479 : Ncnm_val_denom_0 ) &
480 : result( Ncnm )
481 :
482 : ! Description:
483 : ! The overall mean of simplified cloud nuclei concentration, <Ncn>, is
484 : ! calculated from the overall mean of cloud droplet concentration, <Nc>, and
485 : ! some of the PDF parameters.
486 : !
487 : ! At any point, cloud droplet concentration, Nc, is given by:
488 : !
489 : ! Nc = Ncn * H(chi);
490 : !
491 : ! where extended liquid water mixing ratio, chi, is equal to cloud water
492 : ! ratio, rc, when positive. When the atmosphere is saturated at this point,
493 : ! cloud water is found, and Nc = Ncn. Otherwise, only clear air is found,
494 : ! and Nc = 0.
495 : !
496 : ! The overall mean of cloud droplet concentration, <Nc>, is found by
497 : ! integrating over the PDF of chi and Ncn, such that:
498 : !
499 : ! <Nc> = INT(-inf:inf) INT(0:inf) Ncn * H(chi) * P(chi,Ncn) dNcn dchi;
500 : !
501 : ! which can also be written as:
502 : !
503 : ! <Nc> = SUM(i=1,n) mixt_frac_i
504 : ! * INT(-inf:inf) INT(0:inf) Ncn * H(chi) * P_i(chi,Ncn) dNcn dchi;
505 : !
506 : ! where n is the number of multivariate joint PDF components, mixt_frac_i is
507 : ! the weight of the ith PDF component, and P_i is the functional form of the
508 : ! multivariate joint PDF in the ith PDF component.
509 : !
510 : ! This equation is rewritten as:
511 : !
512 : ! <Nc> = SUM(i=1,n) mixt_frac_i
513 : ! * INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi.
514 : !
515 : ! When both chi and Ncn vary in the ith PDF component, the integral is
516 : ! evaluated and the result is:
517 : !
518 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi
519 : ! = (1/2) * exp{ mu_Ncn_i_n + (1/2) * sigma_Ncn_i_n^2 }
520 : ! * erfc( - ( 1 / sqrt(2) ) * ( ( mu_chi_i / sigma_chi_i )
521 : ! + rho_chi_Ncn_i_n * sigma_Ncn_i_n ) );
522 : !
523 : ! which can be reduced to:
524 : !
525 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi
526 : ! = (1/2) * mu_Ncn_i
527 : ! * erfc( - ( 1 / sqrt(2) ) * ( ( mu_chi_i / sigma_chi_i )
528 : ! + rho_chi_Ncn_i_n * sigma_Ncn_i_n ) ).
529 : !
530 : ! When chi is constant, but Ncn varies, in the ith PDF component, the
531 : ! integral is evaluated and results in:
532 : !
533 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi = mu_Ncn_i;
534 : !
535 : ! when mu_chi_i > 0; and
536 : !
537 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi = 0;
538 : !
539 : ! when mu_chi_i <= 0.
540 : !
541 : ! When chi varies, but Ncn is constant, in the ith PDF component, the
542 : ! integral is evaluated and results in:
543 : !
544 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi
545 : ! = mu_Ncn_i * (1/2) * erfc( - ( mu_chi_i / ( sqrt(2) * sigma_chi_i ) ) ).
546 : !
547 : ! When both chi and Ncn are constant in the ith PDF component, the integral
548 : ! is evaluated and results in:
549 : !
550 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi = mu_Ncn_i;
551 : !
552 : ! when mu_chi_i > 0; and
553 : !
554 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi = 0;
555 : !
556 : ! when mu_chi_i <= 0.
557 : !
558 : !
559 : ! Solving for <Ncn>
560 : ! =================
561 : !
562 : ! The individual marginal for simplified cloud nuclei concentration, Ncn, is
563 : ! a single lognormal distribution over the entire horizontal domain. In
564 : ! order to accomplish this in a two-component PDF structure, the PDF
565 : ! parameters involving Ncn are set equal between the two components. This
566 : ! results in:
567 : !
568 : ! mu_Ncn_1 = mu_Ncn_2 = mu_Ncn_i = <Ncn>;
569 : ! mu_Ncn_1_n = mu_Ncn_2_n = mu_Ncn_i_n;
570 : ! sigma_Ncn_1 = sigma_Ncn_2 = sigma_Ncn_i = sqrt( <Ncn'^2> );
571 : ! sigma_Ncn_1_n = sigma_Ncn_2_n = sigma_Ncn_i_n;
572 : ! rho_chi_Ncn_1 = rho_chi_Ncn_2 = rho_chi_Ncn_i = rho_chi_Ncn; and
573 : ! rho_chi_Ncn_1_n = rho_chi_Ncn_2_n = rho_chi_Ncn_i_n.
574 : !
575 : ! Additionally, the equation for sigma_Ncn_i_n is:
576 : !
577 : ! sigma_Ncn_i_n = sqrt( ln( 1 + ( sigma_Ncn_i^2 / mu_Ncn_i^2 ) ) );
578 : !
579 : ! and the equation for rho_chi_Ncn_i_n is:
580 : !
581 : ! rho_chi_Ncn_i_n
582 : ! = rho_chi_Ncn_i * sqrt( exp{ sigma_Ncn_i_n^2 } - 1 ) / sigma_Ncn_i_n.
583 : !
584 : ! The product of rho_chi_Ncn_i_n and sigma_Ncn_i_n is:
585 : !
586 : ! rho_chi_Ncn_i_n * sigma_Ncn_i_n
587 : ! = rho_chi_Ncn_i * sqrt( exp{ sigma_Ncn_i_n^2 } - 1 ).
588 : !
589 : ! After substituting for sigma_Ncn_i_n^2, the equation for the product of
590 : ! rho_chi_Ncn_i_n and sigma_Ncn_i_n is:
591 : !
592 : ! rho_chi_Ncn_i_n * sigma_Ncn_i_n
593 : ! = rho_chi_Ncn_i * sqrt( sigma_Ncn_i^2 / mu_Ncn_i^2 );
594 : !
595 : ! which can be rewritten as:
596 : !
597 : ! rho_chi_Ncn_i_n * sigma_Ncn_i_n
598 : ! = rho_chi_Ncn * sqrt( <Ncn'^2> / <Ncn>^2 ).
599 : !
600 : ! Substituting all of this into the equation for <Nc>, the equation for <Nc>
601 : ! becomes:
602 : !
603 : ! <Nc> = <Ncn>
604 : ! * SUM(i=1,n) mixt_frac_i
605 : ! ---
606 : ! | (1/2) * erfc( - ( 1 / sqrt(2) )
607 : ! | * ( ( mu_chi_i / sigma_chi_i )
608 : ! | + rho_chi_Ncn * sqrt(<Ncn'^2>/<Ncn>^2) ) );
609 : ! | where sigma_chi_i > 0 and <Ncn'^2> > 0;
610 : ! |
611 : ! * | (1/2) * erfc( - ( mu_chi_i / ( sqrt(2) * sigma_chi_i ) ) );
612 : ! | where sigma_chi_i > 0 and <Ncn'^2> = 0;
613 : ! |
614 : ! | 1; where sigma_chi_i = 0 and mu_chi_i > 0;
615 : ! |
616 : ! | 0; where sigma_chi_i = 0 and mu_chi_i <= 0.
617 : ! ---
618 : !
619 : ! In order to isolate <Ncn>, the value of <Ncn'^2>/<Ncn>^2 is set to a
620 : ! constant value, const_Ncn. The value of this constant does not depend on
621 : ! <Ncn>. Likewise, the value of rho_chi_Ncn does not depend on <Ncn>.
622 : ! Solving for <Ncn>, the equation becomes:
623 : !
624 : ! <Ncn>
625 : ! = <Nc> / ( SUM(i=1,n) mixt_frac_i
626 : ! ---
627 : ! | (1/2) * erfc( - ( 1 / sqrt(2) )
628 : ! | * ( ( mu_chi_i / sigma_chi_i )
629 : ! | + rho_chi_Ncn * sqrt( const_Ncn ) ) );
630 : ! | where sigma_chi_i > 0 and const_Ncn > 0;
631 : ! |
632 : ! * | (1/2) * erfc( - ( mu_chi_i / ( sqrt(2) * sigma_chi_i ) ) );
633 : ! | where sigma_chi_i > 0 and const_Ncn = 0;
634 : ! |
635 : ! | 1; where sigma_chi_i = 0 and mu_chi_i > 0;
636 : ! |
637 : ! | 0; where sigma_chi_i = 0 and mu_chi_i <= 0 ).
638 : ! ---
639 : !
640 : ! When the denominator term is 0, there is only clear air. Both the
641 : ! numerator (<Nc>) and the denominator have a value of 0, and <Ncn> is set
642 : ! to an appropriate value.
643 :
644 : ! References:
645 : !-----------------------------------------------------------------------
646 :
647 : use constants_clubb, only: &
648 : one, & ! Constant(s)
649 : zero
650 :
651 : use clubb_precision, only: &
652 : core_rknd ! Variable(s)
653 :
654 : implicit none
655 :
656 : ! Input Variables
657 : real( kind = core_rknd ), intent(in) :: &
658 : mu_chi_1, & ! Mean of chi (old s) (1st PDF component) [kg/kg]
659 : mu_chi_2, & ! Mean of chi (old s) (2nd PDF component) [kg/kg]
660 : sigma_chi_1, & ! Standard deviation of chi (1st PDF component) [kg/kg]
661 : sigma_chi_2, & ! Standard deviation of chi (2nd PDF component) [kg/kg]
662 : mixt_frac ! Mixture fraction [-]
663 :
664 : real( kind = core_rknd ), intent(in) :: &
665 : Ncm, & ! Mean cloud droplet conc. (overall) [num/kg]
666 : const_Ncnp2_on_Ncnm2, & ! Prescribed ratio of <Ncn'^2> to <Ncn>^2 [-]
667 : const_corr_chi_Ncn, & ! Prescribed correlation of chi and Ncn [-]
668 : Ncnm_val_denom_0 ! Ncnm value -- denominator in eqn. is 0 [num/kg]
669 :
670 : ! Return Variable
671 : real( kind = core_rknd ) :: &
672 : Ncnm ! Mean simplified cloud nuclei concentration (overall) [num/kg]
673 :
674 : ! Local Variable
675 : real( kind = core_rknd ) :: &
676 : denominator_term ! Denominator in the equation for <Ncn> [-]
677 :
678 :
679 : denominator_term &
680 : = mixt_frac &
681 : * bivar_Ncnm_eqn_comp( mu_chi_1, sigma_chi_1, &
682 : const_Ncnp2_on_Ncnm2, const_corr_chi_Ncn ) &
683 : + ( one - mixt_frac ) &
684 : * bivar_Ncnm_eqn_comp( mu_chi_2, sigma_chi_2, &
685 0 : const_Ncnp2_on_Ncnm2, const_corr_chi_Ncn )
686 :
687 :
688 0 : if ( denominator_term > zero ) then
689 :
690 0 : Ncnm = Ncm / denominator_term
691 :
692 : else ! denominator_term = 0
693 :
694 : ! When the denominator is 0, it is usually because there is only clear
695 : ! air. In that scenario, Ncm should also be 0. Set Ncnm to a value that
696 : ! is usual or typical
697 0 : Ncnm = Ncnm_val_denom_0
698 :
699 : endif ! denominator_term > 0
700 :
701 :
702 : return
703 :
704 : end function Ncm_to_Ncnm
705 :
706 : !=============================================================================
707 0 : elemental function bivar_NL_chi_Ncn_mean( mu_chi_i, mu_Ncn_i, sigma_chi_i, &
708 : sigma_Ncn_i, sigma_Ncn_i_n, &
709 : corr_chi_Ncn_i_n )
710 :
711 : ! Description:
712 : ! The double integral over Ncn * H(chi) multiplied by the
713 : ! bivariate normal-lognormal joint PDF of chi and Ncn is evaluated. The
714 : ! integral is given by:
715 : !
716 : ! INT(-inf:inf) INT(0:inf) Ncn * H(chi) * P_i(chi,Ncn) dNcn dchi;
717 : !
718 : ! which reduces to:
719 : !
720 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi;
721 : !
722 : ! where the individual marginal distribution of chi is normal in the ith PDF
723 : ! component and the individual marginal distribution of Ncn is lognormal in
724 : ! the ith PDF component.
725 : !
726 : ! When both chi and Ncn vary in the ith PDF component, the integral is
727 : ! evaluated and the result is:
728 : !
729 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi
730 : ! = (1/2) * exp{ mu_Ncn_i_n + (1/2) * sigma_Ncn_i_n^2 }
731 : ! * erfc( - ( 1 / sqrt(2) ) * ( ( mu_chi_i / sigma_chi_i )
732 : ! + rho_chi_Ncn_i_n * sigma_Ncn_i_n ) );
733 : !
734 : ! which can be reduced to:
735 : !
736 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi
737 : ! = (1/2) * mu_Ncn_i
738 : ! * erfc( - ( 1 / sqrt(2) ) * ( ( mu_chi_i / sigma_chi_i )
739 : ! + rho_chi_Ncn_i_n * sigma_Ncn_i_n ) ).
740 : !
741 : ! When chi is constant, but Ncn varies, in the ith PDF component, the
742 : ! integral is evaluated and results in:
743 : !
744 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi = mu_Ncn_i;
745 : !
746 : ! when mu_chi_i > 0; and
747 : !
748 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi = 0;
749 : !
750 : ! when mu_chi_i <= 0.
751 : !
752 : ! When chi varies, but Ncn is constant, in the ith PDF component, the
753 : ! integral is evaluated and results in:
754 : !
755 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi
756 : ! = mu_Ncn_i * (1/2) * erfc( - ( mu_chi_i / ( sqrt(2) * sigma_chi_i ) ) ).
757 : !
758 : ! When both chi and Ncn are constant in the ith PDF component, the integral
759 : ! is evaluated and results in:
760 : !
761 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi = mu_Ncn_i;
762 : !
763 : ! when mu_chi_i > 0; and
764 : !
765 : ! INT(0:inf) INT(0:inf) Ncn * P_i(chi,Ncn) dNcn dchi = 0;
766 : !
767 : ! when mu_chi_i <= 0.
768 :
769 : ! References:
770 : !-----------------------------------------------------------------------
771 :
772 : use constants_clubb, only: &
773 : sqrt_2, & ! Constant(s)
774 : one, &
775 : one_half, &
776 : zero, &
777 : chi_tol, &
778 : Ncn_tol
779 :
780 : use clubb_precision, only: &
781 : core_rknd ! Variable(s)
782 :
783 : implicit none
784 :
785 : ! Input Variables
786 : real( kind = core_rknd ), intent(in) :: &
787 : mu_chi_i, & ! Mean of chi (old s) (ith PDF component) [kg/kg]
788 : mu_Ncn_i, & ! Mean of Ncn (ith PDF component) [num/kg]
789 : sigma_chi_i, & ! Standard deviation of chi (ith PDF comp.) [kg/kg]
790 : sigma_Ncn_i, & ! Standard deviation of Ncn (ith PDF comp.) [num/kg]
791 : sigma_Ncn_i_n, & ! Standard deviation of ln Ncn (ith PDF component) [-]
792 : corr_chi_Ncn_i_n ! Correlation of chi and ln Ncn (ith PDF comp.) [-]
793 :
794 : ! Return Variable
795 : real( kind = core_rknd ) :: &
796 : bivar_NL_chi_Ncn_mean
797 :
798 :
799 0 : if ( sigma_chi_i <= chi_tol .and. sigma_Ncn_i <= Ncn_tol ) then
800 :
801 : ! The ith PDF component variances of both chi and Ncn are 0.
802 :
803 0 : if ( mu_chi_i > zero ) then
804 :
805 0 : bivar_NL_chi_Ncn_mean = mu_Ncn_i
806 :
807 : else ! mu_chi_i <= 0
808 :
809 : bivar_NL_chi_Ncn_mean = zero
810 :
811 : endif
812 :
813 :
814 0 : elseif ( sigma_chi_i <= chi_tol ) then
815 :
816 : ! The ith PDF component variance of chi is 0.
817 :
818 0 : if ( mu_chi_i > zero ) then
819 :
820 0 : bivar_NL_chi_Ncn_mean = mu_Ncn_i
821 :
822 : else ! mu_chi_i <= 0
823 :
824 : bivar_NL_chi_Ncn_mean = zero
825 :
826 : endif
827 :
828 :
829 0 : elseif ( sigma_Ncn_i <= Ncn_tol ) then
830 :
831 : ! The ith PDF component variance of Ncn is 0.
832 :
833 : bivar_NL_chi_Ncn_mean &
834 0 : = mu_Ncn_i * one_half * erfc( - ( mu_chi_i / ( sqrt_2 * sigma_chi_i ) ) )
835 :
836 :
837 : else
838 :
839 : ! Both chi and Ncn vary in the ith PDF component.
840 :
841 : bivar_NL_chi_Ncn_mean &
842 : = one_half * mu_Ncn_i &
843 : * erfc( - ( one / sqrt_2 ) &
844 : * ( ( mu_chi_i / sigma_chi_i ) &
845 0 : + corr_chi_Ncn_i_n * sigma_Ncn_i_n ) )
846 :
847 :
848 : endif
849 :
850 :
851 : return
852 :
853 : end function bivar_NL_chi_Ncn_mean
854 :
855 : !=============================================================================
856 0 : elemental function bivar_Ncnm_eqn_comp( mu_chi_i, sigma_chi_i, &
857 : const_Ncnp2_on_Ncnm2, &
858 : const_corr_chi_Ncn )
859 :
860 : ! Description:
861 : ! When <Ncn> is found based on the value of <Nc>, the following equation is
862 : ! used:
863 : !
864 : ! <Ncn>
865 : ! = <Nc> / ( SUM(i=1,n) mixt_frac_i
866 : ! ---
867 : ! | (1/2) * erfc( - ( 1 / sqrt(2) )
868 : ! | * ( ( mu_chi_i / sigma_chi_i )
869 : ! | + rho_chi_Ncn * sqrt( const_Ncn ) ) );
870 : ! | where sigma_chi_i > 0 and const_Ncn > 0;
871 : ! |
872 : ! * | (1/2) * erfc( - ( mu_chi_i / ( sqrt(2) * sigma_chi_i ) ) );
873 : ! | where sigma_chi_i > 0 and const_Ncn = 0;
874 : ! |
875 : ! | 1; where sigma_chi_i = 0 and mu_chi_i > 0;
876 : ! |
877 : ! | 0; where sigma_chi_i = 0 and mu_chi_i <= 0 ).
878 : ! ---
879 : !
880 : ! In the above equation, const_Ncn = <Ncn'^2> / <Ncn>^2. It is a constant,
881 : ! prescribed parameter. Likewise, rho_chi_Ncn is a parameter that is not
882 : ! based on the value of <Ncn>.
883 : !
884 : ! When the denominator term is 0, there is only clear air. Both the
885 : ! numerator (<Nc>) and the denominator have a value of 0, and <Ncn> is set
886 : ! to an appropriate value.
887 : !
888 : ! The contribution of the ith PDF component to the denominator term in the
889 : ! equation is calculated here.
890 :
891 : ! References:
892 : !-----------------------------------------------------------------------
893 :
894 : use constants_clubb, only: &
895 : sqrt_2, & ! Constant(s)
896 : one, &
897 : one_half, &
898 : zero, &
899 : chi_tol
900 :
901 : use clubb_precision, only: &
902 : core_rknd ! Variable(s)
903 :
904 : implicit none
905 :
906 : ! Input Variables
907 : real( kind = core_rknd ), intent(in) :: &
908 : mu_chi_i, & ! Mean of chi (old s) (ith PDF component) [kg/kg]
909 : sigma_chi_i ! Standard deviation of chi (ith PDF component) [kg/kg]
910 :
911 : real( kind = core_rknd ), intent(in) :: &
912 : const_Ncnp2_on_Ncnm2, & ! Prescribed ratio of <Ncn'^2> to <Ncn>^2 [-]
913 : const_corr_chi_Ncn ! Prescribed correlation of chi and Ncn [-]
914 :
915 : ! Return Variable
916 : real( kind = core_rknd ) :: &
917 : bivar_Ncnm_eqn_comp
918 :
919 :
920 0 : if ( sigma_chi_i <= chi_tol ) then
921 :
922 : ! The ith PDF component variances of chi is 0. The value of the ith PDF
923 : ! component variance of Ncn does not matter in this scenario.
924 :
925 0 : if ( mu_chi_i > zero ) then
926 :
927 : bivar_Ncnm_eqn_comp = one
928 :
929 : else ! mu_chi_i <= 0
930 :
931 0 : bivar_Ncnm_eqn_comp = zero
932 :
933 : endif
934 :
935 : else
936 :
937 : ! Both chi and Ncn vary in the ith PDF component.
938 :
939 : bivar_Ncnm_eqn_comp &
940 : = one_half * erfc( - ( one / sqrt_2 ) * ( ( mu_chi_i / sigma_chi_i ) &
941 0 : + const_corr_chi_Ncn * sqrt( const_Ncnp2_on_Ncnm2 ) ) )
942 :
943 :
944 : endif
945 :
946 :
947 : return
948 :
949 : end function bivar_Ncnm_eqn_comp
950 :
951 : !===============================================================================
952 :
953 : end module Nc_Ncn_eqns
|