- Timestamp:
- Mar 25, 2014 12:21:40 PM (11 years ago)
- Location:
- palm/trunk/SOURCE
- Files:
-
- 2 edited
-
init_cloud_physics.f90 (modified) (3 diffs)
-
microphysics.f90 (modified) (18 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/init_cloud_physics.f90
r1323 r1334 20 20 ! Current revisions: 21 21 ! ----------------- 22 ! 22 ! Bugfix: REAL constants provided with KIND-attribute 23 23 ! 24 24 ! Former revisions: … … 118 118 !-- t / pt : ratio of actual and potential temperature (t_d_pt) 119 119 !-- p_0(z) : vertical profile of the hydrostatic pressure (hyp) 120 t_surface = pt_surface * ( surface_pressure / 1000.0 )**0.286 120 t_surface = pt_surface * ( surface_pressure / 1000.0 )**0.286_wp 121 121 DO k = nzb, nzt+1 122 122 ! … … 129 129 hyp(k) = surface_pressure * 100.0 * & 130 130 ( (t_surface - g/cp * zu(k)) / t_surface )**(1.0_wp/0.286_wp) 131 pt_d_t(k) = ( 100000.0 / hyp(k) )**0.286 131 pt_d_t(k) = ( 100000.0 / hyp(k) )**0.286_wp 132 132 t_d_pt(k) = 1.0 / pt_d_t(k) 133 133 hyrho(k) = hyp(k) / ( r_d * t_d_pt(k) * pt_init(k) ) -
palm/trunk/SOURCE/microphysics.f90
r1323 r1334 20 20 ! Current revisions: 21 21 ! ------------------ 22 ! 22 ! Bugfix: REAL constants provided with KIND-attribute 23 23 ! 24 24 ! Former revisions: … … 361 361 !-- t / pt : ratio of actual and potential temperature (t_d_pt) 362 362 !-- p_0(z) : vertical profile of the hydrostatic pressure (hyp) 363 t_surface = pt_surface * ( surface_pressure / 1000.0 )**0.286 363 t_surface = pt_surface * ( surface_pressure / 1000.0 )**0.286_wp 364 364 DO k = nzb, nzt+1 365 365 hyp(k) = surface_pressure * 100.0 * & 366 ( (t_surface - g/cp * zu(k)) / t_surface )**(1.0 /0.286)367 pt_d_t(k) = ( 100000.0 / hyp(k) )**0.286 366 ( (t_surface - g/cp * zu(k)) / t_surface )**(1.0_wp/0.286_wp) 367 pt_d_t(k) = ( 100000.0 / hyp(k) )**0.286_wp 368 368 t_d_pt(k) = 1.0 / pt_d_t(k) 369 369 hyrho(k) = hyp(k) / ( r_d * t_d_pt(k) * pt_init(k) ) … … 511 511 !-- Universal function for autoconversion process 512 512 !-- (Seifert and Beheng, 2006): 513 phi_au = 600.0 * tau_cloud**0.68 * ( 1.0 - tau_cloud**0.68)**3513 phi_au = 600.0 * tau_cloud**0.68_wp * ( 1.0 - tau_cloud**0.68_wp )**3 514 514 ! 515 515 !-- Shape parameter of gamma distribution (Geoffroy et al., 2010): 516 516 !-- (Use constant nu_c = 1.0 instead?) 517 nu_c = 1.0 !MAX( 0.0 , 1580.0 * hyrho(k) * qc(k,j,i) - 0.28)517 nu_c = 1.0 !MAX( 0.0_wp, 1580.0 * hyrho(k) * qc(k,j,i) - 0.28_wp ) 518 518 ! 519 519 !-- Mean weight of cloud droplets: … … 525 525 ! 526 526 !-- Weight averaged radius of cloud droplets: 527 rc = 0.5 * ( xc * dpirho_l )**( 1.0 / 3.0_wp )527 rc = 0.5 * ( xc * dpirho_l )**( 1.0_wp / 3.0_wp ) 528 528 529 529 alpha_cc = ( a_1 + a_2 * nu_c ) / ( 1.0 + a_3 * nu_c ) … … 532 532 ! 533 533 !-- Mixing length (neglecting distance to ground and stratification) 534 l_mix = ( dx * dy * dzu(k) )**( 1.0 / 3.0_wp )534 l_mix = ( dx * dy * dzu(k) )**( 1.0_wp / 3.0_wp ) 535 535 ! 536 536 !-- Limit dissipation rate according to Seifert, Nuijens and 537 537 !-- Stevens (2010) 538 epsilon = MIN( 0.06 , diss(k,j,i) )538 epsilon = MIN( 0.06_wp, diss(k,j,i) ) 539 539 ! 540 540 !-- Compute Taylor-microscale Reynolds number: 541 re_lambda = 6.0 / 11.0 * ( l_mix / c_const )**( 2.0 / 3.0_wp ) * &542 SQRT( 15.0 / kin_vis_air ) * epsilon**( 1.0 / 6.0_wp )541 re_lambda = 6.0 / 11.0 * ( l_mix / c_const )**( 2.0_wp / 3.0_wp ) * & 542 SQRT( 15.0 / kin_vis_air ) * epsilon**( 1.0_wp / 6.0_wp ) 543 543 ! 544 544 !-- The factor of 1.0E4 is needed to convert the dissipation rate 545 545 !-- from m2 s-3 to cm2 s-3. 546 k_au = k_au * ( 1.0 +&547 epsilon * 1.0E4 * ( re_lambda * 1.0E-3 )**0.25 *&546 k_au = k_au * ( 1.0_wp + & 547 epsilon * 1.0E4 * ( re_lambda * 1.0E-3 )**0.25_wp * & 548 548 ( alpha_cc * EXP( -1.0 * ( ( rc - r_cc ) / & 549 549 sigma_cc )**2 ) + beta_cc ) ) … … 552 552 !-- Autoconversion rate (Seifert and Beheng, 2006): 553 553 autocon = k_au * ( nu_c + 2.0 ) * ( nu_c + 4.0 ) / & 554 ( nu_c + 1.0 )**2 * qc_1d(k)**2 * xc**2* &555 ( 1.0 + phi_au / ( 1.0 - tau_cloud )**2 ) * &554 ( nu_c + 1.0 )**2.0_wp * qc_1d(k)**2.0_wp * xc**2.0_wp * & 555 ( 1.0 + phi_au / ( 1.0 - tau_cloud )**2.0_wp ) * & 556 556 rho_surface 557 557 autocon = MIN( autocon, qc_1d(k) / dt_micro ) … … 605 605 !-- (Seifert and Beheng, 2001): 606 606 phi_ac = tau_cloud / ( tau_cloud + 5.0E-5 ) 607 phi_ac = ( phi_ac**2 )**2607 phi_ac = ( phi_ac**2.0_wp )**2.0_wp 608 608 ! 609 609 !-- Parameterized turbulence effects on autoconversion (Seifert, … … 612 612 IF ( turbulence ) THEN 613 613 k_cr = k_cr0 * ( 1.0 + 0.05 * & 614 MIN( 600.0 , diss(k,j,i) * 1.0E4 )**0.25)614 MIN( 600.0_wp, diss(k,j,i) * 1.0E4 )**0.25_wp ) 615 615 ELSE 616 616 k_cr = k_cr0 … … 667 667 ! 668 668 !-- Weight averaged diameter of rain drops: 669 dr = ( hyrho(k) * qr_1d(k) / nr_1d(k) * dpirho_l )**( 1.0 / 3.0_wp )669 dr = ( hyrho(k) * qr_1d(k) / nr_1d(k) * dpirho_l )**( 1.0_wp / 3.0_wp ) 670 670 ! 671 671 !-- Collisional breakup rate (Seifert, 2008): … … 750 750 ! 751 751 !-- Supersaturation: 752 sat = MIN( 0.0 , ( q_1d(k) - qr_1d(k) - qc_1d(k) ) / q_s - 1.0 )752 sat = MIN( 0.0_wp, ( q_1d(k) - qr_1d(k) - qc_1d(k) ) / q_s - 1.0 ) 753 753 ! 754 754 !-- Actual temperature: … … 763 763 ! 764 764 !-- Weight averaged diameter of rain drops: 765 dr = ( xr * dpirho_l )**( 1.0 / 3.0_wp )765 dr = ( xr * dpirho_l )**( 1.0_wp / 3.0_wp ) 766 766 ! 767 767 !-- Compute ventilation factor and intercept parameter … … 775 775 !-- Slope parameter of gamma distribution (Seifert, 2008): 776 776 lambda_r = ( ( mu_r + 3.0 ) * ( mu_r + 2.0 ) * & 777 ( mu_r + 1.0 ) )**( 1.0 / 3.0_wp ) / dr777 ( mu_r + 1.0 ) )**( 1.0_wp / 3.0_wp ) / dr 778 778 779 779 mu_r_2 = mu_r + 2.0 … … 787 787 ( lambda_r / & 788 788 ( c_term + lambda_r ) )**mu_r_5d2 - & 789 0.125 * ( b_term / a_term )**2 *&789 0.125 * ( b_term / a_term )**2.0_wp * & 790 790 ( lambda_r / & 791 791 ( 2.0 * c_term + lambda_r ) )**mu_r_5d2 - & 792 0.0625 * ( b_term / a_term )**3 *&792 0.0625 * ( b_term / a_term )**3.0_wp * & 793 793 ( lambda_r / & 794 794 ( 3.0 * c_term + lambda_r ) )**mu_r_5d2 - & 795 0.0390625 * ( b_term / a_term )**4 *&795 0.0390625 * ( b_term / a_term )**4.0_wp * & 796 796 ( lambda_r / & 797 797 ( 4.0 * c_term + lambda_r ) )**mu_r_5d2 ) … … 852 852 ! 853 853 !-- Sedimentation of cloud droplets (Heus et al., 2010): 854 sed_qc_const = k_st * ( 3.0 / ( 4.0 * pi * rho_l ))**( 2.0 / 3.0_wp ) * &854 sed_qc_const = k_st * ( 3.0 / ( 4.0 * pi * rho_l ))**( 2.0_wp / 3.0_wp ) * & 855 855 EXP( 5.0 * LOG( sigma_gc )**2 ) 856 856 … … 859 859 DO k = nzt, nzb_s_inner(j,i)+1, -1 860 860 IF ( qc_1d(k) > eps_sb ) THEN 861 sed_qc(k) = sed_qc_const * nc_1d(k)**( -2.0 / 3.0_wp ) * &862 ( qc_1d(k) * hyrho(k) )**( 5.0 / 3.0_wp )861 sed_qc(k) = sed_qc_const * nc_1d(k)**( -2.0_wp / 3.0_wp ) * & 862 ( qc_1d(k) * hyrho(k) )**( 5.0_wp / 3.0_wp ) 863 863 ELSE 864 864 sed_qc(k) = 0.0 … … 943 943 ! 944 944 !-- Weight averaged diameter of rain drops: 945 dr = ( hyrho(k) * qr_1d(k) / nr_1d(k) * dpirho_l )**( 1.0 / 3.0_wp )945 dr = ( hyrho(k) * qr_1d(k) / nr_1d(k) * dpirho_l )**( 1.0_wp / 3.0_wp ) 946 946 ! 947 947 !-- Shape parameter of gamma distribution (Milbrandt and Yau, 2005; … … 951 951 !-- Slope parameter of gamma distribution (Seifert, 2008): 952 952 lambda_r = ( ( mu_r + 3.0 ) * ( mu_r + 2.0 ) * & 953 ( mu_r + 1.0 ) )**( 1.0 / 3.0_wp ) / dr954 955 w_nr(k) = MAX( 0.1 , MIN( 20.0, a_term - b_term * ( 1.0 +&953 ( mu_r + 1.0 ) )**( 1.0_wp / 3.0_wp ) / dr 954 955 w_nr(k) = MAX( 0.1_wp, MIN( 20.0_wp, a_term - b_term * ( 1.0 + & 956 956 c_term / lambda_r )**( -1.0 * ( mu_r + 1.0 ) ) ) ) 957 w_qr(k) = MAX( 0.1 , MIN( 20.0, a_term - b_term * ( 1.0 +&957 w_qr(k) = MAX( 0.1_wp, MIN( 20.0_wp, a_term - b_term * ( 1.0 + & 958 958 c_term / lambda_r )**( -1.0 * ( mu_r + 4.0 ) ) ) ) 959 959 ELSE … … 1109 1109 tmp = x_gamm + 5.5 1110 1110 tmp = ( x_gamm + 0.5 ) * LOG( tmp ) - tmp 1111 ser = 1.000000000190015 1111 ser = 1.000000000190015_wp 1112 1112 1113 1113 DO j = 1, 6
Note: See TracChangeset
for help on using the changeset viewer.
