source: palm/trunk/SOURCE/lpm_droplet_condensation.f90 @ 1952

Last change on this file since 1952 was 1891, checked in by hoffmann, 8 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 22.5 KB
Line 
1!> @file lpm_droplet_condensation.f90
2!--------------------------------------------------------------------------------!
3! This file is part of PALM.
4!
5! PALM is free software: you can redistribute it and/or modify it under the terms
6! of the GNU General Public License as published by the Free Software Foundation,
7! either version 3 of the License, or (at your option) any later version.
8!
9! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
10! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
11! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
12!
13! You should have received a copy of the GNU General Public License along with
14! PALM. If not, see <http://www.gnu.org/licenses/>.
15!
16! Copyright 1997-2016 Leibniz Universitaet Hannover
17!--------------------------------------------------------------------------------!
18!
19! Current revisions:
20! ------------------
21!
22!
23! Former revisions:
24! -----------------
25! $Id: lpm_droplet_condensation.f90 1891 2016-04-22 08:53:22Z suehring $
26!
27! 1890 2016-04-22 08:52:11Z hoffmann
28! Some improvements of the Rosenbrock method. If the Rosenbrock method needs more
29! than 40 iterations to find a sufficient time setp, the model is not aborted.
30! This might lead to small erros. But they can be assumend as negligible, since
31! the maximum timesetp of the Rosenbrock method after 40 iterations will be
32! smaller than 10^-16 s. 
33!
34! 1871 2016-04-15 11:46:09Z hoffmann
35! Initialization of aerosols added.
36!
37! 1849 2016-04-08 11:33:18Z hoffmann
38! Interpolation of supersaturation has been removed because it is not in
39! accordance with the release/depletion of latent heat/water vapor in
40! interaction_droplets_ptq.
41! Calculation of particle Reynolds number has been corrected.
42! eps_ros added from modules.
43!
44! 1831 2016-04-07 13:15:51Z hoffmann
45! curvature_solution_effects moved to particle_attributes
46!
47! 1822 2016-04-07 07:49:42Z hoffmann
48! Unused variables removed.
49!
50! 1682 2015-10-07 23:56:08Z knoop
51! Code annotations made doxygen readable
52!
53! 1359 2014-04-11 17:15:14Z hoffmann
54! New particle structure integrated.
55! Kind definition added to all floating point numbers.
56!
57! 1346 2014-03-27 13:18:20Z heinze
58! Bugfix: REAL constants provided with KIND-attribute especially in call of
59! intrinsic function like MAX, MIN, SIGN
60!
61! 1322 2014-03-20 16:38:49Z raasch
62! REAL constants defined as wp-kind
63!
64! 1320 2014-03-20 08:40:49Z raasch
65! ONLY-attribute added to USE-statements,
66! kind-parameters added to all INTEGER and REAL declaration statements,
67! kinds are defined in new module kinds,
68! comment fields (!:) to be used for variable explanations added to
69! all variable declaration statements
70!
71! 1318 2014-03-17 13:35:16Z raasch
72! module interfaces removed
73!
74! 1092 2013-02-02 11:24:22Z raasch
75! unused variables removed
76!
77! 1071 2012-11-29 16:54:55Z franke
78! Ventilation effect for evaporation of large droplets included
79! Check for unreasonable results included in calculation of Rosenbrock method
80! since physically unlikely results were observed and for the same
81! reason the first internal time step in Rosenbrock method should be < 1.0E02 in
82! case of evaporation
83! Unnecessary calculation of ql_int removed
84! Unnecessary calculations in Rosenbrock method (d2rdt2, drdt_m, dt_ros_last)
85! removed
86! Bugfix: factor in calculation of surface tension changed from 0.00155 to
87! 0.000155
88!
89! 1036 2012-10-22 13:43:42Z raasch
90! code put under GPL (PALM 3.9)
91!
92! 849 2012-03-15 10:35:09Z raasch
93! initial revision (former part of advec_particles)
94!
95!
96! Description:
97! ------------
98!> Calculates change in droplet radius by condensation/evaporation, using
99!> either an analytic formula or by numerically integrating the radius growth
100!> equation including curvature and solution effects using Rosenbrocks method
101!> (see Numerical recipes in FORTRAN, 2nd edition, p. 731).
102!> The analytical formula and growth equation follow those given in
103!> Rogers and Yau (A short course in cloud physics, 3rd edition, p. 102/103).
104!------------------------------------------------------------------------------!
105 SUBROUTINE lpm_droplet_condensation (ip,jp,kp)
106 
107
108    USE arrays_3d,                                                             &
109        ONLY:  hyp, pt, q,  ql_c, ql_v
110
111    USE cloud_parameters,                                                      &
112        ONLY:  l_d_rv, l_v, rho_l, r_v
113
114    USE constants,                                                             &
115        ONLY:  pi
116
117    USE control_parameters,                                                    &
118        ONLY:  dt_3d, dz, message_string, molecular_viscosity, rho_surface
119
120    USE cpulog,                                                                &
121        ONLY:  cpu_log, log_point_s
122
123    USE grid_variables,                                                        &
124        ONLY:  dx, dy
125
126    USE lpm_collision_kernels_mod,                                             &
127        ONLY:  rclass_lbound, rclass_ubound
128
129    USE kinds
130
131    USE particle_attributes,                                                   &
132        ONLY:  curvature_solution_effects, hall_kernel,                        &
133               molecular_weight_of_solute, molecular_weight_of_water,          &
134               number_of_particles, particles, radius_classes, rho_s,          &
135               use_kernel_tables, vanthoff, wang_kernel
136
137
138    IMPLICIT NONE
139
140    INTEGER(iwp) :: ip                         !<
141    INTEGER(iwp) :: internal_timestep_count    !<
142    INTEGER(iwp) :: jp                         !<
143    INTEGER(iwp) :: jtry                       !<
144    INTEGER(iwp) :: kp                         !<
145    INTEGER(iwp) :: n                          !<
146    INTEGER(iwp) :: ros_count                  !<
147 
148    INTEGER(iwp), PARAMETER ::  maxtry = 40    !<
149
150    LOGICAL ::  repeat                         !<
151
152    REAL(wp) ::  aa                            !<
153    REAL(wp) ::  afactor                       !< curvature effects
154    REAL(wp) ::  arg                           !<
155    REAL(wp) ::  bfactor                       !< solute effects
156    REAL(wp) ::  ddenom                        !<
157    REAL(wp) ::  delta_r                       !<
158    REAL(wp) ::  diameter                      !< diameter of cloud droplets
159    REAL(wp) ::  diff_coeff_v                  !< diffusivity for water vapor
160    REAL(wp) ::  drdt                          !<
161    REAL(wp) ::  drdt_ini                      !<
162    REAL(wp) ::  dt_ros                        !<
163    REAL(wp) ::  dt_ros_next                   !<
164    REAL(wp) ::  dt_ros_sum                    !<
165    REAL(wp) ::  dt_ros_sum_ini                !<
166    REAL(wp) ::  d2rdtdr                       !<
167    REAL(wp) ::  errmax                        !<
168    REAL(wp) ::  e_a                           !< current vapor pressure
169    REAL(wp) ::  e_s                           !< current saturation vapor pressure
170    REAL(wp) ::  err_ros                       !<
171    REAL(wp) ::  g1                            !<
172    REAL(wp) ::  g2                            !<
173    REAL(wp) ::  g3                            !<
174    REAL(wp) ::  g4                            !<
175    REAL(wp) ::  r_ros                         !<
176    REAL(wp) ::  r_ros_ini                     !<
177    REAL(wp) ::  sigma                         !<
178    REAL(wp) ::  thermal_conductivity_v        !< thermal conductivity for water
179    REAL(wp) ::  t_int                         !< temperature
180    REAL(wp) ::  w_s                           !< terminal velocity of droplets
181    REAL(wp) ::  re_p                          !<
182
183!
184!-- Parameters for Rosenbrock method
185    REAL(wp), PARAMETER ::  a21 = 2.0_wp               !<
186    REAL(wp), PARAMETER ::  a31 = 48.0_wp / 25.0_wp    !<
187    REAL(wp), PARAMETER ::  a32 = 6.0_wp / 25.0_wp     !<
188    REAL(wp), PARAMETER ::  b1 = 19.0_wp / 9.0_wp      !<
189    REAL(wp), PARAMETER ::  b2 = 0.5_wp                !<
190    REAL(wp), PARAMETER ::  b3 = 25.0_wp / 108.0_wp    !<
191    REAL(wp), PARAMETER ::  b4 = 125.0_wp / 108.0_wp   !<
192    REAL(wp), PARAMETER ::  c21 = -8.0_wp              !<
193    REAL(wp), PARAMETER ::  c31 = 372.0_wp / 25.0_wp   !<
194    REAL(wp), PARAMETER ::  c32 = 12.0_wp / 5.0_wp     !<
195    REAL(wp), PARAMETER ::  c41 = -112.0_wp / 125.0_wp !<
196    REAL(wp), PARAMETER ::  c42 = -54.0_wp / 125.0_wp  !<
197    REAL(wp), PARAMETER ::  c43 = -2.0_wp / 5.0_wp     !<
198    REAL(wp), PARAMETER ::  errcon = 0.1296_wp         !<
199    REAL(wp), PARAMETER ::  e1 = 17.0_wp / 54.0_wp     !<
200    REAL(wp), PARAMETER ::  e2 = 7.0_wp / 36.0_wp      !<
201    REAL(wp), PARAMETER ::  e3 = 0.0_wp                !<
202    REAL(wp), PARAMETER ::  e4 = 125.0_wp / 108.0_wp   !<
203    REAL(wp), PARAMETER ::  eps_ros = 1.0E-3_wp        !< accuracy of Rosenbrock method
204    REAL(wp), PARAMETER ::  gam = 0.5_wp               !<
205    REAL(wp), PARAMETER ::  grow = 1.5_wp              !<
206    REAL(wp), PARAMETER ::  pgrow = -0.25_wp           !<
207    REAL(wp), PARAMETER ::  pshrnk = -1.0_wp /3.0_wp   !<
208    REAL(wp), PARAMETER ::  shrnk = 0.5_wp             !<
209
210!
211!-- Parameters for terminal velocity
212    REAL(wp), PARAMETER ::  a_rog = 9.65_wp      !< parameter for fall velocity
213    REAL(wp), PARAMETER ::  b_rog = 10.43_wp     !< parameter for fall velocity
214    REAL(wp), PARAMETER ::  c_rog = 0.6_wp       !< parameter for fall velocity
215    REAL(wp), PARAMETER ::  k_cap_rog = 4.0_wp   !< parameter for fall velocity
216    REAL(wp), PARAMETER ::  k_low_rog = 12.0_wp  !< parameter for fall velocity
217    REAL(wp), PARAMETER ::  d0_rog = 0.745_wp    !< separation diameter
218
219    REAL(wp), DIMENSION(number_of_particles) ::  ventilation_effect     !<
220    REAL(wp), DIMENSION(number_of_particles) ::  new_r                  !<
221
222
223
224    CALL cpu_log( log_point_s(42), 'lpm_droplet_condens', 'start' )
225
226!
227!-- Calculate temperature, saturation vapor pressure and current vapor pressure
228    t_int = pt(kp,jp,ip) * ( hyp(kp) / 100000.0_wp )**0.286_wp
229    e_s   = 611.0_wp * EXP( l_d_rv * ( 3.6609E-3_wp - 1.0_wp / t_int ) )
230    e_a   = q(kp,jp,ip) * hyp(kp) / ( 0.378_wp * q(kp,jp,ip) + 0.622_wp )
231!
232!-- Thermal conductivity for water (from Rogers and Yau, Table 7.1),
233!-- diffusivity for water vapor (after Hall und Pruppacher, 1976)
234    thermal_conductivity_v = 7.94048E-05_wp * t_int + 0.00227011_wp
235    diff_coeff_v           = 0.211E-4_wp * ( t_int / 273.15_wp )**1.94_wp * &
236                             ( 101325.0_wp / hyp(kp) )
237!
238!-- Calculate effects of heat conductivity and diffusion of water vapor on the
239!-- condensation/evaporation process (typically known as 1.0 / (F_k + F_d) )
240    ddenom  = 1.0_wp / ( rho_l * r_v * t_int / ( e_s * diff_coeff_v ) +        &
241                         ( l_v / ( r_v * t_int ) - 1.0_wp ) * rho_l *          &
242                         l_v / ( thermal_conductivity_v * t_int )              &
243                       )
244
245    new_r = 0.0_wp
246
247!
248!-- Determine ventilation effect on evaporation of large drops
249    DO  n = 1, number_of_particles
250
251       IF ( particles(n)%radius >= 4.0E-5_wp  .AND.  e_a / e_s < 1.0_wp )  THEN
252!
253!--       Terminal velocity is computed for vertical direction (Rogers et al.,
254!--       1993, J. Appl. Meteorol.)
255          diameter = particles(n)%radius * 2000.0_wp !diameter in mm
256          IF ( diameter <= d0_rog )  THEN
257             w_s = k_cap_rog * diameter * ( 1.0_wp - EXP( -k_low_rog * diameter ) )
258          ELSE
259             w_s = a_rog - b_rog * EXP( -c_rog * diameter )
260          ENDIF
261!
262!--       First calculate droplet's Reynolds number
263          re_p = 2.0_wp * particles(n)%radius * w_s / molecular_viscosity
264!
265!--       Ventilation coefficient (Rogers and Yau, 1989):
266          IF ( re_p > 2.5_wp )  THEN
267             ventilation_effect(n) = 0.78_wp + 0.28_wp * SQRT( re_p )
268          ELSE
269             ventilation_effect(n) = 1.0_wp + 0.09_wp * re_p
270          ENDIF
271       ELSE
272!
273!--       For small droplets or in supersaturated environments, the ventilation
274!--       effect does not play a role
275          ventilation_effect(n) = 1.0_wp
276       ENDIF
277    ENDDO
278
279!
280!-- Use analytic model for condensational growth
281    IF( .NOT. curvature_solution_effects ) then
282       DO  n = 1, number_of_particles
283          arg = particles(n)%radius**2 + 2.0_wp * dt_3d * ddenom *             &
284                                         ventilation_effect(n) *               &
285                                         ( e_a / e_s - 1.0_wp )
286          arg = MAX( arg, 1.0E-16_wp )
287          new_r(n) = SQRT( arg )
288       ENDDO
289    ENDIF
290
291!
292!-- If selected, use numerical solution of the condensational growth
293!-- equation (e.g., for studying the activation of aerosols).
294!-- Curvature and solutions effects are included in growth equation.
295!-- Change in Radius is calculated with a 4th-order Rosenbrock method
296!-- for stiff o.d.e's with monitoring local truncation error to adjust
297!-- stepsize (see Numerical recipes in FORTRAN, 2nd edition, p. 731).
298    DO  n = 1, number_of_particles
299       IF ( curvature_solution_effects )  THEN
300
301          ros_count = 0
302          repeat = .TRUE.
303!
304!--       Carry out the Rosenbrock algorithm. In case of unreasonable results
305!--       the switch "repeat" will be set true and the algorithm will be carried
306!--       out again with the internal time step set to its initial (small) value.
307!--       Unreasonable results may occur if the external conditions, especially
308!--       the supersaturation, has significantly changed compared to the last
309!--       PALM timestep.
310          DO WHILE ( repeat )
311
312             repeat = .FALSE.
313!
314!--          Curvature effect (afactor) with surface tension parameterization
315!--          by Straka (2009)
316             sigma = 0.0761_wp - 0.000155_wp * ( t_int - 273.15_wp )
317             afactor = 2.0_wp * sigma / ( rho_l * r_v * t_int )
318!
319!--          Solute effect (bfactor)
320             bfactor = vanthoff * rho_s * particles(n)%rvar2**3 *              &
321                       molecular_weight_of_water /                             &
322                       ( rho_l * molecular_weight_of_solute )
323
324             r_ros = particles(n)%radius
325             dt_ros_sum  = 0.0_wp      ! internal integrated time (s)
326             internal_timestep_count = 0
327!
328!--          Take internal time step values from the end of last PALM time step
329             dt_ros_next = particles(n)%rvar1
330
331!
332!--          Internal time step should not be > 1.0E-2 and < 1.0E-6
333             IF ( dt_ros_next > 1.0E-2_wp )  THEN
334                dt_ros_next = 1.0E-2_wp
335             ELSEIF ( dt_ros_next < 1.0E-6_wp )  THEN
336                dt_ros_next = 1.0E-6_wp
337             ENDIF
338
339!
340!--          If calculation of Rosenbrock method is repeated due to unreasonalble
341!--          results during previous try the initial internal time step has to be
342!--          reduced
343             IF ( ros_count > 1 )  THEN
344                dt_ros_next = dt_ros_next * 0.1_wp
345             ELSEIF ( ros_count > 5 )  THEN
346!
347!--             Prevent creation of infinite loop
348                message_string = 'ros_count > 5 in Rosenbrock method'
349                CALL message( 'lpm_droplet_condensation', 'PA0018', 2, 2, &
350                               0, 6, 0 )
351             ENDIF
352
353!
354!--          Internal time step must not be larger than PALM time step
355             dt_ros = MIN( dt_ros_next, dt_3d )
356
357!
358!--          Integrate growth equation in time unless PALM time step is reached
359             DO WHILE ( dt_ros_sum < dt_3d )
360
361                internal_timestep_count = internal_timestep_count + 1
362
363!
364!--             Derivative at starting value
365                drdt = ddenom * ventilation_effect(n) * ( e_a / e_s - 1.0_wp - &
366                                                          afactor / r_ros +    &
367                                                          bfactor / r_ros**3   &
368                                                        ) / r_ros
369
370                drdt_ini       = drdt
371                dt_ros_sum_ini = dt_ros_sum
372                r_ros_ini      = r_ros
373
374!
375!--             Calculate radial derivative of dr/dt
376                d2rdtdr = ddenom * ventilation_effect(n) *                     &
377                                       ( ( 1.0_wp - e_a / e_s ) / r_ros**2 +   &
378                                         2.0_wp * afactor / r_ros**3 -         &
379                                         4.0_wp * bfactor / r_ros**5           &
380                                       )
381!
382!--             Adjust stepsize unless required accuracy is reached
383                DO  jtry = 1, maxtry+1
384
385                   IF ( jtry == maxtry+1 )  THEN
386                      message_string = 'maxtry > 40 in Rosenbrock method'
387                      CALL message( 'lpm_droplet_condensation', 'PA0347', 0,   &
388                                    1, 0, 6, 0 )
389                   ENDIF
390
391                   aa    = 1.0_wp / ( gam * dt_ros ) - d2rdtdr
392                   g1    = drdt_ini / aa
393                   r_ros = r_ros_ini + a21 * g1
394                   drdt  = ddenom * ventilation_effect(n) * ( e_a / e_s - 1.0_wp - &
395                                                              afactor / r_ros +    &
396                                                              bfactor / r_ros**3   &
397                                                            ) / r_ros
398
399                   g2    = ( drdt + c21 * g1 / dt_ros )&
400                           / aa
401                   r_ros = r_ros_ini + a31 * g1 + a32 * g2
402                   drdt  = ddenom * ventilation_effect(n) * ( e_a / e_s - 1.0_wp - &
403                                                              afactor / r_ros +    &
404                                                              bfactor / r_ros**3   &
405                                                            ) / r_ros
406
407                   g3    = ( drdt +  &
408                             ( c31 * g1 + c32 * g2 ) / dt_ros ) / aa
409                   g4    = ( drdt +  &
410                             ( c41 * g1 + c42 * g2 + c43 * g3 ) / dt_ros ) / aa
411                   r_ros = r_ros_ini + b1 * g1 + b2 * g2 + b3 * g3 + b4 * g4
412
413                   dt_ros_sum = dt_ros_sum_ini + dt_ros
414
415                   IF ( dt_ros_sum == dt_ros_sum_ini )  THEN
416                      message_string = 'zero stepsize in Rosenbrock method'
417                      CALL message( 'lpm_droplet_condensation', 'PA0348', 2,   &
418                                    2, 0, 6, 0 )
419                   ENDIF
420!
421!--                Calculate error
422                   err_ros = e1 * g1 + e2 * g2 + e3 * g3 + e4 * g4
423                   errmax  = 0.0_wp
424                   errmax  = MAX( errmax, ABS( err_ros / r_ros_ini ) ) / eps_ros
425!
426!--                Leave loop if accuracy is sufficient, otherwise try again
427!--                with a reduced stepsize
428                   IF ( errmax <= 1.0_wp )  THEN
429                      EXIT
430                   ELSE
431                      dt_ros = MAX( ABS( 0.9_wp * dt_ros * errmax**pshrnk ),   &
432                                    shrnk * ABS( dt_ros ) )
433                   ENDIF
434
435                ENDDO  ! loop for stepsize adjustment
436
437!
438!--             Calculate next internal time step
439                IF ( errmax > errcon )  THEN
440                   dt_ros_next = 0.9_wp * dt_ros * errmax**pgrow
441                ELSE
442                   dt_ros_next = grow * dt_ros
443                ENDIF
444
445!
446!--             Estimated time step is reduced if the PALM time step is exceeded
447                IF ( ( dt_ros_next + dt_ros_sum ) >= dt_3d )  THEN
448                   dt_ros = dt_3d - dt_ros_sum
449                ELSE
450                   dt_ros = dt_ros_next
451                ENDIF
452
453             ENDDO
454!
455!--          Store internal time step value for next PALM step
456             particles(n)%rvar1 = dt_ros_next
457
458!
459!--          Radius should not fall below 1E-8 because Rosenbrock method may
460!--          lead to errors otherwise
461             new_r(n) = MAX( r_ros, particles(n)%rvar2 )
462!
463!--          Check if calculated droplet radius change is reasonable since in
464!--          case of droplet evaporation the Rosenbrock method may lead to
465!--          secondary solutions which are physically unlikely.
466!--          Due to the solution effect the droplets may grow for relative
467!--          humidities below 100%, but change of radius should not be too
468!--          large. In case of unreasonable droplet growth the Rosenbrock
469!--          method is recalculated using a smaller initial time step.
470!--          Limiting values are tested for droplets down to 1.0E-7
471             IF ( new_r(n) - particles(n)%radius >= 3.0E-7_wp  .AND.  &
472                  e_a / e_s < 0.97_wp )  THEN
473                ros_count = ros_count + 1
474                repeat = .TRUE.
475             ENDIF
476
477          ENDDO    ! Rosenbrock method
478
479       ENDIF
480
481       delta_r = new_r(n) - particles(n)%radius
482
483!
484!--    Sum up the change in volume of liquid water for the respective grid
485!--    volume (this is needed later in lpm_calc_liquid_water_content for
486!--    calculating the release of latent heat)
487       ql_c(kp,jp,ip) = ql_c(kp,jp,ip) + particles(n)%weight_factor *          &
488                                   rho_l * 1.33333333_wp * pi *                &
489                                   ( new_r(n)**3 - particles(n)%radius**3 ) /  &
490                                   ( rho_surface * dx * dy * dz )
491       IF ( ql_c(kp,jp,ip) > 100.0_wp )  THEN
492          WRITE( message_string, * ) 'k=',kp,' j=',jp,' i=',ip,      &
493                       ' ql_c=',ql_c(kp,jp,ip), ' &part(',n,')%wf=', &
494                       particles(n)%weight_factor,' delta_r=',delta_r
495          CALL message( 'lpm_droplet_condensation', 'PA0143', 2, 2, -1, 6, 1 )
496       ENDIF
497
498!
499!--    Change the droplet radius
500       IF ( ( new_r(n) - particles(n)%radius ) < 0.0_wp  .AND.        &
501            new_r(n) < 0.0_wp )  THEN
502          WRITE( message_string, * ) '#1 k=',kp,' j=',jp,' i=',ip,    &
503                       ' e_s=',e_s, ' e_a=',e_a,' t_int=',t_int,      &
504                       ' &delta_r=',delta_r,                          &
505                       ' particle_radius=',particles(n)%radius
506          CALL message( 'lpm_droplet_condensation', 'PA0144', 2, 2, -1, 6, 1 )
507       ENDIF
508
509!
510!--    Sum up the total volume of liquid water (needed below for
511!--    re-calculating the weighting factors)
512       ql_v(kp,jp,ip) = ql_v(kp,jp,ip) + particles(n)%weight_factor * new_r(n)**3
513
514       particles(n)%radius = new_r(n)
515
516!
517!--    Determine radius class of the particle needed for collision
518       IF ( ( hall_kernel  .OR.  wang_kernel )  .AND.  use_kernel_tables )     &
519       THEN
520          particles(n)%class = ( LOG( new_r(n) ) - rclass_lbound ) /           &
521                               ( rclass_ubound - rclass_lbound ) *             &
522                               radius_classes
523          particles(n)%class = MIN( particles(n)%class, radius_classes )
524          particles(n)%class = MAX( particles(n)%class, 1 )
525       ENDIF
526
527    ENDDO
528
529    CALL cpu_log( log_point_s(42), 'lpm_droplet_condens', 'stop' )
530
531
532 END SUBROUTINE lpm_droplet_condensation
Note: See TracBrowser for help on using the repository browser.