[1682] | 1 | !> @file lpm_droplet_collision.f90 |
---|
[1036] | 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 | ! |
---|
[1818] | 16 | ! Copyright 1997-2016 Leibniz Universitaet Hannover |
---|
[1036] | 17 | !--------------------------------------------------------------------------------! |
---|
| 18 | ! |
---|
[849] | 19 | ! Current revisions: |
---|
| 20 | ! ------------------ |
---|
[1823] | 21 | ! |
---|
[1885] | 22 | ! |
---|
[1321] | 23 | ! Former revisions: |
---|
| 24 | ! ----------------- |
---|
| 25 | ! $Id: lpm_droplet_collision.f90 1885 2016-04-21 11:12:26Z gronemeier $ |
---|
| 26 | ! |
---|
[1885] | 27 | ! 1884 2016-04-21 11:11:40Z hoffmann |
---|
| 28 | ! Conservation of mass should only be checked if collisions took place. |
---|
| 29 | ! |
---|
[1861] | 30 | ! 1860 2016-04-13 13:21:28Z hoffmann |
---|
| 31 | ! Interpolation of dissipation rate adjusted to more reasonable values. |
---|
| 32 | ! |
---|
[1823] | 33 | ! 1822 2016-04-07 07:49:42Z hoffmann |
---|
| 34 | ! Integration of a new collision algortithm based on Shima et al. (2009) and |
---|
| 35 | ! Soelch and Kaercher (2010) called all_or_nothing. The previous implemented |
---|
| 36 | ! collision algorithm is called average_impact. Moreover, both algorithms are |
---|
| 37 | ! now positive definit due to their construction, i.e., no negative weighting |
---|
| 38 | ! factors should occur. |
---|
| 39 | ! |
---|
[1683] | 40 | ! 1682 2015-10-07 23:56:08Z knoop |
---|
| 41 | ! Code annotations made doxygen readable |
---|
| 42 | ! |
---|
[1360] | 43 | ! 1359 2014-04-11 17:15:14Z hoffmann |
---|
| 44 | ! New particle structure integrated. |
---|
| 45 | ! Kind definition added to all floating point numbers. |
---|
| 46 | ! |
---|
[1323] | 47 | ! 1322 2014-03-20 16:38:49Z raasch |
---|
| 48 | ! REAL constants defined as wp_kind |
---|
| 49 | ! |
---|
[1321] | 50 | ! 1320 2014-03-20 08:40:49Z raasch |
---|
[1320] | 51 | ! ONLY-attribute added to USE-statements, |
---|
| 52 | ! kind-parameters added to all INTEGER and REAL declaration statements, |
---|
| 53 | ! kinds are defined in new module kinds, |
---|
| 54 | ! revision history before 2012 removed, |
---|
| 55 | ! comment fields (!:) to be used for variable explanations added to |
---|
| 56 | ! all variable declaration statements |
---|
[1072] | 57 | ! |
---|
[1093] | 58 | ! 1092 2013-02-02 11:24:22Z raasch |
---|
| 59 | ! unused variables removed |
---|
| 60 | ! |
---|
[1072] | 61 | ! 1071 2012-11-29 16:54:55Z franke |
---|
[1071] | 62 | ! Calculation of Hall and Wang kernel now uses collision-coalescence formulation |
---|
| 63 | ! proposed by Wang instead of the continuous collection equation (for more |
---|
| 64 | ! information about new method see PALM documentation) |
---|
| 65 | ! Bugfix: message identifiers added |
---|
[849] | 66 | ! |
---|
[1037] | 67 | ! 1036 2012-10-22 13:43:42Z raasch |
---|
| 68 | ! code put under GPL (PALM 3.9) |
---|
| 69 | ! |
---|
[850] | 70 | ! 849 2012-03-15 10:35:09Z raasch |
---|
| 71 | ! initial revision (former part of advec_particles) |
---|
[849] | 72 | ! |
---|
[850] | 73 | ! |
---|
[849] | 74 | ! Description: |
---|
| 75 | ! ------------ |
---|
[1682] | 76 | !> Calculates change in droplet radius by collision. Droplet collision is |
---|
| 77 | !> calculated for each grid box seperately. Collision is parameterized by |
---|
[1822] | 78 | !> using collision kernels. Two different kernels are available: |
---|
[1682] | 79 | !> Hall kernel: Kernel from Hall (1980, J. Atmos. Sci., 2486-2507), which |
---|
| 80 | !> considers collision due to pure gravitational effects. |
---|
| 81 | !> Wang kernel: Beside gravitational effects (treated with the Hall-kernel) also |
---|
| 82 | !> the effects of turbulence on the collision are considered using |
---|
| 83 | !> parameterizations of Ayala et al. (2008, New J. Phys., 10, |
---|
| 84 | !> 075015) and Wang and Grabowski (2009, Atmos. Sci. Lett., 10, |
---|
| 85 | !> 1-8). This kernel includes three possible effects of turbulence: |
---|
| 86 | !> the modification of the relative velocity between the droplets, |
---|
| 87 | !> the effect of preferential concentration, and the enhancement of |
---|
| 88 | !> collision efficiencies. |
---|
[849] | 89 | !------------------------------------------------------------------------------! |
---|
[1682] | 90 | SUBROUTINE lpm_droplet_collision (i,j,k) |
---|
| 91 | |
---|
[849] | 92 | |
---|
[1359] | 93 | |
---|
[1320] | 94 | USE arrays_3d, & |
---|
[1822] | 95 | ONLY: diss, ql_v, ql_vp |
---|
[849] | 96 | |
---|
[1320] | 97 | USE cloud_parameters, & |
---|
[1822] | 98 | ONLY: rho_l |
---|
[1320] | 99 | |
---|
| 100 | USE constants, & |
---|
| 101 | ONLY: pi |
---|
| 102 | |
---|
| 103 | USE control_parameters, & |
---|
[1822] | 104 | ONLY: dt_3d, message_string, dz |
---|
[1320] | 105 | |
---|
| 106 | USE cpulog, & |
---|
| 107 | ONLY: cpu_log, log_point_s |
---|
| 108 | |
---|
| 109 | USE grid_variables, & |
---|
[1822] | 110 | ONLY: dx, dy |
---|
[1320] | 111 | |
---|
| 112 | USE kinds |
---|
| 113 | |
---|
| 114 | USE lpm_collision_kernels_mod, & |
---|
[1822] | 115 | ONLY: ckernel, recalculate_kernel |
---|
[1320] | 116 | |
---|
| 117 | USE particle_attributes, & |
---|
[1822] | 118 | ONLY: all_or_nothing, average_impact, dissipation_classes, & |
---|
| 119 | hall_kernel, iran_part, number_of_particles, particles, & |
---|
| 120 | particle_type, prt_count, use_kernel_tables, wang_kernel |
---|
[1320] | 121 | |
---|
[1822] | 122 | USE random_function_mod, & |
---|
| 123 | ONLY: random_function |
---|
| 124 | |
---|
[1359] | 125 | USE pegrid |
---|
| 126 | |
---|
[849] | 127 | IMPLICIT NONE |
---|
| 128 | |
---|
[1682] | 129 | INTEGER(iwp) :: eclass !< |
---|
| 130 | INTEGER(iwp) :: i !< |
---|
| 131 | INTEGER(iwp) :: j !< |
---|
| 132 | INTEGER(iwp) :: k !< |
---|
| 133 | INTEGER(iwp) :: n !< |
---|
[1822] | 134 | INTEGER(iwp) :: m !< |
---|
[1682] | 135 | INTEGER(iwp) :: rclass_l !< |
---|
| 136 | INTEGER(iwp) :: rclass_s !< |
---|
[849] | 137 | |
---|
[1822] | 138 | REAL(wp) :: collection_probability !< probability for collection |
---|
| 139 | REAL(wp) :: ddV !< inverse grid box volume |
---|
| 140 | REAL(wp) :: epsilon !< dissipation rate |
---|
| 141 | REAL(wp) :: factor_volume_to_mass !< 4.0 / 3.0 * pi * rho_l |
---|
| 142 | REAL(wp) :: xm !< mean mass of droplet m |
---|
| 143 | REAL(wp) :: xn !< mean mass of droplet n |
---|
[1359] | 144 | |
---|
[1822] | 145 | REAL(wp), DIMENSION(:), ALLOCATABLE :: weight !< weighting factor |
---|
| 146 | REAL(wp), DIMENSION(:), ALLOCATABLE :: mass !< total mass of super droplet |
---|
[1359] | 147 | |
---|
[849] | 148 | CALL cpu_log( log_point_s(43), 'lpm_droplet_coll', 'start' ) |
---|
| 149 | |
---|
[1822] | 150 | number_of_particles = prt_count(k,j,i) |
---|
| 151 | factor_volume_to_mass = 4.0_wp / 3.0_wp * pi * rho_l |
---|
| 152 | ddV = 1 / ( dx * dy * dz ) |
---|
[849] | 153 | ! |
---|
[1822] | 154 | !-- Collision requires at least one super droplet inside the box |
---|
| 155 | IF ( number_of_particles > 0 ) THEN |
---|
[849] | 156 | |
---|
| 157 | ! |
---|
[1359] | 158 | !-- Now apply the different kernels |
---|
[1822] | 159 | IF ( use_kernel_tables ) THEN |
---|
[849] | 160 | ! |
---|
[1822] | 161 | !-- Fast method with pre-calculated collection kernels for |
---|
[1359] | 162 | !-- discrete radius- and dissipation-classes. |
---|
| 163 | !-- |
---|
| 164 | !-- Determine dissipation class index of this gridbox |
---|
| 165 | IF ( wang_kernel ) THEN |
---|
[1860] | 166 | eclass = INT( diss(k,j,i) * 1.0E4_wp / 600.0_wp * & |
---|
[1359] | 167 | dissipation_classes ) + 1 |
---|
| 168 | epsilon = diss(k,j,i) |
---|
| 169 | ELSE |
---|
| 170 | epsilon = 0.0_wp |
---|
| 171 | ENDIF |
---|
| 172 | IF ( hall_kernel .OR. epsilon * 1.0E4_wp < 0.001_wp ) THEN |
---|
| 173 | eclass = 0 ! Hall kernel is used |
---|
| 174 | ELSE |
---|
| 175 | eclass = MIN( dissipation_classes, eclass ) |
---|
| 176 | ENDIF |
---|
| 177 | |
---|
[849] | 178 | ! |
---|
[1359] | 179 | !-- Droplet collision are calculated using collision-coalescence |
---|
| 180 | !-- formulation proposed by Wang (see PALM documentation) |
---|
[1822] | 181 | !-- Temporary fields for total mass of super-droplet and weighting factors |
---|
| 182 | !-- are allocated. |
---|
| 183 | ALLOCATE(mass(1:number_of_particles), weight(1:number_of_particles)) |
---|
[849] | 184 | |
---|
[1822] | 185 | mass(1:number_of_particles) = particles(1:number_of_particles)%weight_factor * & |
---|
| 186 | particles(1:number_of_particles)%radius**3 * & |
---|
| 187 | factor_volume_to_mass |
---|
| 188 | weight(1:number_of_particles) = particles(1:number_of_particles)%weight_factor |
---|
[849] | 189 | |
---|
[1822] | 190 | IF ( average_impact ) THEN ! select collision algorithm |
---|
[1071] | 191 | |
---|
[1822] | 192 | DO n = 1, number_of_particles |
---|
[1071] | 193 | |
---|
[1822] | 194 | rclass_l = particles(n)%class |
---|
| 195 | xn = mass(n) / weight(n) |
---|
[849] | 196 | |
---|
[1822] | 197 | DO m = n, number_of_particles |
---|
[849] | 198 | |
---|
[1822] | 199 | rclass_s = particles(m)%class |
---|
| 200 | xm = mass(m) / weight(m) |
---|
[849] | 201 | |
---|
[1822] | 202 | IF ( xm .LT. xn ) THEN |
---|
| 203 | |
---|
| 204 | ! |
---|
| 205 | !-- Particle n collects smaller particle m |
---|
| 206 | collection_probability = ckernel(rclass_l,rclass_s,eclass) * & |
---|
| 207 | weight(n) * ddV * dt_3d |
---|
[849] | 208 | |
---|
[1822] | 209 | mass(n) = mass(n) + mass(m) * collection_probability |
---|
| 210 | weight(m) = weight(m) - weight(m) * collection_probability |
---|
| 211 | mass(m) = mass(m) - mass(m) * collection_probability |
---|
| 212 | ELSEIF ( xm .GT. xn ) THEN |
---|
[849] | 213 | ! |
---|
[1822] | 214 | !-- Particle m collects smaller particle n |
---|
| 215 | collection_probability = ckernel(rclass_l,rclass_s,eclass) * & |
---|
| 216 | weight(m) * ddV * dt_3d |
---|
[849] | 217 | |
---|
[1822] | 218 | mass(m) = mass(m) + mass(n) * collection_probability |
---|
| 219 | weight(n) = weight(n) - weight(n) * collection_probability |
---|
| 220 | mass(n) = mass(n) - mass(n) * collection_probability |
---|
| 221 | ELSE |
---|
[1071] | 222 | ! |
---|
[1822] | 223 | !-- Same-size collections. If n = m, weight is reduced by the |
---|
| 224 | !-- number of possible same-size collections; the total mass |
---|
| 225 | !-- is not changed during same-size collection. |
---|
| 226 | !-- Same-size collections of different |
---|
| 227 | !-- particles ( n /= m ) are treated as same-size collections |
---|
| 228 | !-- of ONE partilce with weight = weight(n) + weight(m) and |
---|
| 229 | !-- mass = mass(n) + mass(m). |
---|
| 230 | !-- Accordingly, each particle loses the same number of |
---|
| 231 | !-- droplets to the other particle, but this has no effect on |
---|
| 232 | !-- total mass mass, since the exchanged droplets have the |
---|
| 233 | !-- same radius. |
---|
[849] | 234 | |
---|
[1822] | 235 | !-- Note: For m = n this equation is an approximation only |
---|
| 236 | !-- valid for weight >> 1 (which is usually the case). The |
---|
| 237 | !-- approximation is weight(n)-1 = weight(n). |
---|
| 238 | weight(n) = weight(n) - 0.5_wp * weight(n) * & |
---|
| 239 | ckernel(rclass_l,rclass_s,eclass) * & |
---|
| 240 | weight(m) * ddV * dt_3d |
---|
| 241 | IF ( n .NE. m ) THEN |
---|
| 242 | weight(m) = weight(m) - 0.5_wp * weight(m) * & |
---|
| 243 | ckernel(rclass_l,rclass_s,eclass) * & |
---|
| 244 | weight(n) * ddV * dt_3d |
---|
| 245 | ENDIF |
---|
| 246 | ENDIF |
---|
[1071] | 247 | |
---|
[1822] | 248 | ENDDO |
---|
[1071] | 249 | |
---|
[1822] | 250 | ql_vp(k,j,i) = ql_vp(k,j,i) + mass(n) / factor_volume_to_mass |
---|
| 251 | |
---|
[1359] | 252 | ENDDO |
---|
[849] | 253 | |
---|
[1822] | 254 | ELSEIF ( all_or_nothing ) THEN ! select collision algorithm |
---|
[849] | 255 | |
---|
[1822] | 256 | DO n = 1, number_of_particles |
---|
[849] | 257 | |
---|
[1822] | 258 | rclass_l = particles(n)%class |
---|
| 259 | xn = mass(n) / weight(n) ! mean mass of droplet n |
---|
[849] | 260 | |
---|
[1822] | 261 | DO m = n, number_of_particles |
---|
[849] | 262 | |
---|
[1822] | 263 | rclass_s = particles(m)%class |
---|
| 264 | xm = mass(m) / weight(m) ! mean mass of droplet m |
---|
[849] | 265 | |
---|
[1822] | 266 | IF ( weight(n) .LT. weight(m) ) THEN |
---|
| 267 | ! |
---|
| 268 | !-- Particle n collects weight(n) droplets of particle m |
---|
| 269 | collection_probability = ckernel(rclass_l,rclass_s,eclass) * & |
---|
| 270 | weight(m) * ddV * dt_3d |
---|
[1071] | 271 | |
---|
[1822] | 272 | IF ( collection_probability .GT. random_function( iran_part ) ) THEN |
---|
| 273 | mass(n) = mass(n) + weight(n) * xm |
---|
| 274 | weight(m) = weight(m) - weight(n) |
---|
| 275 | mass(m) = mass(m) - weight(n) * xm |
---|
| 276 | ENDIF |
---|
| 277 | |
---|
| 278 | ELSEIF ( weight(m) .LT. weight(n) ) THEN |
---|
[849] | 279 | ! |
---|
[1822] | 280 | !-- Particle m collects weight(m) droplets of particle n |
---|
| 281 | collection_probability = ckernel(rclass_l,rclass_s,eclass) * & |
---|
| 282 | weight(n) * ddV * dt_3d |
---|
| 283 | |
---|
| 284 | IF ( collection_probability .GT. random_function( iran_part ) ) THEN |
---|
| 285 | mass(m) = mass(m) + weight(m) * xn |
---|
| 286 | weight(n) = weight(n) - weight(m) |
---|
| 287 | mass(n) = mass(n) - weight(m) * xn |
---|
| 288 | ENDIF |
---|
| 289 | ELSE |
---|
[849] | 290 | ! |
---|
[1822] | 291 | !-- Collisions of particles of the same weighting factor. |
---|
| 292 | !-- Particle n collects 1/2 weight(n) droplets of particle m, |
---|
| 293 | !-- particle m collects 1/2 weight(m) droplets of particle n. |
---|
| 294 | !-- The total mass mass changes accordingly. |
---|
| 295 | !-- If n = m, the first half of the droplets coalesces with the |
---|
| 296 | !-- second half of the droplets; mass is unchanged because |
---|
| 297 | !-- xm = xn for n = m. |
---|
[849] | 298 | |
---|
[1822] | 299 | !-- Note: For m = n this equation is an approximation only |
---|
| 300 | !-- valid for weight >> 1 (which is usually the case). The |
---|
| 301 | !-- approximation is weight(n)-1 = weight(n). |
---|
| 302 | collection_probability = ckernel(rclass_l,rclass_s,eclass) * & |
---|
| 303 | weight(n) * ddV * dt_3d |
---|
[849] | 304 | |
---|
[1822] | 305 | IF ( collection_probability .GT. random_function( iran_part ) ) THEN |
---|
| 306 | mass(n) = mass(n) + 0.5_wp * weight(n) * ( xm - xn ) |
---|
| 307 | mass(m) = mass(m) + 0.5_wp * weight(m) * ( xn - xm ) |
---|
| 308 | weight(n) = weight(n) - 0.5_wp * weight(m) |
---|
| 309 | weight(m) = weight(n) |
---|
| 310 | ENDIF |
---|
| 311 | ENDIF |
---|
[849] | 312 | |
---|
[1822] | 313 | ENDDO |
---|
[849] | 314 | |
---|
[1822] | 315 | ql_vp(k,j,i) = ql_vp(k,j,i) + mass(n) / factor_volume_to_mass |
---|
[849] | 316 | |
---|
[1822] | 317 | ENDDO |
---|
[849] | 318 | |
---|
[1822] | 319 | ENDIF |
---|
[849] | 320 | |
---|
| 321 | |
---|
| 322 | |
---|
| 323 | |
---|
[1822] | 324 | IF ( ANY(weight < 0.0_wp) ) THEN |
---|
| 325 | WRITE( message_string, * ) 'negative weighting' |
---|
| 326 | CALL message( 'lpm_droplet_collision', 'PA0028', & |
---|
| 327 | 2, 2, -1, 6, 1 ) |
---|
| 328 | ENDIF |
---|
[849] | 329 | |
---|
[1822] | 330 | particles(1:number_of_particles)%radius = ( mass(1:number_of_particles) / & |
---|
| 331 | ( weight(1:number_of_particles) & |
---|
| 332 | * factor_volume_to_mass & |
---|
| 333 | ) & |
---|
| 334 | )**0.33333333333333_wp |
---|
[849] | 335 | |
---|
[1822] | 336 | particles(1:number_of_particles)%weight_factor = weight(1:number_of_particles) |
---|
[849] | 337 | |
---|
[1822] | 338 | DEALLOCATE(weight, mass) |
---|
[849] | 339 | |
---|
[1822] | 340 | ELSEIF ( .NOT. use_kernel_tables ) THEN |
---|
| 341 | ! |
---|
| 342 | !-- Collection kernels are calculated for every new |
---|
| 343 | !-- grid box. First, allocate memory for kernel table. |
---|
| 344 | !-- Third dimension is 1, because table is re-calculated for |
---|
| 345 | !-- every new dissipation value. |
---|
| 346 | ALLOCATE( ckernel(1:number_of_particles,1:number_of_particles,1:1) ) |
---|
| 347 | ! |
---|
| 348 | !-- Now calculate collection kernel for this box. Note that |
---|
| 349 | !-- the kernel is based on the previous time step |
---|
| 350 | CALL recalculate_kernel( i, j, k ) |
---|
| 351 | ! |
---|
| 352 | !-- Droplet collision are calculated using collision-coalescence |
---|
| 353 | !-- formulation proposed by Wang (see PALM documentation) |
---|
| 354 | !-- Temporary fields for total mass of super-droplet and weighting factors |
---|
| 355 | !-- are allocated. |
---|
| 356 | ALLOCATE(mass(1:number_of_particles), weight(1:number_of_particles)) |
---|
[849] | 357 | |
---|
[1822] | 358 | mass(1:number_of_particles) = particles(1:number_of_particles)%weight_factor * & |
---|
| 359 | particles(1:number_of_particles)%radius**3 * & |
---|
| 360 | factor_volume_to_mass |
---|
[849] | 361 | |
---|
[1822] | 362 | weight(1:number_of_particles) = particles(1:number_of_particles)%weight_factor |
---|
[849] | 363 | |
---|
[1822] | 364 | IF ( average_impact ) THEN ! select collision algorithm |
---|
[849] | 365 | |
---|
[1822] | 366 | DO n = 1, number_of_particles |
---|
[849] | 367 | |
---|
[1822] | 368 | xn = mass(n) / weight(n) ! mean mass of droplet n |
---|
[849] | 369 | |
---|
[1822] | 370 | DO m = n, number_of_particles |
---|
[849] | 371 | |
---|
[1822] | 372 | xm = mass(m) / weight(m) !mean mass of droplet m |
---|
[849] | 373 | |
---|
[1822] | 374 | IF ( xm .LT. xn ) THEN |
---|
[849] | 375 | ! |
---|
[1822] | 376 | !-- Particle n collects smaller particle m |
---|
| 377 | collection_probability = ckernel(n,m,1) * weight(n) * & |
---|
| 378 | ddV * dt_3d |
---|
| 379 | |
---|
| 380 | mass(n) = mass(n) + mass(m) * collection_probability |
---|
| 381 | weight(m) = weight(m) - weight(m) * collection_probability |
---|
| 382 | mass(m) = mass(m) - mass(m) * collection_probability |
---|
| 383 | ELSEIF ( xm .GT. xn ) THEN |
---|
[849] | 384 | ! |
---|
[1822] | 385 | !-- Particle m collects smaller particle n |
---|
| 386 | collection_probability = ckernel(n,m,1) * weight(m) * & |
---|
| 387 | ddV * dt_3d |
---|
[849] | 388 | |
---|
[1822] | 389 | mass(m) = mass(m) + mass(n) * collection_probability |
---|
| 390 | weight(n) = weight(n) - weight(n) * collection_probability |
---|
| 391 | mass(n) = mass(n) - mass(n) * collection_probability |
---|
| 392 | ELSE |
---|
[849] | 393 | ! |
---|
[1822] | 394 | !-- Same-size collections. If n = m, weight is reduced by the |
---|
| 395 | !-- number of possible same-size collections; the total mass |
---|
| 396 | !-- mass is not changed during same-size collection. |
---|
| 397 | !-- Same-size collections of different |
---|
| 398 | !-- particles ( n /= m ) are treated as same-size collections |
---|
| 399 | !-- of ONE partilce with weight = weight(n) + weight(m) and |
---|
| 400 | !-- mass = mass(n) + mass(m). |
---|
| 401 | !-- Accordingly, each particle loses the same number of |
---|
| 402 | !-- droplets to the other particle, but this has no effect on |
---|
| 403 | !-- total mass mass, since the exchanged droplets have the |
---|
| 404 | !-- same radius. |
---|
| 405 | !-- |
---|
| 406 | !-- Note: For m = n this equation is an approximation only |
---|
| 407 | !-- valid for weight >> 1 (which is usually the case). The |
---|
| 408 | !-- approximation is weight(n)-1 = weight(n). |
---|
| 409 | weight(n) = weight(n) - 0.5_wp * weight(n) * & |
---|
| 410 | ckernel(n,m,1) * weight(m) * & |
---|
| 411 | ddV * dt_3d |
---|
| 412 | IF ( n .NE. m ) THEN |
---|
| 413 | weight(m) = weight(m) - 0.5_wp * weight(m) * & |
---|
| 414 | ckernel(n,m,1) * weight(n) * & |
---|
| 415 | ddV * dt_3d |
---|
| 416 | ENDIF |
---|
| 417 | ENDIF |
---|
[849] | 418 | |
---|
| 419 | |
---|
[1359] | 420 | ENDDO |
---|
[849] | 421 | |
---|
[1822] | 422 | ql_vp(k,j,i) = ql_vp(k,j,i) + mass(n) / factor_volume_to_mass |
---|
[849] | 423 | |
---|
[1822] | 424 | ENDDO |
---|
[849] | 425 | |
---|
[1822] | 426 | ELSEIF ( all_or_nothing ) THEN ! select collision algorithm |
---|
[849] | 427 | |
---|
[1822] | 428 | DO n = 1, number_of_particles |
---|
[849] | 429 | |
---|
[1822] | 430 | xn = mass(n) / weight(n) ! mean mass of droplet n |
---|
[849] | 431 | |
---|
[1822] | 432 | DO m = n, number_of_particles |
---|
[849] | 433 | |
---|
[1822] | 434 | xm = mass(m) / weight(m) !mean mass of droplet m |
---|
[849] | 435 | |
---|
[1822] | 436 | IF ( weight(n) .LT. weight(m) ) THEN |
---|
| 437 | ! |
---|
| 438 | !-- Particle n collects smaller particle m |
---|
| 439 | collection_probability = ckernel(n,m,1) * weight(m) * & |
---|
| 440 | ddV * dt_3d |
---|
[1071] | 441 | |
---|
[1822] | 442 | IF ( collection_probability .GT. random_function( iran_part ) ) THEN |
---|
| 443 | mass(n) = mass(n) + weight(n) * xm |
---|
| 444 | weight(m) = weight(m) - weight(n) |
---|
| 445 | mass(m) = mass(m) - weight(n) * xm |
---|
| 446 | ENDIF |
---|
[1359] | 447 | |
---|
[1822] | 448 | ELSEIF ( weight(m) .LT. weight(n) ) THEN |
---|
[1071] | 449 | ! |
---|
[1822] | 450 | !-- Particle m collects smaller particle n |
---|
| 451 | collection_probability = ckernel(n,m,1) * weight(n) * & |
---|
| 452 | ddV * dt_3d |
---|
[1071] | 453 | |
---|
[1822] | 454 | IF ( collection_probability .GT. random_function( iran_part ) ) THEN |
---|
| 455 | mass(m) = mass(m) + weight(m) * xn |
---|
| 456 | weight(n) = weight(n) - weight(m) |
---|
| 457 | mass(n) = mass(n) - weight(m) * xn |
---|
| 458 | ENDIF |
---|
| 459 | ELSE |
---|
| 460 | ! |
---|
| 461 | !-- Collisions of particles of the same weighting factor. |
---|
| 462 | !-- Particle n collects 1/2 weight(n) droplets of particle m, |
---|
| 463 | !-- particle m collects 1/2 weight(m) droplets of particle n. |
---|
| 464 | !-- The total mass mass changes accordingly. |
---|
| 465 | !-- If n = m, the first half of the droplets coalesces with the |
---|
| 466 | !-- second half of the droplets; mass is unchanged because |
---|
| 467 | !-- xm = xn for n = m. |
---|
| 468 | !-- |
---|
| 469 | !-- Note: For m = n this equation is an approximation only |
---|
| 470 | !-- valid for weight >> 1 (which is usually the case). The |
---|
| 471 | !-- approximation is weight(n)-1 = weight(n). |
---|
| 472 | collection_probability = ckernel(n,m,1) * weight(n) * & |
---|
| 473 | ddV * dt_3d |
---|
[1071] | 474 | |
---|
[1822] | 475 | IF ( collection_probability .GT. random_function( iran_part ) ) THEN |
---|
| 476 | mass(n) = mass(n) + 0.5_wp * weight(n) * ( xm - xn ) |
---|
| 477 | mass(m) = mass(m) + 0.5_wp * weight(m) * ( xn - xm ) |
---|
| 478 | weight(n) = weight(n) - 0.5_wp * weight(m) |
---|
| 479 | weight(m) = weight(n) |
---|
| 480 | ENDIF |
---|
| 481 | ENDIF |
---|
[1071] | 482 | |
---|
| 483 | |
---|
[1822] | 484 | ENDDO |
---|
[1071] | 485 | |
---|
[1822] | 486 | ql_vp(k,j,i) = ql_vp(k,j,i) + mass(n) / factor_volume_to_mass |
---|
[1071] | 487 | |
---|
[1822] | 488 | ENDDO |
---|
[1071] | 489 | |
---|
[1822] | 490 | ENDIF |
---|
[1071] | 491 | |
---|
[1822] | 492 | IF ( ANY(weight < 0.0_wp) ) THEN |
---|
| 493 | WRITE( message_string, * ) 'negative weighting' |
---|
| 494 | CALL message( 'lpm_droplet_collision', 'PA0028', & |
---|
| 495 | 2, 2, -1, 6, 1 ) |
---|
| 496 | ENDIF |
---|
| 497 | |
---|
| 498 | particles(1:number_of_particles)%radius = ( mass(1:number_of_particles) / & |
---|
| 499 | ( weight(1:number_of_particles) & |
---|
| 500 | * factor_volume_to_mass & |
---|
| 501 | ) & |
---|
| 502 | )**0.33333333333333_wp |
---|
| 503 | |
---|
| 504 | particles(1:number_of_particles)%weight_factor = weight(1:number_of_particles) |
---|
| 505 | |
---|
| 506 | DEALLOCATE( weight, mass, ckernel ) |
---|
| 507 | |
---|
| 508 | ENDIF |
---|
| 509 | |
---|
[849] | 510 | ! |
---|
[1884] | 511 | !-- Check if LWC is conserved during collision process |
---|
| 512 | IF ( ql_v(k,j,i) /= 0.0_wp ) THEN |
---|
| 513 | IF ( ql_vp(k,j,i) / ql_v(k,j,i) >= 1.0001_wp .OR. & |
---|
| 514 | ql_vp(k,j,i) / ql_v(k,j,i) <= 0.9999_wp ) THEN |
---|
| 515 | WRITE( message_string, * ) ' LWC is not conserved during', & |
---|
| 516 | ' collision! ', & |
---|
| 517 | ' LWC after condensation: ', ql_v(k,j,i), & |
---|
| 518 | ' LWC after collision: ', ql_vp(k,j,i) |
---|
| 519 | CALL message( 'lpm_droplet_collision', 'PA0040', 2, 2, -1, 6, 1 ) |
---|
| 520 | ENDIF |
---|
[1359] | 521 | ENDIF |
---|
[1884] | 522 | |
---|
[1359] | 523 | ENDIF |
---|
[1884] | 524 | |
---|
[849] | 525 | CALL cpu_log( log_point_s(43), 'lpm_droplet_coll', 'stop' ) |
---|
| 526 | |
---|
[1359] | 527 | END SUBROUTINE lpm_droplet_collision |
---|