forked from ESCOMP/CLUBB_CESM
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Skx_module.F90
210 lines (151 loc) · 6.4 KB
/
Skx_module.F90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
!-------------------------------------------------------------------------
!$Id$
!===============================================================================
module Skx_module
implicit none
private ! Default Scope
public :: Skx_func, &
LG_2005_ansatz, &
xp3_LG_2005_ansatz
contains
!-----------------------------------------------------------------------------
elemental function Skx_func( xp2, xp3, x_tol ) &
result( Skx )
! Description:
! Calculate the skewness of x
! References:
! None
!-----------------------------------------------------------------------
use constants_clubb, only: &
three_halves ! 3/2
use parameters_tunable, only: &
Skw_denom_coef, & ! Variable(s)
Skw_max_mag ! Max magnitude of skewness
use clubb_precision, only: &
core_rknd ! Variable(s)
implicit none
! External
intrinsic :: min, max
! Parameter Constants
! Whether to apply clipping to the final result
logical, parameter :: &
l_clipping_kluge = .false.
! Input Variables
real( kind = core_rknd ), intent(in) :: &
xp2, & ! <x'^2> [(x units)^2]
xp3, & ! <x'^3> [(x units)^3]
x_tol ! x tolerance value [(x units)]
! Output Variable
real( kind = core_rknd ) :: &
Skx ! Skewness of x [-]
! ---- Begin Code ----
!Skx = xp3 / ( max( xp2, x_tol**two ) )**three_halves
! Calculation of skewness to help reduce the sensitivity of this value to
! small values of xp2.
Skx = xp3 / ( xp2 + Skw_denom_coef * x_tol**2 )**three_halves
! This is no longer needed since clipping is already
! imposed on wp2 and wp3 elsewhere in the code
! I turned clipping on in this local copy since thlp3 and rtp3 are not clipped
if ( l_clipping_kluge ) then
Skx = min( max( Skx, -Skw_max_mag ), Skw_max_mag )
end if
return
end function Skx_func
!-----------------------------------------------------------------------------
elemental function LG_2005_ansatz( Skw, wpxp, wp2, &
xp2, beta, sigma_sqd_w, x_tol ) &
result( Skx )
! Description:
! Calculate the skewness of x using the diagnostic ansatz of Larson and
! Golaz (2005).
! References:
! Vincent E. Larson and Jean-Christophe Golaz, 2005: Using Probability
! Density Functions to Derive Consistent Closure Relationships among
! Higher-Order Moments. Mon. Wea. Rev., 133, 1023–1042.
!-----------------------------------------------------------------------
use constants_clubb, only: &
three_halves, & ! Variable(s)
one, &
w_tol_sqd
use clubb_precision, only: &
core_rknd ! Variable(s)
implicit none
! External
intrinsic :: sqrt
! Input Variables
real( kind = core_rknd ), intent(in) :: &
Skw, & ! Skewness of w [-]
wpxp, & ! Turbulent flux of x [m/s (x units)]
wp2, & ! Variance of w [m^2/s^2]
xp2, & ! Variance of x [(x units)^2]
beta, & ! Tunable parameter [-]
sigma_sqd_w, & ! Normalized variance of w [-]
x_tol ! Minimum tolerance of x [(x units)]
! Output Variable
real( kind = core_rknd ) :: &
Skx ! Skewness of x [-]
! Local Variables
real( kind = core_rknd ) :: &
nrmlzd_corr_wx, & ! Normalized correlation of w and x [-]
nrmlzd_Skw ! Normalized skewness of w [-]
! ---- Begin Code ----
! weberjk, 8-July 2015. Commented this out for now. cgils was failing during some tests.
! Larson and Golaz (2005) eq. 16
nrmlzd_corr_wx &
= ( wpxp / ( sqrt( max( wp2, w_tol_sqd ) ) &
* sqrt( max( xp2, x_tol**2 ) ) ) ) &
/ sqrt( one - sigma_sqd_w )
! Larson and Golaz (2005) eq. 11
nrmlzd_Skw = Skw * ( one - sigma_sqd_w)**(-three_halves)
! Larson and Golaz (2005) eq. 33
Skx = nrmlzd_Skw * nrmlzd_corr_wx &
* ( beta + ( one - beta ) * nrmlzd_corr_wx**2 )
return
end function LG_2005_ansatz
!-----------------------------------------------------------------------------
function xp3_LG_2005_ansatz( Skw_zt, wpxp_zt, wp2_zt, &
xp2_zt, sigma_sqd_w_zt, x_tol ) &
result( xp3 )
! Description:
! Calculate <x'^3> after calculating the skewness of x using the ansatz of
! Larson and Golaz (2005).
! References:
!-----------------------------------------------------------------------
use grid_class, only: &
gr ! Variable Type
use constants_clubb, only: &
three_halves ! Variable(s)
use parameters_tunable, only: &
beta, & ! Variable(s)
Skw_denom_coef
use clubb_precision, only: &
core_rknd ! Variable(s)
implicit none
! External
intrinsic :: sqrt
! Input Variables
real( kind = core_rknd ), dimension(gr%nz), intent(in) :: &
Skw_zt, & ! Skewness of w on thermodynamic levels [-]
wpxp_zt, & ! Flux of x (interp. to t-levs.) [m/s(x units)]
wp2_zt, & ! Variance of w (interp. to t-levs.) [m^2/s^2]
xp2_zt, & ! Variance of x (interp. to t-levs.) [(x units)^2]
sigma_sqd_w_zt ! Normalized variance of w (interp. to t-levs.) [-]
real( kind = core_rknd ), intent(in) :: &
x_tol ! Minimum tolerance of x [(x units)]
! Return Variable
real( kind = core_rknd ), dimension(gr%nz) :: &
xp3 ! <x'^3> (thermodynamic levels) [(x units)^3]
! Local Variable
real( kind = core_rknd ), dimension(gr%nz) :: &
Skx_zt ! Skewness of x on thermodynamic levels [-]
! Calculate skewness of x using the ansatz of LG05.
Skx_zt(1:gr%nz) &
= LG_2005_ansatz( Skw_zt(1:gr%nz), wpxp_zt(1:gr%nz), wp2_zt(1:gr%nz), &
xp2_zt(1:gr%nz), beta, sigma_sqd_w_zt(1:gr%nz), x_tol )
! Calculate <x'^3> using the reverse of the special sensitivity reduction
! formula in function Skx_func above.
xp3 = Skx_zt * ( xp2_zt + Skw_denom_coef * x_tol**2 )**three_halves
return
end function xp3_LG_2005_ansatz
!-----------------------------------------------------------------------------
end module Skx_module