Line data Source code
1 : !-------------------------------------------------------------------------
2 : !$Id$
3 : !===============================================================================
4 : module Skx_module
5 :
6 : implicit none
7 :
8 : private ! Default Scope
9 :
10 : public :: Skx_func, &
11 : LG_2005_ansatz, &
12 : xp3_LG_2005_ansatz
13 :
14 : contains
15 :
16 : !-----------------------------------------------------------------------------
17 4588272 : subroutine Skx_func( nz, ngrdcol, xp2, xp3, &
18 : x_tol, Skw_denom_coef, Skw_max_mag, &
19 4588272 : Skx )
20 :
21 : ! Description:
22 : ! Calculate the skewness of x
23 :
24 : ! References:
25 : ! None
26 : !-----------------------------------------------------------------------
27 :
28 : use clubb_precision, only: &
29 : core_rknd ! Variable(s)
30 :
31 : implicit none
32 :
33 : integer, intent(in) :: &
34 : nz, &
35 : ngrdcol
36 :
37 : ! External
38 : intrinsic :: min, max
39 :
40 : ! Parameter Constants
41 : ! Whether to apply clipping to the final result
42 : logical, parameter :: &
43 : l_clipping_kluge = .false.
44 :
45 : ! Input Variables
46 : real( kind = core_rknd ), dimension(ngrdcol,nz), intent(in) :: &
47 : xp2, & ! <x'^2> [(x units)^2]
48 : xp3 ! <x'^3> [(x units)^3]
49 :
50 : real( kind = core_rknd ), intent(in) :: &
51 : x_tol, & ! x tolerance value [(x units)]
52 : Skw_denom_coef, & ! CLUBB tunable parameter Skw_denom_coef [-]
53 : Skw_max_mag ! Max magnitude of skewness [-]
54 :
55 : ! Output Variable
56 : real( kind = core_rknd ), dimension(ngrdcol,nz), intent(out) :: &
57 : Skx ! Skewness of x [-]
58 :
59 : ! Local Variable
60 : real( kind = core_rknd ) :: &
61 : Skx_denom_tol
62 :
63 : integer :: i, k
64 :
65 : ! ---- Begin Code ----
66 :
67 : !$acc data copyin( xp2, xp3 ) &
68 : !$acc copyout( Skx )
69 :
70 4588272 : Skx_denom_tol = Skw_denom_coef * x_tol**2
71 :
72 : !Skx = xp3 / ( max( xp2, x_tol**two ) )**three_halves
73 : ! Calculation of skewness to help reduce the sensitivity of this value to
74 : ! small values of xp2.
75 : !$acc parallel loop gang vector collapse(2) default(present)
76 394591392 : do k = 1, nz
77 6516733392 : do i = 1, ngrdcol
78 6512145120 : Skx(i,k) = xp3(i,k) * sqrt( xp2(i,k) + Skx_denom_tol )**(-3)
79 : end do
80 : end do
81 : !$acc end parallel loop
82 :
83 : ! This is no longer needed since clipping is already
84 : ! imposed on wp2 and wp3 elsewhere in the code
85 :
86 : ! I turned clipping on in this local copy since thlp3 and rtp3 are not clipped
87 : if ( l_clipping_kluge ) then
88 : !$acc parallel loop gang vector collapse(2) default(present)
89 : do k = 1, nz
90 : do i = 1, ngrdcol
91 : Skx(i,k) = min( max( Skx(i,k), -Skw_max_mag ), Skw_max_mag )
92 : end do
93 : end do
94 : !$acc end parallel loop
95 : end if
96 :
97 : !$acc end data
98 :
99 4588272 : return
100 :
101 : end subroutine Skx_func
102 :
103 : !-----------------------------------------------------------------------------
104 1411776 : subroutine LG_2005_ansatz( nz, ngrdcol, Skw, wpxp, wp2, &
105 1411776 : xp2, beta, sigma_sqd_w, x_tol, &
106 1411776 : Skx )
107 :
108 : ! Description:
109 : ! Calculate the skewness of x using the diagnostic ansatz of Larson and
110 : ! Golaz (2005).
111 :
112 : ! References:
113 : ! Vincent E. Larson and Jean-Christophe Golaz, 2005: Using Probability
114 : ! Density Functions to Derive Consistent Closure Relationships among
115 : ! Higher-Order Moments. Mon. Wea. Rev., 133, 1023–1042.
116 : !-----------------------------------------------------------------------
117 :
118 : use constants_clubb, only: &
119 : one, & ! Variable(s)
120 : w_tol_sqd
121 :
122 : use clubb_precision, only: &
123 : core_rknd ! Variable(s)
124 :
125 : implicit none
126 :
127 : !-------------------------- Input Variables --------------------------
128 : integer, intent(in) :: &
129 : nz, &
130 : ngrdcol
131 :
132 : real( kind = core_rknd ), dimension(ngrdcol,nz), intent(in) :: &
133 : Skw, & ! Skewness of w [-]
134 : wpxp, & ! Turbulent flux of x [m/s (x units)]
135 : wp2, & ! Variance of w [m^2/s^2]
136 : xp2, & ! Variance of x [(x units)^2]
137 : sigma_sqd_w ! Normalized variance of w [-]
138 :
139 : real( kind = core_rknd ), intent(in) :: &
140 : beta, & ! Tunable parameter [-]
141 : x_tol ! Minimum tolerance of x [(x units)]
142 :
143 : !-------------------------- Output Variable --------------------------
144 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
145 : Skx ! Skewness of x [-]
146 :
147 : !-------------------------- Local Variables --------------------------
148 : real( kind = core_rknd ) :: &
149 : nrmlzd_corr_wx, & ! Normalized correlation of w and x [-]
150 : nrmlzd_Skw ! Normalized skewness of w [-]
151 :
152 : integer :: i, k
153 :
154 : !--------------------------Begin Code --------------------------
155 :
156 : ! weberjk, 8-July 2015. Commented this out for now. cgils was failing during some tests.
157 :
158 : ! Larson and Golaz (2005) eq. 16
159 : !$acc parallel loop gang vector collapse(2) default(present)
160 121412736 : do k = 1, nz
161 2005148736 : do i = 1, ngrdcol
162 : nrmlzd_corr_wx = &
163 3767472000 : wpxp(i,k) / sqrt( max( wp2(i,k), w_tol_sqd ) &
164 3767472000 : * max( xp2(i,k), x_tol**2 ) * ( one - sigma_sqd_w(i,k) ) )
165 :
166 : ! Larson and Golaz (2005) eq. 11
167 1883736000 : nrmlzd_Skw = Skw(i,k) / ( ( one - sigma_sqd_w(i,k)) * sqrt( one - sigma_sqd_w(i,k) ) )
168 :
169 : ! Larson and Golaz (2005) eq. 33
170 : Skx(i,k) = nrmlzd_Skw * nrmlzd_corr_wx &
171 2003736960 : * ( beta + ( one - beta ) * nrmlzd_corr_wx**2 )
172 : end do
173 : end do
174 : !$acc end parallel loop
175 :
176 1411776 : return
177 :
178 : end subroutine LG_2005_ansatz
179 :
180 : !-----------------------------------------------------------------------------
181 1411776 : subroutine xp3_LG_2005_ansatz( nz, ngrdcol, Skw_zt, wpxp_zt, wp2_zt, &
182 1411776 : xp2_zt, sigma_sqd_w_zt, &
183 : beta, Skw_denom_coef, x_tol, &
184 1411776 : xp3 )
185 : ! Description:
186 : ! Calculate <x'^3> after calculating the skewness of x using the ansatz of
187 : ! Larson and Golaz (2005).
188 :
189 : ! References:
190 : !-----------------------------------------------------------------------
191 :
192 : use grid_class, only: &
193 : grid ! Type
194 :
195 : use clubb_precision, only: &
196 : core_rknd ! Variable(s)
197 :
198 : implicit none
199 :
200 : !-------------------------- Input Variables--------------------------
201 : integer, intent(in) :: &
202 : nz, &
203 : ngrdcol
204 :
205 : real( kind = core_rknd ), dimension(ngrdcol,nz), intent(in) :: &
206 : Skw_zt, & ! Skewness of w on thermodynamic levels [-]
207 : wpxp_zt, & ! Flux of x (interp. to t-levs.) [m/s(x units)]
208 : wp2_zt, & ! Variance of w (interp. to t-levs.) [m^2/s^2]
209 : xp2_zt, & ! Variance of x (interp. to t-levs.) [(x units)^2]
210 : sigma_sqd_w_zt ! Normalized variance of w (interp. to t-levs.) [-]
211 :
212 : real( kind = core_rknd ), intent(in) :: &
213 : beta, & ! CLUBB tunable parameter beta [-]
214 : Skw_denom_coef, & ! CLUBB tunable parameter Skw_denom_coef [-]
215 : x_tol ! Minimum tolerance of x [(x units)]
216 :
217 : !-------------------------- Return Variable --------------------------
218 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
219 : xp3 ! <x'^3> (thermodynamic levels) [(x units)^3]
220 :
221 : !-------------------------- Local Variable --------------------------
222 : real( kind = core_rknd ), dimension(ngrdcol,nz) :: &
223 2823552 : Skx_zt ! Skewness of x on thermodynamic levels [-]
224 :
225 : real( kind = core_rknd ) :: &
226 : Skx_denom_tol
227 :
228 : integer :: i, k
229 :
230 : !-------------------------- Begin Code --------------------------
231 :
232 : !$acc data create(Skx_zt) &
233 : !$acc copyin( Skw_zt, wpxp_zt, wp2_zt, xp2_zt, sigma_sqd_w_zt ) &
234 : !$acc copyout( xp3 )
235 :
236 : ! Calculate skewness of x using the ansatz of LG05.
237 : call LG_2005_ansatz( nz, ngrdcol, Skw_zt, wpxp_zt, wp2_zt, &
238 : xp2_zt, beta, sigma_sqd_w_zt, x_tol, &
239 1411776 : Skx_zt )
240 :
241 1411776 : Skx_denom_tol = Skw_denom_coef * x_tol**2
242 :
243 : ! Calculate <x'^3> using the reverse of the special sensitivity reduction
244 : ! formula in function Skx_func above.
245 : !$acc parallel loop gang vector collapse(2) default(present)
246 121412736 : do k = 1, nz
247 2005148736 : do i = 1, ngrdcol
248 3767472000 : xp3(i,k) = Skx_zt(i,k) * ( xp2_zt(i,k) + Skx_denom_tol ) &
249 5771208960 : * sqrt( xp2_zt(i,k) + Skx_denom_tol )
250 : end do
251 : end do
252 : !$acc end parallel loop
253 :
254 : !$acc end data
255 :
256 1411776 : return
257 :
258 : end subroutine xp3_LG_2005_ansatz
259 :
260 : !-----------------------------------------------------------------------------
261 :
262 : end module Skx_module
|