updates to MYNN-EDMF and Thompson #120
+1,278
−1,013
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
|
|
@@ -993,6 +993,7 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, & | |
| rainprod, evapprod, & | ||
| #endif | ||
| refl_10cm, diagflag, do_radar_ref, & | ||
| max_hail_diam_sfc, & | ||
| vt_dbz_wt, first_time_step, & | ||
| re_cloud, re_ice, re_snow, & | ||
| has_reqc, has_reqi, has_reqs, & | ||
|
|
@@ -1062,6 +1063,8 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, & | |
| GRAUPELNC, GRAUPELNCV | ||
| REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT):: & | ||
| refl_10cm | ||
| REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT):: & | ||
| max_hail_diam_sfc | ||
| REAL, DIMENSION(ims:ime, kms:kme, jms:jme), OPTIONAL, INTENT(INOUT):: & | ||
| vt_dbz_wt | ||
| LOGICAL, INTENT(IN) :: first_time_step | ||
|
|
@@ -1416,14 +1419,14 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, & | |
| qcten1(k) = 0. | ||
| endif initialize_extended_diagnostics | ||
| enddo | ||
| lsml = lsm(i,j) | ||
| if (is_aerosol_aware .or. merra2_aerosol_aware) then | ||
| do k = kts, kte | ||
| nc1d(k) = nc(i,k,j) | ||
| nwfa1d(k) = nwfa(i,k,j) | ||
| nifa1d(k) = nifa(i,k,j) | ||
| enddo | ||
| else | ||
| do k = kts, kte | ||
| if(lsml == 1) then | ||
| nc1d(k) = Nt_c_l/rho(k) | ||
|
|
@@ -1679,6 +1682,8 @@ SUBROUTINE mp_gt_driver(qv, qc, qr, qi, qs, qg, ni, nr, nc, & | |
| (nsteps>1 .and. istep==nsteps) .or. & | ||
| (nsteps==1 .and. ndt==1)) THEN | ||
|
|
||
| max_hail_diam_sfc(i,j) = hail_mass_99th_percentile(kts, kte, qg1d, t1d, p1d, qv1d) | ||
|
|
||
| !> - Call calc_refl10cm() | ||
|
|
||
| diagflag_present: IF ( PRESENT (diagflag) ) THEN | ||
|
|
@@ -2464,17 +2469,17 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, & | |
| !+---+-----------------------------------------------------------------+ | ||
| !> - Calculate y-intercept, slope values for graupel. | ||
| !+---+-----------------------------------------------------------------+ | ||
| ! do k = kte, kts, -1 | ||
| ! ygra1 = alog10(max(1.E-9, rg(k))) | ||
| ! zans1 = 3.4 + 2./7.*(ygra1+8.) + rand1 | ||
| ! N0_exp = 10.**(zans1) | ||
| ! N0_exp = MAX(DBLE(gonv_min), MIN(N0_exp, DBLE(gonv_max))) | ||
| ! lam_exp = (N0_exp*am_g*cgg(1)/rg(k))**oge1 | ||
| ! lamg = lam_exp * (cgg(3)*ogg2*ogg1)**obmg | ||
| ! ilamg(k) = 1./lamg | ||
| ! N0_g(k) = N0_exp/(cgg(2)*lam_exp) * lamg**cge(2) | ||
| ! enddo | ||
| call graupel_psd_parameters(kts, kte, rand1, rg, ilamg, N0_g) | ||
| endif | ||
|
|
||
| !+---+-----------------------------------------------------------------+ | ||
|
|
@@ -3541,17 +3546,17 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, & | |
| !+---+-----------------------------------------------------------------+ | ||
| !> - Calculate y-intercept, slope values for graupel. | ||
| !+---+-----------------------------------------------------------------+ | ||
| ! do k = kte, kts, -1 | ||
| ! ygra1 = alog10(max(1.E-9, rg(k))) | ||
| ! zans1 = 3.4 + 2./7.*(ygra1+8.) + rand1 | ||
| ! N0_exp = 10.**(zans1) | ||
| ! N0_exp = MAX(DBLE(gonv_min), MIN(N0_exp, DBLE(gonv_max))) | ||
| ! lam_exp = (N0_exp*am_g*cgg(1)/rg(k))**oge1 | ||
| ! lamg = lam_exp * (cgg(3)*ogg2*ogg1)**obmg | ||
| ! ilamg(k) = 1./lamg | ||
| ! N0_g(k) = N0_exp/(cgg(2)*lam_exp) * lamg**cge(2) | ||
| ! enddo | ||
| call graupel_psd_parameters(kts, kte, rand1, rg, ilamg, N0_g) | ||
| endif | ||
|
|
||
| !+---+-----------------------------------------------------------------+ | ||
|
|
@@ -3589,7 +3594,7 @@ subroutine mp_thompson (qv1d, qc1d, qi1d, qr1d, qs1d, qg1d, ni1d, & | |
| !+---+-----------------------------------------------------------------+ ! DROPLET NUCLEATION | ||
| if (clap .gt. eps) then | ||
| if (is_aerosol_aware .or. merra2_aerosol_aware) then | ||
| xnc = MAX(2., activ_ncloud(temp(k), w1d(k)+rand3, nwfa(k), lsml)) | ||
| else | ||
| if(lsml == 1) then | ||
| xnc = Nt_c_l | ||
|
|
@@ -5366,14 +5371,15 @@ end subroutine table_ccnAct | |
| ! TO_DO ITEM: For radiation cooling producing fog, in which case the | ||
| !.. updraft velocity could easily be negative, we could use the temp | ||
| !.. and its tendency to diagnose a pretend postive updraft velocity. | ||
| real function activ_ncloud(Tt, Ww, NCCN, lsm_in) | ||
|
|
||
| implicit none | ||
| REAL, INTENT(IN):: Tt, Ww, NCCN | ||
| INTEGER, INTENT(IN):: lsm_in | ||
| REAL:: n_local, w_local | ||
| INTEGER:: i, j, k, l, m, n | ||
| REAL:: A, B, C, D, t, u, x1, x2, y1, y2, nx, wy, fraction | ||
| REAL:: lower_lim_nuc_frac | ||
|
|
||
| ! ta_Na = (/10.0, 31.6, 100.0, 316.0, 1000.0, 3160.0, 10000.0/) ntb_arc | ||
| ! ta_Ww = (/0.01, 0.0316, 0.1, 0.316, 1.0, 3.16, 10.0, 31.6, 100.0/) ntb_arw | ||
|
|
@@ -5420,6 +5426,14 @@ real function activ_ncloud(Tt, Ww, NCCN) | |
| l = 3 | ||
| m = 2 | ||
|
|
||
| if (lsm_in .eq. 1) then ! land | ||
| lower_lim_nuc_frac = 0. | ||
| else if (lsm_in .eq. 0) then ! water | ||
| lower_lim_nuc_frac = 0.15 | ||
| else | ||
| lower_lim_nuc_frac = 0.15 ! catch-all for anything else | ||
| endif | ||
|
|
||
| A = tnccn_act(i-1,j-1,k,l,m) | ||
| B = tnccn_act(i,j-1,k,l,m) | ||
| C = tnccn_act(i,j,k,l,m) | ||
|
|
@@ -5434,7 +5448,8 @@ real function activ_ncloud(Tt, Ww, NCCN) | |
| ! u = (w_local-ta_Ww(j-1))/(ta_Ww(j)-ta_Ww(j-1)) | ||
|
|
||
| fraction = (1.0-t)*(1.0-u)*A + t*(1.0-u)*B + t*u*C + (1.0-t)*u*D | ||
| fraction = MAX(fraction, lower_lim_nuc_frac) | ||
|
|
||
| ! if (NCCN*fraction .gt. 0.75*Nt_c_max) then | ||
| ! write(*,*) ' DEBUG-GT ', n_local, w_local, Tt, i, j, k | ||
| ! endif | ||
|
|
@@ -5854,6 +5869,7 @@ subroutine calc_effectRad (t1d, p1d, qv1d, qc1d, nc1d, qi1d, ni1d, qs1d, & | |
| endif | ||
| lamc = (nc(k)*am_r*g_ratio(inu_c)/rc(k))**obmr | ||
| re_qc1d(k) = SNGL(0.5D0 * DBLE(3.+inu_c)/lamc) | ||
| if (lsml .ne. 1) re_qc1d(k) = max(re_qc1d(k), 7.0E-6) | ||
joeolson42 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
| enddo | ||
| endif | ||
|
|
||
|
|
@@ -6085,16 +6101,17 @@ subroutine calc_refl10cm (qv1d, qc1d, qr1d, nr1d, qs1d, qg1d, & | |
| !+---+-----------------------------------------------------------------+ | ||
|
|
||
| if (ANY(L_qg .eqv. .true.)) then | ||
| ! do k = kte, kts, -1 | ||
| ! ygra1 = alog10(max(1.E-9, rg(k))) | ||
| ! zans1 = 3.4 + 2./7.*(ygra1+8.) + rand1 | ||
| ! N0_exp = 10.**(zans1) | ||
| ! N0_exp = MAX(DBLE(gonv_min), MIN(N0_exp, DBLE(gonv_max))) | ||
| ! lam_exp = (N0_exp*am_g*cgg(1)/rg(k))**oge1 | ||
| ! lamg = lam_exp * (cgg(3)*ogg2*ogg1)**obmg | ||
| ! ilamg(k) = 1./lamg | ||
| ! N0_g(k) = N0_exp/(cgg(2)*lam_exp) * lamg**cge(2) | ||
| ! enddo | ||
| call graupel_psd_parameters(kts, kte, rand1, rg, ilamg, N0_g) | ||
| endif | ||
|
|
||
| !+---+-----------------------------------------------------------------+ | ||
|
|
@@ -6471,6 +6488,86 @@ SUBROUTINE semi_lagrange_sedim(km,dzl,wwl,rql,precip,pfsan,dt,R1) | |
|
|
||
| END SUBROUTINE semi_lagrange_sedim | ||
|
|
||
| !>\ingroup aathompson | ||
| !! @brief Calculates graupel size distribution parameters | ||
| !! | ||
| !! Calculates graupel intercept and slope parameters for | ||
| !! for a vertical column | ||
| !! | ||
| !! @param[in] kts integer start index for vertical column | ||
| !! @param[in] kte integer end index for vertical column | ||
| !! @param[in] rand1 real random number for stochastic physics | ||
| !! @param[in] rg real array, size(kts:kte) for graupel mass concentration [kg m^3] | ||
| !! @param[out] ilamg double array, size(kts:kte) for inverse graupel slope parameter [m] | ||
| !! @param[out] N0_g double array, size(kts:kte) for graupel intercept paramter [m-4] | ||
| subroutine graupel_psd_parameters(kts, kte, rand1, rg, ilamg, N0_g) | ||
|
|
||
| implicit none | ||
|
|
||
| integer, intent(in) :: kts, kte | ||
| real, intent(in) :: rand1 | ||
| real, intent(in) :: rg(:) | ||
| double precision, intent(out) :: ilamg(:), N0_g(:) | ||
|
|
||
| integer :: k | ||
| real :: ygra1, zans1 | ||
| double precision :: N0_exp, lam_exp, lamg | ||
|
|
||
| do k = kte, kts, -1 | ||
| ygra1 = alog10(max(1.e-9, rg(k))) | ||
| zans1 = 3.4 + 2./7.*(ygra1+8.) + rand1 | ||
| N0_exp = 10.**(zans1) | ||
| N0_exp = max1(dble(gonv_min), min(N0_exp, dble(gonv_max))) | ||
|
There was a problem hiding this comment. @AndersJensen-NOAA We settled on changing this to "max" or "dmax1", right? "max1" should be returning an integer that is then cast to double precision, which I don't think is what you want, right? There was a problem hiding this comment. Nevermind, I see that you already did this and this is outdated. |
||
| lam_exp = (N0_exp*am_g*cgg(1)/rg(k))**oge1 | ||
| lamg = lam_exp * (cgg(3)*ogg2*ogg1)**obmg | ||
| ilamg(k) = 1./lamg | ||
| N0_g(k) = N0_exp/(cgg(2)*lam_exp) * lamg**cge(2) | ||
| enddo | ||
|
|
||
| end subroutine graupel_psd_parameters | ||
|
|
||
| !>\ingroup aathompson | ||
| !! @brief Calculates graupel/hail maximum diameter | ||
| !! | ||
| !! Calculates graupel/hail maximum diameter (currently the 99th percentile of mass distribtuion) | ||
| !! for a vertical column | ||
| !! | ||
| !! @param[in] kts integer start index for vertical column | ||
| !! @param[in] kte integer end index for vertical column | ||
| !! @param[in] qg real array, size(kts:kte) for graupel mass mixing ratio [kg kg^-1] | ||
| !! @param[in] temperature double array, size(kts:kte) temperature [K] | ||
| !! @param[in] pressure double array, size(kts:kte) pressure [Pa] | ||
| !! @param[in] qv real array, size(kts:kte) water vapor mixing ratio [kg kg^-1] | ||
| !! @param[out] max_hail_diam real maximum hail diameter [m] | ||
| function hail_mass_99th_percentile(kts, kte, qg, temperature, pressure, qv) result(max_hail_diam) | ||
|
|
||
| implicit none | ||
|
|
||
| integer, intent(in) :: kts, kte | ||
| real, intent(in) :: qg(:), temperature(:), pressure(:), qv(:) | ||
| real :: max_hail_diam | ||
|
|
||
| integer :: k | ||
| real :: rho(kts:kte), rg(kts:kte), max_hail_column(kts:kte) | ||
| double precision :: ilamg(kts:kte), N0_g(kts:kte) | ||
| real, parameter :: random_number = 0. | ||
|
|
||
| max_hail_column = 0. | ||
| rg = 0. | ||
| do k = kts, kte | ||
| rho(k) = 0.622*pressure(k)/(R*temperature(k)*(max(1.e-10, qv(k))+0.622)) | ||
| if (qg(k) .gt. R1) then | ||
| rg(k) = qg(k)*rho(k) | ||
| endif | ||
| enddo | ||
|
|
||
| call graupel_psd_parameters(kts, kte, random_number, rg, ilamg, N0_g) | ||
|
|
||
| where(rg .gt. 1.e-9) max_hail_column = 10.05 * ilamg | ||
| max_hail_diam = max_hail_column(kts) | ||
|
|
||
| end function hail_mass_99th_percentile | ||
|
|
||
| !+---+-----------------------------------------------------------------+ | ||
| !+---+-----------------------------------------------------------------+ | ||
| END MODULE module_mp_thompson | ||
|
|
||
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Oops, something went wrong.
Oops, something went wrong.
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The comment lines 2472-2481, 3549-3558, and 6104-6113 can be removed