[1682] | 1 | !> @file lpm_droplet_condensation.f90 |
---|
[2000] | 2 | !------------------------------------------------------------------------------! |
---|
[2696] | 3 | ! This file is part of the PALM model system. |
---|
[1036] | 4 | ! |
---|
[2000] | 5 | ! PALM is free software: you can redistribute it and/or modify it under the |
---|
| 6 | ! terms of the GNU General Public License as published by the Free Software |
---|
| 7 | ! Foundation, either version 3 of the License, or (at your option) any later |
---|
| 8 | ! version. |
---|
[1036] | 9 | ! |
---|
| 10 | ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY |
---|
| 11 | ! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR |
---|
| 12 | ! A PARTICULAR PURPOSE. See the GNU General Public License for more details. |
---|
| 13 | ! |
---|
| 14 | ! You should have received a copy of the GNU General Public License along with |
---|
| 15 | ! PALM. If not, see <http://www.gnu.org/licenses/>. |
---|
| 16 | ! |
---|
[2718] | 17 | ! Copyright 1997-2018 Leibniz Universitaet Hannover |
---|
[2000] | 18 | !------------------------------------------------------------------------------! |
---|
[1036] | 19 | ! |
---|
[849] | 20 | ! Current revisions: |
---|
| 21 | ! ------------------ |
---|
[2375] | 22 | ! |
---|
| 23 | ! |
---|
[1891] | 24 | ! Former revisions: |
---|
| 25 | ! ----------------- |
---|
| 26 | ! $Id: lpm_droplet_condensation.f90 2718 2018-01-02 08:49:38Z suehring $ |
---|
[2716] | 27 | ! Corrected "Former revisions" section |
---|
| 28 | ! |
---|
| 29 | ! 2696 2017-12-14 17:12:51Z kanani |
---|
| 30 | ! Change in file header (GPL part) |
---|
| 31 | ! |
---|
| 32 | ! 2608 2017-11-13 14:04:26Z schwenkel |
---|
[2608] | 33 | ! Calculation of magnus equation in external module (diagnostic_quantities_mod). |
---|
| 34 | ! |
---|
| 35 | ! 2375 2017-08-29 14:10:28Z schwenkel |
---|
[2375] | 36 | ! Changed ONLY-dependencies |
---|
| 37 | ! |
---|
| 38 | ! 2312 2017-07-14 20:26:51Z hoffmann |
---|
[2312] | 39 | ! Rosenbrock scheme improved. Gas-kinetic effect added. |
---|
[1891] | 40 | ! |
---|
[2001] | 41 | ! 2000 2016-08-20 18:09:15Z knoop |
---|
| 42 | ! Forced header and separation lines into 80 columns |
---|
[2312] | 43 | ! |
---|
[1891] | 44 | ! 1890 2016-04-22 08:52:11Z hoffmann |
---|
[1890] | 45 | ! Some improvements of the Rosenbrock method. If the Rosenbrock method needs more |
---|
| 46 | ! than 40 iterations to find a sufficient time setp, the model is not aborted. |
---|
[2312] | 47 | ! This might lead to small erros. But they can be assumend as negligible, since |
---|
| 48 | ! the maximum timesetp of the Rosenbrock method after 40 iterations will be |
---|
| 49 | ! smaller than 10^-16 s. |
---|
| 50 | ! |
---|
[1872] | 51 | ! 1871 2016-04-15 11:46:09Z hoffmann |
---|
| 52 | ! Initialization of aerosols added. |
---|
| 53 | ! |
---|
[1851] | 54 | ! 1849 2016-04-08 11:33:18Z hoffmann |
---|
[1852] | 55 | ! Interpolation of supersaturation has been removed because it is not in |
---|
| 56 | ! accordance with the release/depletion of latent heat/water vapor in |
---|
[1849] | 57 | ! interaction_droplets_ptq. |
---|
| 58 | ! Calculation of particle Reynolds number has been corrected. |
---|
[1852] | 59 | ! eps_ros added from modules. |
---|
[1849] | 60 | ! |
---|
[1832] | 61 | ! 1831 2016-04-07 13:15:51Z hoffmann |
---|
| 62 | ! curvature_solution_effects moved to particle_attributes |
---|
| 63 | ! |
---|
[1823] | 64 | ! 1822 2016-04-07 07:49:42Z hoffmann |
---|
| 65 | ! Unused variables removed. |
---|
| 66 | ! |
---|
[1683] | 67 | ! 1682 2015-10-07 23:56:08Z knoop |
---|
[2312] | 68 | ! Code annotations made doxygen readable |
---|
| 69 | ! |
---|
[1360] | 70 | ! 1359 2014-04-11 17:15:14Z hoffmann |
---|
[2312] | 71 | ! New particle structure integrated. |
---|
[1360] | 72 | ! Kind definition added to all floating point numbers. |
---|
| 73 | ! |
---|
[1347] | 74 | ! 1346 2014-03-27 13:18:20Z heinze |
---|
[2312] | 75 | ! Bugfix: REAL constants provided with KIND-attribute especially in call of |
---|
[1347] | 76 | ! intrinsic function like MAX, MIN, SIGN |
---|
| 77 | ! |
---|
[1323] | 78 | ! 1322 2014-03-20 16:38:49Z raasch |
---|
| 79 | ! REAL constants defined as wp-kind |
---|
| 80 | ! |
---|
[1321] | 81 | ! 1320 2014-03-20 08:40:49Z raasch |
---|
[1320] | 82 | ! ONLY-attribute added to USE-statements, |
---|
| 83 | ! kind-parameters added to all INTEGER and REAL declaration statements, |
---|
| 84 | ! kinds are defined in new module kinds, |
---|
| 85 | ! comment fields (!:) to be used for variable explanations added to |
---|
| 86 | ! all variable declaration statements |
---|
[1072] | 87 | ! |
---|
[1319] | 88 | ! 1318 2014-03-17 13:35:16Z raasch |
---|
| 89 | ! module interfaces removed |
---|
| 90 | ! |
---|
[1093] | 91 | ! 1092 2013-02-02 11:24:22Z raasch |
---|
| 92 | ! unused variables removed |
---|
| 93 | ! |
---|
[1072] | 94 | ! 1071 2012-11-29 16:54:55Z franke |
---|
[1071] | 95 | ! Ventilation effect for evaporation of large droplets included |
---|
| 96 | ! Check for unreasonable results included in calculation of Rosenbrock method |
---|
| 97 | ! since physically unlikely results were observed and for the same |
---|
| 98 | ! reason the first internal time step in Rosenbrock method should be < 1.0E02 in |
---|
| 99 | ! case of evaporation |
---|
| 100 | ! Unnecessary calculation of ql_int removed |
---|
| 101 | ! Unnecessary calculations in Rosenbrock method (d2rdt2, drdt_m, dt_ros_last) |
---|
| 102 | ! removed |
---|
| 103 | ! Bugfix: factor in calculation of surface tension changed from 0.00155 to |
---|
| 104 | ! 0.000155 |
---|
[849] | 105 | ! |
---|
[1037] | 106 | ! 1036 2012-10-22 13:43:42Z raasch |
---|
| 107 | ! code put under GPL (PALM 3.9) |
---|
| 108 | ! |
---|
[850] | 109 | ! 849 2012-03-15 10:35:09Z raasch |
---|
| 110 | ! initial revision (former part of advec_particles) |
---|
[849] | 111 | ! |
---|
[850] | 112 | ! |
---|
[849] | 113 | ! Description: |
---|
| 114 | ! ------------ |
---|
[1682] | 115 | !> Calculates change in droplet radius by condensation/evaporation, using |
---|
| 116 | !> either an analytic formula or by numerically integrating the radius growth |
---|
| 117 | !> equation including curvature and solution effects using Rosenbrocks method |
---|
| 118 | !> (see Numerical recipes in FORTRAN, 2nd edition, p. 731). |
---|
| 119 | !> The analytical formula and growth equation follow those given in |
---|
| 120 | !> Rogers and Yau (A short course in cloud physics, 3rd edition, p. 102/103). |
---|
[849] | 121 | !------------------------------------------------------------------------------! |
---|
[1682] | 122 | SUBROUTINE lpm_droplet_condensation (ip,jp,kp) |
---|
[849] | 123 | |
---|
[2312] | 124 | |
---|
[1320] | 125 | USE arrays_3d, & |
---|
[2312] | 126 | ONLY: hyp, pt, q, ql_c, ql_v |
---|
[849] | 127 | |
---|
[1320] | 128 | USE cloud_parameters, & |
---|
[2375] | 129 | ONLY: l_d_rv, l_v, molecular_weight_of_solute, & |
---|
| 130 | molecular_weight_of_water, rho_l, rho_s, r_v, vanthoff |
---|
[849] | 131 | |
---|
[1320] | 132 | USE constants, & |
---|
| 133 | ONLY: pi |
---|
[849] | 134 | |
---|
[1320] | 135 | USE control_parameters, & |
---|
[1822] | 136 | ONLY: dt_3d, dz, message_string, molecular_viscosity, rho_surface |
---|
| 137 | |
---|
[1320] | 138 | USE cpulog, & |
---|
| 139 | ONLY: cpu_log, log_point_s |
---|
[849] | 140 | |
---|
[2608] | 141 | USE diagnostic_quantities_mod, & |
---|
| 142 | ONLY: magnus |
---|
| 143 | |
---|
[1320] | 144 | USE grid_variables, & |
---|
[1822] | 145 | ONLY: dx, dy |
---|
[1071] | 146 | |
---|
[1320] | 147 | USE lpm_collision_kernels_mod, & |
---|
| 148 | ONLY: rclass_lbound, rclass_ubound |
---|
[849] | 149 | |
---|
[1320] | 150 | USE kinds |
---|
| 151 | |
---|
| 152 | USE particle_attributes, & |
---|
[2375] | 153 | ONLY: curvature_solution_effects, hall_kernel, number_of_particles, & |
---|
| 154 | particles, radius_classes, use_kernel_tables, wang_kernel |
---|
[1320] | 155 | |
---|
| 156 | IMPLICIT NONE |
---|
| 157 | |
---|
[1682] | 158 | INTEGER(iwp) :: ip !< |
---|
| 159 | INTEGER(iwp) :: jp !< |
---|
| 160 | INTEGER(iwp) :: kp !< |
---|
| 161 | INTEGER(iwp) :: n !< |
---|
[1320] | 162 | |
---|
[1849] | 163 | REAL(wp) :: afactor !< curvature effects |
---|
[1682] | 164 | REAL(wp) :: arg !< |
---|
[1849] | 165 | REAL(wp) :: bfactor !< solute effects |
---|
[1682] | 166 | REAL(wp) :: ddenom !< |
---|
| 167 | REAL(wp) :: delta_r !< |
---|
[1849] | 168 | REAL(wp) :: diameter !< diameter of cloud droplets |
---|
[2312] | 169 | REAL(wp) :: diff_coeff !< diffusivity for water vapor |
---|
[1682] | 170 | REAL(wp) :: drdt !< |
---|
| 171 | REAL(wp) :: dt_ros !< |
---|
| 172 | REAL(wp) :: dt_ros_sum !< |
---|
| 173 | REAL(wp) :: d2rdtdr !< |
---|
[1849] | 174 | REAL(wp) :: e_a !< current vapor pressure |
---|
| 175 | REAL(wp) :: e_s !< current saturation vapor pressure |
---|
[2312] | 176 | REAL(wp) :: error !< local truncation error in Rosenbrock |
---|
| 177 | REAL(wp) :: k1 !< |
---|
| 178 | REAL(wp) :: k2 !< |
---|
| 179 | REAL(wp) :: r_err !< First order estimate of Rosenbrock radius |
---|
| 180 | REAL(wp) :: r_ros !< Rosenbrock radius |
---|
| 181 | REAL(wp) :: r_ros_ini !< initial Rosenbrock radius |
---|
| 182 | REAL(wp) :: r0 !< gas-kinetic lengthscale |
---|
| 183 | REAL(wp) :: sigma !< surface tension of water |
---|
| 184 | REAL(wp) :: thermal_conductivity !< thermal conductivity for water |
---|
[1849] | 185 | REAL(wp) :: t_int !< temperature |
---|
| 186 | REAL(wp) :: w_s !< terminal velocity of droplets |
---|
[2312] | 187 | REAL(wp) :: re_p !< particle Reynolds number |
---|
[1849] | 188 | ! |
---|
[2312] | 189 | !-- Parameters for Rosenbrock method (see Verwer et al., 1999) |
---|
| 190 | REAL(wp), PARAMETER :: prec = 1.0E-3_wp !< precision of Rosenbrock solution |
---|
| 191 | REAL(wp), PARAMETER :: q_increase = 1.5_wp !< increase factor in timestep |
---|
| 192 | REAL(wp), PARAMETER :: q_decrease = 0.9_wp !< decrease factor in timestep |
---|
| 193 | REAL(wp), PARAMETER :: gamma = 0.292893218814_wp !< = 1.0 - 1.0 / SQRT(2.0) |
---|
[849] | 194 | ! |
---|
[1849] | 195 | !-- Parameters for terminal velocity |
---|
| 196 | REAL(wp), PARAMETER :: a_rog = 9.65_wp !< parameter for fall velocity |
---|
| 197 | REAL(wp), PARAMETER :: b_rog = 10.43_wp !< parameter for fall velocity |
---|
| 198 | REAL(wp), PARAMETER :: c_rog = 0.6_wp !< parameter for fall velocity |
---|
| 199 | REAL(wp), PARAMETER :: k_cap_rog = 4.0_wp !< parameter for fall velocity |
---|
| 200 | REAL(wp), PARAMETER :: k_low_rog = 12.0_wp !< parameter for fall velocity |
---|
| 201 | REAL(wp), PARAMETER :: d0_rog = 0.745_wp !< separation diameter |
---|
[849] | 202 | |
---|
[1849] | 203 | REAL(wp), DIMENSION(number_of_particles) :: ventilation_effect !< |
---|
| 204 | REAL(wp), DIMENSION(number_of_particles) :: new_r !< |
---|
[849] | 205 | |
---|
[1849] | 206 | CALL cpu_log( log_point_s(42), 'lpm_droplet_condens', 'start' ) |
---|
[849] | 207 | |
---|
| 208 | ! |
---|
[2312] | 209 | !-- Absolute temperature |
---|
[1849] | 210 | t_int = pt(kp,jp,ip) * ( hyp(kp) / 100000.0_wp )**0.286_wp |
---|
[849] | 211 | ! |
---|
[2312] | 212 | !-- Saturation vapor pressure (Eq. 10 in Bolton, 1980) |
---|
[2608] | 213 | e_s = magnus( t_int ) |
---|
[1849] | 214 | ! |
---|
[2312] | 215 | !-- Current vapor pressure |
---|
| 216 | e_a = q(kp,jp,ip) * hyp(kp) / ( q(kp,jp,ip) + 0.622_wp ) |
---|
| 217 | ! |
---|
| 218 | !-- Thermal conductivity for water (from Rogers and Yau, Table 7.1) |
---|
| 219 | thermal_conductivity = 7.94048E-05_wp * t_int + 0.00227011_wp |
---|
| 220 | ! |
---|
| 221 | !-- Moldecular diffusivity of water vapor in air (Hall und Pruppacher, 1976) |
---|
| 222 | diff_coeff = 0.211E-4_wp * ( t_int / 273.15_wp )**1.94_wp * & |
---|
| 223 | ( 101325.0_wp / hyp(kp) ) |
---|
| 224 | ! |
---|
| 225 | !-- Lengthscale for gas-kinetic effects (from Mordy, 1959, p. 23): |
---|
| 226 | r0 = diff_coeff / 0.036_wp * SQRT( 2.0_wp * pi / ( r_v * t_int ) ) |
---|
| 227 | ! |
---|
| 228 | !-- Calculate effects of heat conductivity and diffusion of water vapor on the |
---|
| 229 | !-- diffusional growth process (usually known as 1.0 / (F_k + F_d) ) |
---|
| 230 | ddenom = 1.0_wp / ( rho_l * r_v * t_int / ( e_s * diff_coeff ) + & |
---|
[1849] | 231 | ( l_v / ( r_v * t_int ) - 1.0_wp ) * rho_l * & |
---|
[2312] | 232 | l_v / ( thermal_conductivity * t_int ) & |
---|
[1849] | 233 | ) |
---|
[1359] | 234 | new_r = 0.0_wp |
---|
[1849] | 235 | ! |
---|
| 236 | !-- Determine ventilation effect on evaporation of large drops |
---|
[1359] | 237 | DO n = 1, number_of_particles |
---|
[1849] | 238 | |
---|
| 239 | IF ( particles(n)%radius >= 4.0E-5_wp .AND. e_a / e_s < 1.0_wp ) THEN |
---|
[849] | 240 | ! |
---|
[1849] | 241 | !-- Terminal velocity is computed for vertical direction (Rogers et al., |
---|
| 242 | !-- 1993, J. Appl. Meteorol.) |
---|
| 243 | diameter = particles(n)%radius * 2000.0_wp !diameter in mm |
---|
| 244 | IF ( diameter <= d0_rog ) THEN |
---|
| 245 | w_s = k_cap_rog * diameter * ( 1.0_wp - EXP( -k_low_rog * diameter ) ) |
---|
| 246 | ELSE |
---|
| 247 | w_s = a_rog - b_rog * EXP( -c_rog * diameter ) |
---|
| 248 | ENDIF |
---|
[849] | 249 | ! |
---|
[2312] | 250 | !-- Calculate droplet's Reynolds number |
---|
[1849] | 251 | re_p = 2.0_wp * particles(n)%radius * w_s / molecular_viscosity |
---|
[1071] | 252 | ! |
---|
[1359] | 253 | !-- Ventilation coefficient (Rogers and Yau, 1989): |
---|
| 254 | IF ( re_p > 2.5_wp ) THEN |
---|
[1849] | 255 | ventilation_effect(n) = 0.78_wp + 0.28_wp * SQRT( re_p ) |
---|
[1071] | 256 | ELSE |
---|
[1849] | 257 | ventilation_effect(n) = 1.0_wp + 0.09_wp * re_p |
---|
[1071] | 258 | ENDIF |
---|
[1849] | 259 | ELSE |
---|
[1071] | 260 | ! |
---|
[2312] | 261 | !-- For small droplets or in supersaturated environments, the ventilation |
---|
[1849] | 262 | !-- effect does not play a role |
---|
| 263 | ventilation_effect(n) = 1.0_wp |
---|
[849] | 264 | ENDIF |
---|
[1359] | 265 | ENDDO |
---|
[849] | 266 | |
---|
[2312] | 267 | IF( .NOT. curvature_solution_effects ) then |
---|
[849] | 268 | ! |
---|
[2312] | 269 | !-- Use analytic model for diffusional growth including gas-kinetic |
---|
| 270 | !-- effects (Mordy, 1959) but without the impact of aerosols. |
---|
[1849] | 271 | DO n = 1, number_of_particles |
---|
[2312] | 272 | arg = ( particles(n)%radius + r0 )**2 + 2.0_wp * dt_3d * ddenom * & |
---|
| 273 | ventilation_effect(n) * & |
---|
| 274 | ( e_a / e_s - 1.0_wp ) |
---|
| 275 | arg = MAX( arg, ( 0.01E-6 + r0 )**2 ) |
---|
| 276 | new_r(n) = SQRT( arg ) - r0 |
---|
[1849] | 277 | ENDDO |
---|
[1359] | 278 | |
---|
[2312] | 279 | ELSE |
---|
[1849] | 280 | ! |
---|
[2312] | 281 | !-- Integrate the diffusional growth including gas-kinetic (Mordy, 1959), |
---|
| 282 | !-- as well as curvature and solute effects (e.g., Köhler, 1936). |
---|
[849] | 283 | ! |
---|
[2312] | 284 | !-- Curvature effect (afactor) with surface tension (sigma) by Straka (2009) |
---|
| 285 | sigma = 0.0761_wp - 0.000155_wp * ( t_int - 273.15_wp ) |
---|
[1071] | 286 | ! |
---|
[2312] | 287 | !-- Solute effect (afactor) |
---|
| 288 | afactor = 2.0_wp * sigma / ( rho_l * r_v * t_int ) |
---|
[849] | 289 | |
---|
[2312] | 290 | DO n = 1, number_of_particles |
---|
[849] | 291 | ! |
---|
[2312] | 292 | !-- Solute effect (bfactor) |
---|
| 293 | bfactor = vanthoff * rho_s * particles(n)%aux1**3 * & |
---|
| 294 | molecular_weight_of_water / ( rho_l * molecular_weight_of_solute ) |
---|
[1071] | 295 | |
---|
[2312] | 296 | dt_ros = particles(n)%aux2 ! use previously stored Rosenbrock timestep |
---|
| 297 | dt_ros_sum = 0.0_wp |
---|
[1871] | 298 | |
---|
[2312] | 299 | r_ros = particles(n)%radius ! initialize Rosenbrock particle radius |
---|
| 300 | r_ros_ini = r_ros |
---|
[849] | 301 | ! |
---|
[2312] | 302 | !-- Integrate growth equation using a 2nd-order Rosenbrock method |
---|
| 303 | !-- (see Verwer et al., 1999, Eq. (3.2)). The Rosenbrock method adjusts |
---|
| 304 | !-- its with internal timestep to minimize the local truncation error. |
---|
| 305 | DO WHILE ( dt_ros_sum < dt_3d ) |
---|
[1071] | 306 | |
---|
[2312] | 307 | dt_ros = MIN( dt_ros, dt_3d - dt_ros_sum ) |
---|
[1871] | 308 | |
---|
[2312] | 309 | DO |
---|
[849] | 310 | |
---|
[2312] | 311 | drdt = ddenom * ventilation_effect(n) * ( e_a / e_s - 1.0 - & |
---|
[1849] | 312 | afactor / r_ros + & |
---|
| 313 | bfactor / r_ros**3 & |
---|
[2312] | 314 | ) / ( r_ros + r0 ) |
---|
[1849] | 315 | |
---|
[2312] | 316 | d2rdtdr = -ddenom * ventilation_effect(n) * ( & |
---|
| 317 | (e_a / e_s - 1.0) * r_ros**4 - & |
---|
| 318 | afactor * r0 * r_ros**2 - & |
---|
| 319 | 2.0 * afactor * r_ros**3 + & |
---|
| 320 | 3.0 * bfactor * r0 + & |
---|
| 321 | 4.0 * bfactor * r_ros & |
---|
| 322 | ) & |
---|
| 323 | / ( r_ros**4 * ( r_ros + r0 )**2 ) |
---|
[849] | 324 | |
---|
[2312] | 325 | k1 = drdt / ( 1.0 - gamma * dt_ros * d2rdtdr ) |
---|
[849] | 326 | |
---|
[2312] | 327 | r_ros = MAX(r_ros_ini + k1 * dt_ros, particles(n)%aux1) |
---|
| 328 | r_err = r_ros |
---|
[849] | 329 | |
---|
[2312] | 330 | drdt = ddenom * ventilation_effect(n) * ( e_a / e_s - 1.0 - & |
---|
| 331 | afactor / r_ros + & |
---|
| 332 | bfactor / r_ros**3 & |
---|
| 333 | ) / ( r_ros + r0 ) |
---|
[849] | 334 | |
---|
[2312] | 335 | k2 = ( drdt - dt_ros * 2.0 * gamma * d2rdtdr * k1 ) / & |
---|
| 336 | ( 1.0 - dt_ros * gamma * d2rdtdr ) |
---|
[849] | 337 | |
---|
[2312] | 338 | r_ros = MAX(r_ros_ini + dt_ros * ( 1.5 * k1 + 0.5 * k2), particles(n)%aux1) |
---|
| 339 | ! |
---|
| 340 | !-- Check error of the solution, and reduce dt_ros if necessary. |
---|
| 341 | error = ABS(r_err - r_ros) / r_ros |
---|
| 342 | IF ( error .GT. prec ) THEN |
---|
| 343 | dt_ros = SQRT( q_decrease * prec / error ) * dt_ros |
---|
| 344 | r_ros = r_ros_ini |
---|
[849] | 345 | ELSE |
---|
[2312] | 346 | dt_ros_sum = dt_ros_sum + dt_ros |
---|
| 347 | dt_ros = q_increase * dt_ros |
---|
| 348 | r_ros_ini = r_ros |
---|
| 349 | EXIT |
---|
[849] | 350 | ENDIF |
---|
| 351 | |
---|
[2312] | 352 | END DO |
---|
[849] | 353 | |
---|
[2312] | 354 | END DO !Rosenbrock loop |
---|
[849] | 355 | ! |
---|
[2312] | 356 | !-- Store new particle radius |
---|
| 357 | new_r(n) = r_ros |
---|
[849] | 358 | ! |
---|
[2312] | 359 | !-- Store internal time step value for next PALM step |
---|
| 360 | particles(n)%aux2 = dt_ros |
---|
[849] | 361 | |
---|
[2312] | 362 | ENDDO !Particle loop |
---|
[849] | 363 | |
---|
[2312] | 364 | ENDIF |
---|
[849] | 365 | |
---|
[2312] | 366 | DO n = 1, number_of_particles |
---|
[849] | 367 | ! |
---|
[2312] | 368 | !-- Sum up the change in liquid water for the respective grid |
---|
| 369 | !-- box for the computation of the release/depletion of water vapor |
---|
| 370 | !-- and heat. |
---|
[1890] | 371 | ql_c(kp,jp,ip) = ql_c(kp,jp,ip) + particles(n)%weight_factor * & |
---|
[1359] | 372 | rho_l * 1.33333333_wp * pi * & |
---|
| 373 | ( new_r(n)**3 - particles(n)%radius**3 ) / & |
---|
[849] | 374 | ( rho_surface * dx * dy * dz ) |
---|
[2312] | 375 | ! |
---|
| 376 | !-- Check if the increase in liqid water is not too big. If this is the case, |
---|
| 377 | !-- the model timestep might be too long. |
---|
[1890] | 378 | IF ( ql_c(kp,jp,ip) > 100.0_wp ) THEN |
---|
| 379 | WRITE( message_string, * ) 'k=',kp,' j=',jp,' i=',ip, & |
---|
| 380 | ' ql_c=',ql_c(kp,jp,ip), ' &part(',n,')%wf=', & |
---|
[849] | 381 | particles(n)%weight_factor,' delta_r=',delta_r |
---|
| 382 | CALL message( 'lpm_droplet_condensation', 'PA0143', 2, 2, -1, 6, 1 ) |
---|
| 383 | ENDIF |
---|
| 384 | ! |
---|
[2312] | 385 | !-- Check if the change in the droplet radius is not too big. If this is the |
---|
| 386 | !-- case, the model timestep might be too long. |
---|
| 387 | delta_r = new_r(n) - particles(n)%radius |
---|
| 388 | IF ( delta_r < 0.0_wp .AND. new_r(n) < 0.0_wp ) THEN |
---|
[1890] | 389 | WRITE( message_string, * ) '#1 k=',kp,' j=',jp,' i=',ip, & |
---|
[1849] | 390 | ' e_s=',e_s, ' e_a=',e_a,' t_int=',t_int, & |
---|
[1890] | 391 | ' &delta_r=',delta_r, & |
---|
[849] | 392 | ' particle_radius=',particles(n)%radius |
---|
| 393 | CALL message( 'lpm_droplet_condensation', 'PA0144', 2, 2, -1, 6, 1 ) |
---|
| 394 | ENDIF |
---|
| 395 | ! |
---|
| 396 | !-- Sum up the total volume of liquid water (needed below for |
---|
| 397 | !-- re-calculating the weighting factors) |
---|
[1890] | 398 | ql_v(kp,jp,ip) = ql_v(kp,jp,ip) + particles(n)%weight_factor * new_r(n)**3 |
---|
[849] | 399 | ! |
---|
| 400 | !-- Determine radius class of the particle needed for collision |
---|
[2312] | 401 | IF ( use_kernel_tables ) THEN |
---|
[1359] | 402 | particles(n)%class = ( LOG( new_r(n) ) - rclass_lbound ) / & |
---|
| 403 | ( rclass_ubound - rclass_lbound ) * & |
---|
[849] | 404 | radius_classes |
---|
| 405 | particles(n)%class = MIN( particles(n)%class, radius_classes ) |
---|
| 406 | particles(n)%class = MAX( particles(n)%class, 1 ) |
---|
| 407 | ENDIF |
---|
[2312] | 408 | ! |
---|
| 409 | !-- Store new radius to particle features |
---|
| 410 | particles(n)%radius = new_r(n) |
---|
[849] | 411 | |
---|
| 412 | ENDDO |
---|
| 413 | |
---|
| 414 | CALL cpu_log( log_point_s(42), 'lpm_droplet_condens', 'stop' ) |
---|
| 415 | |
---|
| 416 | |
---|
| 417 | END SUBROUTINE lpm_droplet_condensation |
---|