Line data Source code
1 : !-------------------------------------------------------------------------
2 : ! $Id$
3 : !===============================================================================
4 : module sigma_sqd_w_module
5 :
6 : implicit none
7 :
8 : public :: compute_sigma_sqd_w
9 :
10 : private ! Default scope
11 :
12 : contains
13 :
14 : !=============================================================================
15 17870112 : subroutine compute_sigma_sqd_w( nz, ngrdcol, &
16 17870112 : gamma_Skw_fnc, wp2, thlp2, rtp2, &
17 17870112 : up2, vp2, wpthlp, wprtp, upwp, vpwp, &
18 : l_predict_upwp_vpwp, &
19 17870112 : sigma_sqd_w )
20 :
21 : ! Description:
22 : ! Compute the variable sigma_sqd_w (PDF width parameter).
23 : !
24 : ! The value of sigma_sqd_w is restricted in the ADG1 PDF in order to keep
25 : ! the marginal PDFs of all responder variables (variables that do not set
26 : ! the mixture fraction) valid. The limits on sigma_sqd_w in order to keep
27 : ! the PDF of a responder variable, x, valid are:
28 : !
29 : ! 0 <= sigma_sqd_w <= 1 - <w'x'>^2 / ( <w'^2> * <x'^2> ).
30 : !
31 : ! The overall limits on sigma_sqd_w must be applied based on the most
32 : ! restrictive case so that all Double Gaussian PDF responder variables, x,
33 : ! have realizable PDFs. The overall limits on sigma_sqd_w are:
34 : !
35 : ! 0 <= sigma_sqd_w <= 1 - max( <w'x'>^2 / ( <w'^2> * <x'^2> ), for all x).
36 : !
37 : ! The equation used for sigma_sqd_w is:
38 : !
39 : ! sigma_sqd_w = gamma_Skw_fnc
40 : ! * ( 1 - max( <w'x'>^2 / ( <w'^2> * <x'^2> ), for all x) );
41 : !
42 : ! where 0 <= gamma_Skw_fnc <= 1.
43 :
44 : ! References:
45 : ! Eqn 22 in ``Equations for CLUBB''
46 : !-----------------------------------------------------------------------
47 :
48 : use constants_clubb, only: &
49 : one, & ! Constant(s)
50 : w_tol, &
51 : rt_tol, &
52 : thl_tol, &
53 : one_hundred, &
54 : w_tol_sqd
55 :
56 : use clubb_precision, only: &
57 : core_rknd ! Variable(s)
58 :
59 : implicit none
60 :
61 : ! Input Variables
62 : integer, intent(in) :: &
63 : nz, &
64 : ngrdcol
65 :
66 : real( kind = core_rknd ), dimension(ngrdcol,nz), intent(in) :: &
67 : gamma_Skw_fnc, & ! Gamma as a function of skewness [-]
68 : wp2, & ! Variance of vertical velocity [m^2/s^2]
69 : thlp2, & ! Variance of liquid water potential temp. [K^2]
70 : rtp2, & ! Variance of total water mixing ratio [kg^2/kg^2]
71 : up2, & ! Variance of west-east horizontal velocity [m^2/s^2]
72 : vp2, & ! Variance of south-north horizontal velocity [m^2/s^2]
73 : wpthlp, & ! Flux of liquid water potential temp. [m/s K]
74 : wprtp, & ! Flux of total water mixing ratio [m/s kg/kg]
75 : upwp, & ! Flux of west-east horizontal velocity [m^2/s^2]
76 : vpwp ! Flux of south-north horizontal velocity [m^2/s^2]
77 :
78 : logical, intent(in) :: &
79 : l_predict_upwp_vpwp ! Flag to predict <u'w'> and <v'w'> along with <u> and <v> alongside the
80 : ! advancement of <rt>, <w'rt'>, <thl>, <wpthlp>, <sclr>, and <w'sclr'> in
81 : ! subroutine advance_xm_wpxp. Otherwise, <u'w'> and <v'w'> are still
82 : ! approximated by eddy diffusivity when <u> and <v> are advanced in
83 : ! subroutine advance_windm_edsclrm.
84 :
85 : ! Output Variable
86 : real( kind = core_rknd ), dimension(ngrdcol,nz), intent(out) :: &
87 : sigma_sqd_w ! PDF width parameter [-]
88 :
89 : ! Local Variable
90 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
91 35740224 : max_corr_w_x_sqd ! Max. val. of wpxp^2/(wp2*xp2) for all vars. x [-]
92 :
93 : integer :: i, k
94 :
95 : ! ---- Begin Code ----
96 :
97 : !$acc enter data create( max_corr_w_x_sqd )
98 :
99 : !----------------------------------------------------------------
100 : ! Compute sigma_sqd_w with new formula from Vince
101 : !----------------------------------------------------------------
102 :
103 : ! Find the maximum value of <w'x'>^2 / ( <w'^2> * <x'^2> ) for all
104 : ! variables x that are Double Gaussian PDF responder variables. This
105 : ! includes rt and theta-l. When l_predict_upwp_vpwp is enabled, u and v are
106 : ! also calculated as part of the PDF, and they are included as well.
107 : ! Additionally, when sclr_dim > 0, passive scalars (sclr) are also included.
108 : !$acc parallel loop gang vector collapse(2) default(present)
109 1536829632 : do k = 1, nz
110 25380961632 : do i = 1, ngrdcol
111 47688264000 : max_corr_w_x_sqd(i,k) = max( ( wpthlp(i,k) / ( sqrt( wp2(i,k) * thlp2(i,k) ) &
112 : + one_hundred * w_tol * thl_tol ) )**2, &
113 : ( wprtp(i,k) / ( sqrt( wp2(i,k) * rtp2(i,k) ) &
114 73051355520 : + one_hundred * w_tol * rt_tol ) )**2 )
115 : end do
116 : end do
117 : !$acc end parallel loop
118 :
119 17870112 : if ( l_predict_upwp_vpwp ) then
120 : !$acc parallel loop gang vector collapse(2) default(present)
121 1536829632 : do k = 1, nz
122 25380961632 : do i = 1, ngrdcol
123 47688264000 : max_corr_w_x_sqd(i,k) = max( max_corr_w_x_sqd(i,k), &
124 : ( upwp(i,k) / ( sqrt( up2(i,k) * wp2(i,k) ) &
125 : + one_hundred * w_tol_sqd ) )**2, &
126 : ( vpwp(i,k) / ( sqrt( vp2(i,k) * wp2(i,k) ) &
127 73051355520 : + one_hundred * w_tol_sqd ) )**2 )
128 : end do
129 : end do
130 : !$acc end parallel loop
131 : endif ! l_predict_upwp_vpwp
132 :
133 : ! Calculate the value of sigma_sqd_w
134 : !$acc parallel loop gang vector collapse(2) default(present)
135 1536829632 : do k = 1, nz
136 25380961632 : do i = 1, ngrdcol
137 25363091520 : sigma_sqd_w(i,k) = gamma_Skw_fnc(i,k) * ( one - min( max_corr_w_x_sqd(i,k), one ) )
138 : end do
139 : end do
140 : !$acc end parallel loop
141 :
142 : !$acc exit data delete( max_corr_w_x_sqd )
143 :
144 17870112 : return
145 :
146 : end subroutine compute_sigma_sqd_w
147 :
148 : end module sigma_sqd_w_module
|