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

Last change on this file since 1883 was 1872, checked in by hoffmann, 8 years ago

last commit documented

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