Changeset 1065 for palm/trunk
- Timestamp:
- Nov 22, 2012 5:42:36 PM (12 years ago)
- Location:
- palm/trunk/SOURCE
- Files:
-
- 11 edited
-
check_parameters.f90 (modified) (2 diffs)
-
data_output_2d.f90 (modified) (3 diffs)
-
diffusion_e.f90 (modified) (7 diffs)
-
init_3d_model.f90 (modified) (2 diffs)
-
init_cloud_physics.f90 (modified) (2 diffs)
-
microphysics.f90 (modified) (11 diffs, 1 prop)
-
modules.f90 (modified) (6 diffs)
-
parin.f90 (modified) (4 diffs)
-
read_var_list.f90 (modified) (6 diffs)
-
time_integration.f90 (modified) (2 diffs)
-
write_var_list.f90 (modified) (6 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/check_parameters.f90
r1061 r1065 20 20 ! Current revisions: 21 21 ! ----------------- 22 ! Bugfix: It is not allowed to use cloud_scheme = seifert_beheng without 23 ! precipitation in order to save computational resources. 22 24 ! 23 25 ! … … 886 888 'loop_optimization = cache' 887 889 CALL message( 'check_parameters', 'PA0362', 1, 2, 0, 6, 0 ) 890 ENDIF 891 892 IF ( cloud_physics .AND. icloud_scheme == 0 .AND. & 893 .NOT. precipitation ) THEN 894 message_string = 'cloud_scheme = seifert_beheng requires ' // & 895 'precipitation = .TRUE.' 896 CALL message( 'check_parameters', 'PA0363', 1, 2, 0, 6, 0 ) 888 897 ENDIF 889 898 -
palm/trunk/SOURCE/data_output_2d.f90
r1054 r1065 20 20 ! Current revisions: 21 21 ! ----------------- 22 ! Bugfix: Output of cross sections of ql 22 23 ! 23 24 ! Former revisions: … … 488 489 DO i = nxlg, nxrg 489 490 DO j = nysg, nyng 490 DO k = nzb, nz _do3d491 DO k = nzb, nzt+1 491 492 local_pf(i,j,k) = ql(k,j,i) + qr(k,j,i) 492 493 ENDDO … … 501 502 DO i = nxlg, nxrg 502 503 DO j = nysg, nyng 503 DO k = nzb, nz _do3d504 DO k = nzb, nzt+1 504 505 local_pf(i,j,k) = ql_av(k,j,i) + qr_av(k,j,i) 505 506 ENDDO -
palm/trunk/SOURCE/diffusion_e.f90
r1037 r1065 20 20 ! Current revisions: 21 21 ! ----------------- 22 ! 22 ! Enabled the claculation of diss in case of turbulence = .TRUE. (parameterized 23 ! effects of turbulence on autoconversion and accretion in two-moments cloud 24 ! physics scheme). 23 25 ! 24 26 ! Former revisions: … … 183 185 !-- Store dissipation if needed for calculating the sgs particle 184 186 !-- velocities 185 IF ( use_sgs_for_particles .OR. wang_kernel ) THEN 187 IF ( use_sgs_for_particles .OR. wang_kernel .OR. & 188 turbulence ) THEN 186 189 DO j = nys, nyn 187 190 DO k = nzb_s_inner(j,i)+1, nzt … … 254 257 !-- Store dissipation if needed for calculating the sgs particle 255 258 !-- velocities 256 IF ( use_sgs_for_particles .OR. wang_kernel ) THEN 259 IF ( use_sgs_for_particles .OR. wang_kernel .OR. & 260 turbulence ) THEN 257 261 DO j = nys, nyn 258 262 DO k = nzb_s_inner(j,i)+1, nzt … … 268 272 ! 269 273 !-- Boundary condition for dissipation 270 IF ( use_sgs_for_particles .OR. wang_kernel ) THEN274 IF ( use_sgs_for_particles .OR. wang_kernel .OR. turbulence ) THEN 271 275 DO i = nxl, nxr 272 276 DO j = nys, nyn … … 433 437 !-- Store dissipation if needed for calculating the sgs 434 438 !-- particle velocities 435 IF ( use_sgs_for_particles .OR. wang_kernel ) THEN 439 IF ( use_sgs_for_particles .OR. wang_kernel .OR. & 440 turbulence ) THEN 436 441 diss(k,j,i) = dissipation 437 442 ENDIF … … 448 453 ! 449 454 !-- Boundary condition for dissipation 450 IF ( use_sgs_for_particles .OR. wang_kernel ) THEN455 IF ( use_sgs_for_particles .OR. wang_kernel .OR. turbulence ) THEN 451 456 !$acc kernels present( diss, nzb_s_inner ) 452 457 !$acc loop … … 539 544 ! 540 545 !-- Store dissipation if needed for calculating the sgs particle velocities 541 IF ( use_sgs_for_particles .OR. wang_kernel ) THEN546 IF ( use_sgs_for_particles .OR. wang_kernel .OR. turbulence ) THEN 542 547 DO k = nzb_s_inner(j,i)+1, nzt 543 548 diss(k,j,i) = dissipation(k) -
palm/trunk/SOURCE/init_3d_model.f90
r1054 r1065 23 23 ! Current revisions: 24 24 ! ------------------ 25 ! allocation of diss (dissipation rate) in case of turbulence = .TRUE. added 25 26 ! 26 27 ! Former revisions: … … 467 468 !-- 3D-array for storing the dissipation, needed for calculating the sgs 468 469 !-- particle velocities 469 IF ( use_sgs_for_particles .OR. wang_kernel ) THEN470 IF ( use_sgs_for_particles .OR. wang_kernel .OR. turbulence ) THEN 470 471 ALLOCATE ( diss(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 471 472 ELSE -
palm/trunk/SOURCE/init_cloud_physics.f90
r1054 r1065 20 20 ! Current revisions: 21 21 ! ----------------- 22 ! The Courant number of sedimentation can be controlled with c_sedimentation. 22 23 ! 23 24 ! Former revisions: … … 88 89 !-- Calculate timestep according to precipitation 89 90 IF ( icloud_scheme == 0 .AND. precipitation ) THEN 90 dt_precipitation = MINVAL( dzu(nzb+2:nzt) ) / w_precipitation 91 dt_precipitation = c_sedimentation * MINVAL( dzu(nzb+2:nzt) ) / & 92 w_precipitation 91 93 ENDIF 92 94 ! -
palm/trunk/SOURCE/microphysics.f90
- Property svn:keywords set to Id
r1054 r1065 4 4 ! Current revisions: 5 5 ! ----------------- 6 ! Sedimentation process implemented according to Stevens and Seifert (2008). 7 ! Turbulence effects on autoconversion and accretion added (Seifert, Nuijens 8 ! and Stevens, 2010). 6 9 ! 7 10 ! Former revisions: … … 241 244 242 245 INTEGER :: i, j, k 243 REAL :: dr_min = 2.0E-6, dr_max = 1.0E-3244 246 245 247 DO k = nzb_2d(j,i)+1, nzt 246 IF ( ( qr(k,j,i) > eps_sb ) .AND. ( nr(k,j,i) > eps_sb ) ) THEN 248 249 IF ( qr(k,j,i) <= eps_sb ) THEN 250 qr(k,j,i) = 0.0 251 ELSE 252 ! 253 !-- Adjust number of raindrops to avoid nonlinear effects in 254 !-- sedimentation and evaporation of rain drops due to too small or 255 !-- too big weights of rain drops (Stevens and Seifert, 2008). 256 IF ( nr(k,j,i) * xrmin > qr(k,j,i) * hyrho(k) ) THEN 257 nr(k,j,i) = qr(k,j,i) * hyrho(k) / xrmin 258 ELSEIF ( nr(k,j,i) * xrmax < qr(k,j,i) * hyrho(k) ) THEN 259 nr(k,j,i) = qr(k,j,i) * hyrho(k) / xrmax 260 ENDIF 261 xr(k) = hyrho(k) * qr(k,j,i) / nr(k,j,i) 247 262 ! 248 263 !-- Weight averaged diameter of rain drops: … … 250 265 dpirho_l )**( 1.0 / 3.0 ) 251 266 ! 252 !-- Adjust number of raindrops to avoid nonlinear effects in253 !-- sedimentation and evaporation of rain drops due to too small or254 !-- big diameters of rain drops (Stevens and Seifert, 2008).255 IF ( dr(k) < dr_min ) THEN256 nr(k,j,i) = qr(k,j,i) * hyrho(k) / dr_min**3 * dpirho_l257 dr(k) = dr_min258 ELSEIF ( dr(k) > dr_max ) THEN259 nr(k,j,i) = qr(k,j,i) * hyrho(k) / dr_max**3 * dpirho_l260 dr(k) = dr_max261 ENDIF262 !263 !-- Mean weight of rain drops (Seifert and Beheng, 2006):264 xr(k) = MIN( MAX( hyrho(k) * qr(k,j,i) / nr(k,j,i), xrmin ), xrmax)265 !266 267 !-- Shape parameter of gamma distribution (Milbrandt and Yau, 2005; 267 268 !-- Stevens and Seifert, 2008): 268 IF ( .NOT. mu_constant ) THEN 269 mu_r(k) = 10.0 * ( 1.0 + TANH( 1.2E3 * ( dr(k) - 1.4E-3 ) ) ) 270 ELSE 271 mu_r(k) = mu_constant_value 272 ENDIF 269 mu_r(k) = 10.0 * ( 1.0 + TANH( 1.2E3 * ( dr(k) - 1.4E-3 ) ) ) 273 270 ! 274 271 !-- Slope parameter of gamma distribution (Seifert, 2008): … … 288 285 USE control_parameters 289 286 USE statistics 287 USE grid_variables 290 288 291 289 IMPLICIT NONE 292 290 293 291 INTEGER :: i, j, k 294 REAL :: k_au, autocon, phi_au, tau_cloud, xc, nu_c 292 REAL :: k_au, autocon, phi_au, tau_cloud, xc, nu_c, rc, & 293 l_mix, re_lambda, alpha_cc, r_cc, sigma_cc, epsilon 295 294 296 295 k_au = k_cc / ( 20.0 * x0 ) … … 310 309 !-- Shape parameter of gamma distribution (Geoffroy et al., 2010): 311 310 !-- (Use constant nu_c = 1.0 instead?) 312 nu_c = 1 580.0 * hyrho(k) * ql(k,j,i) - 0.28311 nu_c = 1.0 !MAX( 0.0, 1580.0 * hyrho(k) * ql(k,j,i) - 0.28 ) 313 312 ! 314 313 !-- Mean weight of cloud droplets: 315 314 xc = hyrho(k) * ql(k,j,i) / nc 316 315 ! 316 !-- Parameterized turbulence effects on autoconversion (Seifert, 317 !-- Nuijens and Stevens, 2010) 318 IF ( turbulence ) THEN 319 ! 320 !-- Weight averaged radius of cloud droplets: 321 rc = 0.5 * ( xc * dpirho_l )**( 1.0 / 3.0 ) 322 323 alpha_cc = ( a_1 + a_2 * nu_c ) / ( 1.0 + a_3 * nu_c ) 324 r_cc = ( b_1 + b_2 * nu_c ) / ( 1.0 + b_3 * nu_c ) 325 sigma_cc = ( c_1 + c_2 * nu_c ) / ( 1.0 + c_3 * nu_c ) 326 ! 327 !-- Mixing length (neglecting distance to ground and stratification) 328 l_mix = ( dx * dy * dzu(k) )**( 1.0 / 3.0 ) 329 ! 330 !-- Limit dissipation rate according to Seifert, Nuijens and 331 !-- Stevens (2010) 332 epsilon = MIN( 0.06, diss(k,j,i) ) 333 ! 334 !-- Compute Taylor-microscale Reynolds number: 335 re_lambda = 6.0 / 11.0 * ( l_mix / c_const )**( 2.0 / 3.0 ) * & 336 SQRT( 15.0 / kin_vis_air ) * epsilon**( 1.0 / 6.0 ) 337 ! 338 !-- The factor of 1.0E4 is needed to convert the dissipation rate 339 !-- from m2 s-3 to cm2 s-3. 340 k_au = k_au * ( 1.0 + & 341 epsilon * 1.0E4 * ( re_lambda * 1.0E-3 )**0.25 * & 342 ( alpha_cc * EXP( -1.0 * ( ( rc - r_cc ) / & 343 sigma_cc )**2 ) + beta_cc ) ) 344 ENDIF 345 ! 317 346 !-- Autoconversion rate (Seifert and Beheng, 2006): 318 autocon = k_au * ( nu_c + 2.0 ) * ( nu_c + 4.0 ) / &319 ( nu_c + 1.0 )**2 * ql(k,j,i)**2 * xc**2 * &320 ( 1.0 + phi_au / ( 1.0 - tau_cloud )**2 ) * &347 autocon = k_au * ( nu_c + 2.0 ) * ( nu_c + 4.0 ) / & 348 ( nu_c + 1.0 )**2 * ql(k,j,i)**2 * xc**2 * & 349 ( 1.0 + phi_au / ( 1.0 - tau_cloud )**2 ) * & 321 350 rho_surface 322 autocon = MIN( autocon, ql(k,j,i) / ( dt_3d * &351 autocon = MIN( autocon, ql(k,j,i) / ( dt_3d * & 323 352 weight_substep(intermediate_timestep_count) ) ) 324 353 ! … … 346 375 347 376 INTEGER :: i, j, k 348 REAL :: accr, phi_ac, tau_cloud 377 REAL :: accr, phi_ac, tau_cloud, k_cr 349 378 350 379 DO k = nzb_2d(j,i)+1, nzt … … 356 385 !-- Universal function for accretion process 357 386 !-- (Seifert and Beheng, 2001): 358 phi_ac = tau_cloud / ( tau_cloud + 5.0E-5 ) 359 phi_ac = ( phi_ac**2 )**2 387 phi_ac = tau_cloud / ( tau_cloud + 5.0E-5 ) 388 phi_ac = ( phi_ac**2 )**2 389 ! 390 !-- Parameterized turbulence effects on autoconversion (Seifert, 391 !-- Nuijens and Stevens, 2010). The factor of 1.0E4 is needed to 392 !-- convert the dissipation (diss) from m2 s-3 to cm2 s-3. 393 IF ( turbulence ) THEN 394 k_cr = k_cr0 * ( 1.0 + 0.05 * & 395 MIN( 600.0, diss(k,j,i) * 1.0E4 )**0.25 ) 396 ELSE 397 k_cr = k_cr0 398 ENDIF 360 399 ! 361 400 !-- Accretion rate (Seifert and Beheng, 2006): 362 accr 363 364 accr 365 401 accr = k_cr * ql(k,j,i) * qr(k,j,i) * phi_ac * & 402 SQRT( rho_surface * hyrho(k) ) 403 accr = MIN( accr, ql(k,j,i) / ( dt_3d * & 404 weight_substep(intermediate_timestep_count) ) ) 366 405 ! 367 406 !-- Tendencies for q, qr, pt: … … 390 429 391 430 DO k = nzb_2d(j,i)+1, nzt 392 IF ( ( qr(k,j,i) > eps_sb ) .AND. ( nr(k,j,i) > eps_sb )) THEN431 IF ( qr(k,j,i) > eps_sb ) THEN 393 432 ! 394 433 !-- Selfcollection rate (Seifert and Beheng, 2006): 395 434 !-- pirho_l**( 1.0 / 3.0 ) is necessary to convert [lambda_r] = m-1 to 396 435 !-- kg**( 1.0 / 3.0 ). 397 phi_sc = ( 1.0 + kappa_rr / lambda_r(k) * &398 pirho_l**( 1.0 / 3.0 ) )**( -9 )436 phi_sc = 1.0 !( 1.0 + kappa_rr / lambda_r(k) * & 437 !pirho_l**( 1.0 / 3.0 ) )**( -9 ) 399 438 400 439 selfcoll = k_rr * nr(k,j,i) * qr(k,j,i) * phi_sc * & … … 438 477 439 478 DO k = nzb_2d(j,i)+1, nzt 440 IF ( ( qr(k,j,i) > eps_sb ) .AND. ( nr(k,j,i) > eps_sb )) THEN479 IF ( qr(k,j,i) > eps_sb ) THEN 441 480 ! 442 481 !-- Actual liquid water temperature: … … 555 594 IMPLICIT NONE 556 595 557 INTEGER :: i, j, k, n, n_substep 558 REAL :: sed_nr_tend, sed_qr_tend 596 INTEGER :: i, j, k, k_run, n, n_substep 597 REAL :: c_run, d_max, d_mean, d_min, dt_sedi, flux, mean, z_run 598 599 REAL, DIMENSION(nzb:nzt) :: c_nr, c_qr, d_nr, d_qr, nr_slope, qr_slope, & 600 w_nr, w_qr 601 ! 602 !-- Computation of sedimentation flux. Implementation according to Stevens 603 !-- and Seifert (2008). 559 604 560 605 IF ( intermediate_timestep_count == 1 ) prr(:,j,i) = 0.0 561 606 562 sed_nr = 0.0 563 sed_qr = 0.0 564 565 DO k = nzb_2d(j,i)+1, nzt 566 IF ( ( qr(k,j,i) > eps_sb ) .AND. ( nr(k,j,i) > eps_sb ) ) THEN 567 ! 568 !-- Sedimentation of rain water content and rain drop concentration 569 !-- according to Stevens and Seifert (2008): 570 sed_nr(k) = MIN( 20.0, MAX( 0.1, a_term - b_term * ( 1.0 + & 571 c_term / lambda_r(k) )**( -1.0 * & 572 ( mu_r(k) + 1.0 ) ) ) ) * nr(k,j,i) 573 sed_qr(k) = MIN( 20.0, MAX( 0.1, a_term - b_term * ( 1.0 + & 574 c_term / lambda_r(k) )**( -1.0 * & 575 ( mu_r(k) + 4.0 ) ) ) ) * qr(k,j,i) * hyrho(k) 576 ! 577 !-- Computation of rain rate 578 prr(k,j,i) = prr(k,j,i) + sed_qr(k) / hyrho(k) * & 579 weight_substep(intermediate_timestep_count) 580 ENDIF 581 ENDDO 582 583 DO k = nzb_2d(j,i)+1, nzt 584 sed_nr_tend = MAX( ( sed_nr(k+1) - sed_nr(k) ) * ddzu(k+1), & 585 -nr(k,j,i) / ( dt_3d * & 586 weight_substep(intermediate_timestep_count) ) ) 587 sed_qr_tend = MAX( ( sed_qr(k+1) - sed_qr(k) ) * ddzu(k+1) / & 588 hyrho(k), & 589 -qr(k,j,i) / ( dt_3d * & 590 weight_substep(intermediate_timestep_count) ) ) 591 592 tend_nr(k,j,i) = tend_nr(k,j,i) + sed_nr_tend 593 tend_qr(k,j,i) = tend_qr(k,j,i) + sed_qr_tend 594 ENDDO 607 dt_sedi = dt_3d * weight_substep(intermediate_timestep_count) 608 609 w_nr = 0.0 610 w_qr = 0.0 611 ! 612 !-- Compute velocities 613 DO k = nzb_s_inner(j,i)+1, nzt 614 IF ( qr(k,j,i) > eps_sb ) THEN 615 w_nr(k) = MAX( 0.1, MIN( 20.0, a_term - b_term * ( 1.0 + & 616 c_term / lambda_r(k) )**( -1.0 * ( mu_r(k) + 1.0 ) ) ) ) 617 w_qr(k) = MAX( 0.1, MIN( 20.0, a_term - b_term * ( 1.0 + & 618 c_term / lambda_r(k) )**( -1.0 * ( mu_r(k) + 4.0 ) ) ) ) 619 ELSE 620 w_nr(k) = 0.0 621 w_qr(k) = 0.0 622 ENDIF 623 ENDDO 624 ! 625 !-- Adjust boundary values 626 w_nr(nzb_2d(j,i)) = w_nr(nzb_2d(j,i)+1) 627 w_qr(nzb_2d(j,i)) = w_qr(nzb_2d(j,i)+1) 628 w_nr(nzt) = w_nr(nzt-1) 629 w_qr(nzt) = w_qr(nzt-1) 630 ! 631 !-- Compute Courant number 632 DO k = nzb_s_inner(j,i)+1, nzt-1 633 c_nr(k) = 0.25 * ( w_nr(k-1) + 2.0 * w_nr(k) + w_nr(k+1) ) * & 634 dt_sedi * ddzu(k) 635 c_qr(k) = 0.25 * ( w_qr(k-1) + 2.0 * w_qr(k) + w_qr(k+1) ) * & 636 dt_sedi * ddzu(k) 637 ENDDO 638 ! 639 !-- Limit slopes with monotonized centered (MC) limiter (van Leer, 1977): 640 IF ( limiter_sedimentation ) THEN 641 642 qr(nzb_s_inner(j,i),j,i) = 0.0 643 nr(nzb_s_inner(j,i),j,i) = 0.0 644 qr(nzt,j,i) = 0.0 645 nr(nzt,j,i) = 0.0 646 647 DO k = nzb_s_inner(j,i)+1, nzt-1 648 d_mean = 0.5 * ( qr(k+1,j,i) + qr(k-1,j,i) ) 649 d_min = qr(k,j,i) - MIN( qr(k+1,j,i), qr(k,j,i), qr(k-1,j,i) ) 650 d_max = MAX( qr(k+1,j,i), qr(k,j,i), qr(k-1,j,i) ) - qr(k,j,i) 651 652 qr_slope(k) = SIGN(1.0, d_mean) * MIN ( 2.0 * d_min, 2.0 * d_max, & 653 ABS( d_mean ) ) 654 655 d_mean = 0.5 * ( nr(k+1,j,i) + nr(k-1,j,i) ) 656 d_min = nr(k,j,i) - MIN( nr(k+1,j,i), nr(k,j,i), nr(k-1,j,i) ) 657 d_max = MAX( nr(k+1,j,i), nr(k,j,i), nr(k-1,j,i) ) - nr(k,j,i) 658 659 nr_slope(k) = SIGN(1.0, d_mean) * MIN ( 2.0 * d_min, 2.0 * d_max, & 660 ABS( d_mean ) ) 661 ENDDO 662 663 ELSE 664 nr_slope = 0.0 665 qr_slope = 0.0 666 ENDIF 667 ! 668 !-- Compute sedimentation flux 669 DO k = nzt-2, nzb_s_inner(j,i)+1, -1 670 ! 671 !-- Sum up all rain drop number densities which contribute to the flux 672 !-- through k-1/2 673 flux = 0.0 674 z_run = 0.0 ! height above z(k) 675 k_run = k 676 c_run = MIN( 1.0, c_nr(k) ) 677 DO WHILE ( c_run > 0.0 .AND. k_run <= nzt-1 ) 678 flux = flux + hyrho(k_run) * & 679 ( nr(k_run,j,i) + nr_slope(k_run) * ( 1.0 - c_run ) * & 680 0.5 ) * c_run * dzu(k_run) 681 z_run = z_run + dzu(k_run) 682 k_run = k_run + 1 683 c_run = MIN( 1.0, c_nr(k_run) - z_run * ddzu(k_run) ) 684 ENDDO 685 ! 686 !-- It is not allowed to sediment more rain drop number density than 687 !-- available 688 flux = MIN( flux, & 689 hyrho(k) * dzu(k) * nr(k,j,i) + sed_nr(k+1) * dt_sedi ) 690 691 sed_nr(k) = flux / dt_sedi 692 tend_nr(k,j,i) = tend_nr(k,j,i) + ( sed_nr(k+1) - sed_nr(k) ) * & 693 ddzu(k+1) / hyrho(k) 694 ! 695 !-- Sum up all rain water content which contributes to the flux 696 !-- through k-1/2 697 flux = 0.0 698 z_run = 0.0 ! height above z(k) 699 k_run = k 700 c_run = MIN( 1.0, c_qr(k) ) 701 DO WHILE ( c_run > 0.0 .AND. k_run <= nzt-1 ) 702 flux = flux + hyrho(k_run) * & 703 ( qr(k_run,j,i) + qr_slope(k_run) * ( 1.0 - c_run ) * & 704 0.5 ) * c_run * dzu(k_run) 705 z_run = z_run + dzu(k_run) 706 k_run = k_run + 1 707 c_run = MIN( 1.0, c_qr(k_run) - z_run * ddzu(k_run) ) 708 ENDDO 709 ! 710 !-- It is not allowed to sediment more rain water content than available 711 flux = MIN( flux, & 712 hyrho(k) * dzu(k) * qr(k,j,i) + sed_qr(k+1) * dt_sedi ) 713 714 sed_qr(k) = flux / dt_sedi 715 tend_qr(k,j,i) = tend_qr(k,j,i) + ( sed_qr(k+1) - sed_qr(k) ) * & 716 ddzu(k+1) / hyrho(k) 717 ! 718 !-- Compute the rain rate 719 prr(k,j,i) = prr(k,j,i) + sed_qr(k) / hyrho(k) * & 720 weight_substep(intermediate_timestep_count) 721 ENDDO 595 722 ! 596 723 !-- Precipitation amount … … 616 743 IMPLICIT NONE 617 744 618 REAL :: gamm, xx, & 619 ser, tmp, x_gamm, y_gamm 745 REAL :: gamm, ser, tmp, x_gamm, xx, y_gamm 620 746 INTEGER :: j 621 747 -
palm/trunk/SOURCE/modules.f90
r1054 r1065 20 20 ! Current revisions: 21 21 ! ----------------- 22 ! + c_sedimentation, limiter_sedimentation, turbulence, a_1, a_2, a_3, b_1, b_2, 23 ! + b_3, c_1, c_2, c_3, beta_cc 24 ! 25 ! bottom boundary condition of qr, nr changed from Dirichlet to Neumann 22 26 ! 23 27 ! Former revisions: … … 468 472 469 473 LOGICAL :: curvature_solution_effects = .FALSE., & 470 ventilation_effect = .FALSE., & 471 mu_constant = .FALSE. 472 473 REAL :: a_vent = 0.78, & ! coef. for ventilation effect 474 a_term = 9.65, & ! coef. for terminal velocity (m s-1) 475 b_vent = 0.308, & ! coef. for ventilation effect 476 b_term = 9.8, & ! coef. for terminal velocity (m s-1) 477 bfactor, & 478 c_evap = 0.7, & ! constant in evaporation 479 c_term = 600.0, & ! coef. for terminal velocity (m-1) 474 limiter_sedimentation = .TRUE., & 475 ventilation_effect = .FALSE. 476 477 478 REAL :: a_1 = 8.69E-4, & ! coef. in turb. parametrization (cm-2 s3) 479 a_2 = -7.38E-5, & ! coef. in turb. parametrization (cm-2 s3) 480 a_3 = -1.40E-2, & ! coef. in turb. parametrization 481 a_term = 9.65, & ! coef. for terminal velocity (m s-1) 482 a_vent = 0.78, & ! coef. for ventilation effect 483 b_1 = 11.45E-6, & ! coef. in turb. parametrization (m) 484 b_2 = 9.68E-6, & ! coef. in turb. parametrization (m) 485 b_3 = 0.62, & ! coef. in turb. parametrization 486 b_term = 9.8, & ! coef. for terminal velocity (m s-1) 487 b_vent = 0.308, & ! coef. for ventilation effect 488 beta_cc = 3.09E-4, & ! coef. in turb. parametrization (cm-2 s3) 489 bfactor, & 490 c_1 = 4.82E-6, & ! coef. in turb. parametrization (m) 491 c_2 = 4.8E-6, & ! coef. in turb. parametrization (m) 492 c_3 = 0.76, & ! coef. in turb. parametrization 493 c_const = 0.93, & ! const. in Taylor-microscale Reynolds number 494 c_evap = 0.7, & ! constant in evaporation 495 c_sedimentation = 2.0, & ! Courant number of sedimentation process 496 c_term = 600.0, & ! coef. for terminal velocity (m-1) 480 497 cof(6) = (/ 76.18009172947146, & ! coefficients in the 481 498 -86.50532032941677, & ! numerical … … 484 501 0.1208650973866179E-2, & 485 502 -0.5395239384953E-5 /), & 486 cp = 1005.0, 503 cp = 1005.0, & ! heat capacity of dry air (J kg-1 K-1) 487 504 diff_coeff_l = 0.23E-4, & ! diffusivity of water vapor (m2 s-1) 488 505 effective_coll_efficiency, & 489 eps_ros = 1.0E-4, 490 eps_sb = 1.0E-20, 491 k_cc = 4.44E09,& ! const. rain-rain kernel (m3 kg-2 s-1)492 k_cr = 5.25,& ! const. cloud-rain kernel (m3 kg-1 s-1)493 k_rr = 7.12, 494 k_br = 1000., 495 kappa_rr = 60.7, 506 eps_ros = 1.0E-4, & ! accuracy of Rosenbrock method 507 eps_sb = 1.0E-20, & ! threshold in two-moments scheme 508 k_cc = 9.44E09, & ! const. rain-rain kernel (m3 kg-2 s-1) 509 k_cr0 = 4.33, & ! const. cloud-rain kernel (m3 kg-1 s-1) 510 k_rr = 7.12, & ! const. rain-rain kernel (m3 kg-1 s-1) 511 k_br = 1000., & ! const. in breakup parametrization (m-1) 512 kappa_rr = 60.7, & ! const. in collision kernel (kg-1/3) 496 513 kin_vis_air = 1.4086E-5, & ! kin. viscosity of air (m2 s-1) 497 l_v = 2.5E+06, 514 l_v = 2.5E+06, & ! latent heat of vaporization (J kg-1) 498 515 l_d_cp, l_d_r, l_d_rv, & ! l_v / cp, l_v / r_d, l_v / r_v 499 516 mass_of_solute = 1.0E-17, & ! soluted NaCl (kg) 500 517 molecular_weight_of_solute = 0.05844, & ! mol. m. NaCl (kg mol-1) 501 518 molecular_weight_of_water = 0.01801528, & ! mol. m. H2O (kg mol-1) 502 mu_constant_value = 0.0, & ! shape param. of gamma distribution 503 nc = 70.0E6, & ! cloud droplet concentration 519 nc = 70.0E6, & ! cloud droplet concentration 504 520 prec_time_const = 0.001, & !coef. in Kessler scheme 505 pirho_l, dpirho_l, 506 rho_l = 1.0E3, 507 ql_crit = 0.0005, 508 r_d = 287.0, 509 r_v = 461.51, 510 schmidt = 0.71, 511 schmidt_p_1d3, 521 pirho_l, dpirho_l, & ! pi * rho_l / 6.0; 6.0 / ( pi * rho_l ) 522 rho_l = 1.0E3, & ! density of water (kg m-3) 523 ql_crit = 0.0005, & ! coef. in Kessler scheme 524 r_d = 287.0, & ! sp. gas const. dry air (J kg-1 K-1) 525 r_v = 461.51, & ! sp. gas const. water vapor (J kg-1 K-1) 526 schmidt = 0.71, & ! Schmidt number 527 schmidt_p_1d3, & ! schmidt**( 1.0 / 3.0 ) 512 528 stp = 2.5066282746310005, & ! parameter in gamma function 513 529 thermal_conductivity_l = 2.43E-2, & ! therm. cond. air (J m-1 s-1 K-1) 514 vanthoff = 2.0, 515 x0 = 2.6E-10, 516 xrmin = 2.6E-10, 517 xrmax = 5.0E-6, 530 vanthoff = 2.0, & ! van't Hoff factor for NaCl 531 x0 = 2.6E-10, & ! separating drop mass (kg) 532 xrmin = 2.6E-10, & ! minimum rain drop size (kg) 533 xrmax = 5.0E-6, & ! maximum rain drop site (kg) 518 534 dt_precipitation = 100.0, & ! timestep precipitation (s) 519 535 w_precipitation = 9.65 ! maximum terminal velocity (m s-1) … … 595 611 scalar_advec = 'ws-scheme' 596 612 CHARACTER (LEN=20) :: bc_e_b = 'neumann', bc_lr = 'cyclic', & 597 bc_nr_b = ' dirichlet', bc_nr_t = 'neumann', &613 bc_nr_b = 'neumann', bc_nr_t = 'neumann', & 598 614 bc_ns = 'cyclic', bc_p_b = 'neumann', & 599 615 bc_p_t = 'dirichlet', bc_pt_b = 'dirichlet', & 600 616 bc_pt_t = 'initial_gradient', & 601 617 bc_q_b = 'dirichlet', bc_q_t = 'neumann', & 602 bc_qr_b = ' dirichlet', bc_qr_t = 'neumann',&618 bc_qr_b = 'neumann', bc_qr_t = 'neumann',& 603 619 bc_s_b = 'dirichlet', bc_s_t = 'neumann', & 604 620 bc_sa_t = 'neumann', & … … 710 726 do3d_compress = .FALSE., do_sum = .FALSE., & 711 727 dp_external = .FALSE., dp_smooth = .FALSE., & 712 drizzle = . TRUE., dt_fixed = .FALSE., &728 drizzle = .FALSE., dt_fixed = .FALSE., & 713 729 dt_3d_reached, dt_3d_reached_l, exchange_mg = .FALSE., & 714 730 first_call_lpm = .TRUE., & … … 728 744 scalar_rayleigh_damping = .TRUE., sloping_surface = .FALSE., & 729 745 stop_dt = .FALSE., synchronous_exchange = .FALSE., & 730 terminate_run = .FALSE., turbulent_inflow = .FALSE., & 746 terminate_run = .FALSE., turbulence = .FALSE., & 747 turbulent_inflow = .FALSE., & 731 748 use_prescribed_profile_data = .FALSE., use_reference = .FALSE.,& 732 749 use_surface_fluxes = .FALSE., use_top_fluxes = .FALSE., & -
palm/trunk/SOURCE/parin.f90
r1054 r1065 20 20 ! Current revisions: 21 21 ! ----------------- 22 ! +nc, c_sedimentation, limiter_sedimentation, turbulence 23 ! -mu_constant, mu_constant_value 22 24 ! 23 25 ! Former revisions: … … 199 201 call_psolver_at_all_substeps, canopy_mode, canyon_height, & 200 202 canyon_width_x, canyon_width_y, canyon_wall_left, & 201 canyon_wall_south, cfl_factor, cloud_droplets, cloud_physics, & 202 cloud_scheme, collective_wait, conserve_volume_flow, & 203 canyon_wall_south, c_sedimentation, cfl_factor, cloud_droplets, & 204 cloud_physics, cloud_scheme, collective_wait, & 205 conserve_volume_flow, & 203 206 conserve_volume_flow_mode, coupling_start_time, cthf, & 204 207 curvature_solution_effects, cycle_mg, damp_level_1d, & … … 212 215 initializing_actions, km_constant, lad_surface, & 213 216 lad_vertical_gradient, lad_vertical_gradient_level, & 214 leaf_surface_concentration, &217 leaf_surface_concentration, limiter_sedimentation, & 215 218 loop_optimization, masking_method, mg_cycles, & 216 219 mg_switch_to_pe0_level, mixing_length_1d, momentum_advec, & 217 mu_constant, mu_constant_value, & 218 netcdf_precision, neutral, ngsrb, nr_surface, & 220 nc, netcdf_precision, neutral, ngsrb, nr_surface, & 219 221 nr_surface_initial_change, nr_vertical_gradient, & 220 222 nr_vertical_gradient_level, nsor, & … … 241 243 topography, topography_grid_convention, top_heatflux, & 242 244 top_momentumflux_u, top_momentumflux_v, top_salinityflux, & 243 turbulen t_inflow, ug_surface, ug_vertical_gradient, &245 turbulence, turbulent_inflow, ug_surface, ug_vertical_gradient, & 244 246 ug_vertical_gradient_level, use_surface_fluxes, & 245 247 use_top_fluxes, use_ug_for_galilei_tr, use_upstream_for_tke, & -
palm/trunk/SOURCE/read_var_list.f90
r1054 r1065 20 20 ! Current revisions: 21 21 ! ------------------ 22 ! +nc, c_sedimentation, limiter_sedimentation, turbulence 23 ! -mu_constant, mu_constant_value 22 24 ! 23 25 ! Former revisions: … … 179 181 !-- Make version number check first 180 182 READ ( 13 ) version_on_file 181 binary_version = '3. 6'183 binary_version = '3.8' 182 184 IF ( TRIM( version_on_file ) /= TRIM( binary_version ) ) THEN 183 185 WRITE( message_string, * ) 'version mismatch concerning control ', & … … 342 344 CASE ( 'canyon_wall_south' ) 343 345 READ ( 13 ) canyon_wall_south 346 CASE ( 'c_sedimentation' ) 347 READ ( 13 ) c_sedimentation 344 348 CASE ( 'cfl_factor' ) 345 349 READ ( 13 ) cfl_factor … … 447 451 CASE ( 'leaf_surface_concentration' ) 448 452 READ ( 13 ) leaf_surface_concentration 453 CASE ( 'limiter_sedimentation' ) 454 READ ( 13 ) limiter_sedimentation 449 455 CASE ( 'loop_optimization' ) 450 456 READ ( 13 ) loop_optimization … … 464 470 CASE ( 'momentum_advec' ) 465 471 READ ( 13 ) momentum_advec 466 CASE ( 'mu_constant' ) 467 READ ( 13 ) mu_constant 468 CASE ( 'mu_constant_value' ) 469 READ ( 13 ) mu_constant_value 472 CASE ( 'nc' ) 473 READ ( 13 ) nc 470 474 CASE ( 'netcdf_precision' ) 471 475 READ ( 13 ) netcdf_precision … … 678 682 CASE ( 'tsc' ) 679 683 READ ( 13 ) tsc 684 CASE ( 'turbulence' ) 685 READ ( 13 ) turbulence 680 686 CASE ( 'turbulent_inflow' ) 681 687 READ ( 13 ) turbulent_inflow -
palm/trunk/SOURCE/time_integration.f90
r1054 r1065 20 20 ! Current revisions: 21 21 ! ----------------- 22 ! exchange of diss (dissipation rate) in case of turbulence = .TRUE. added 22 23 ! 23 24 ! Former revisions: … … 255 256 CALL exchange_horiz( ql_vp, nbgp ) 256 257 ENDIF 257 IF ( wang_kernel ) CALL exchange_horiz( diss, nbgp )258 IF ( wang_kernel .OR. turbulence ) CALL exchange_horiz( diss, nbgp ) 258 259 259 260 CALL cpu_log( log_point(26), 'exchange-horiz-progn', 'stop' ) -
palm/trunk/SOURCE/write_var_list.f90
r1054 r1065 20 20 ! Current revisions: 21 21 ! ----------------- 22 ! +nc, c_sedimentation, turbulence, limiter_sedimentation 23 ! -mu_constant, mu_constant_value 22 24 ! 23 25 ! Former revisions: … … 161 163 162 164 163 binary_version = '3. 6'165 binary_version = '3.8' 164 166 165 167 WRITE ( 14 ) binary_version … … 264 266 WRITE ( 14 ) 'canyon_wall_south ' 265 267 WRITE ( 14 ) canyon_wall_south 268 WRITE ( 14 ) 'c_sedimentation ' 269 WRITE ( 14 ) c_sedimentation 266 270 WRITE ( 14 ) 'cfl_factor ' 267 271 WRITE ( 14 ) cfl_factor … … 368 372 WRITE ( 14 ) 'leaf_surface_concentration ' 369 373 WRITE ( 14 ) leaf_surface_concentration 374 WRITE ( 14 ) 'limiter_sedimentation ' 375 WRITE ( 14 ) limiter_sedimentation 370 376 WRITE ( 14 ) 'loop_optimization ' 371 377 WRITE ( 14 ) loop_optimization … … 384 390 WRITE ( 14 ) 'momentum_advec ' 385 391 WRITE ( 14 ) momentum_advec 386 WRITE ( 14 ) 'mu_constant ' 387 WRITE ( 14 ) mu_constant 388 WRITE ( 14 ) 'mu_constant_value ' 389 WRITE ( 14 ) mu_constant_value 392 WRITE ( 14 ) 'nc ' 393 WRITE ( 14 ) nc 390 394 WRITE ( 14 ) 'netcdf_precision ' 391 395 WRITE ( 14 ) netcdf_precision … … 596 600 WRITE ( 14 ) 'tsc ' 597 601 WRITE ( 14 ) tsc 602 WRITE ( 14 ) 'turbulence ' 603 WRITE ( 14 ) turbulence 598 604 WRITE ( 14 ) 'turbulent_inflow ' 599 605 WRITE ( 14 ) turbulent_inflow
Note: See TracChangeset
for help on using the changeset viewer.
