GCC Code Coverage Report


Directory: ./
File: phys/wx_pbl_mod.f90
Date: 2022-01-11 19:19:34
Exec Total Coverage
Lines: 0 348 0.0%
Branches: 0 176 0.0%

Line Branch Exec Source
1 MODULE wx_pbl_mod
2 !
3 ! Split Planetary Boundary Layer
4 !
5 ! This module manages the splitting of the boundary layer between two regions; the (w)
6 ! region (inside cold pools) and the (x) region (outside cold pools)
7 !
8 USE dimphy
9
10 IMPLICIT NONE
11
12 CONTAINS
13 !
14 !****************************************************************************************
15 !
16 SUBROUTINE wx_pbl0_merge(knon, ypplay, ypaprs, &
17 sigw, dTs_forcing, dqs_forcing, &
18 yt_x, yt_w, yq_x, yq_w, &
19 yu_x, yu_w, yv_x, yv_w, &
20 ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w, &
21 ycdragm_x, ycdragm_w, &
22 AcoefT_x, AcoefT_w, AcoefQ_x, AcoefQ_w, &
23 AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w, &
24 BcoefT_x, BcoefT_w, BcoefQ_x, BcoefQ_w, &
25 BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w, &
26 AcoefT, AcoefQ, AcoefU, AcoefV, &
27 BcoefT, BcoefQ, BcoefU, BcoefV, &
28 ycdragh, ycdragq, ycdragm, &
29 yt1, yq1, yu1, yv1 &
30 )
31 !
32
33 USE wx_pbl_var_mod
34
35 USE print_control_mod, ONLY: prt_level,lunout
36 USE indice_sol_mod, ONLY: is_oce
37 !
38 INCLUDE "YOMCST.h"
39 INCLUDE "FCTTRE.h"
40 INCLUDE "YOETHF.h"
41 INCLUDE "clesphys.h"
42 !
43 INTEGER, INTENT(IN) :: knon ! number of grid cells
44 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa)
45 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa)
46 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area
47 REAL, DIMENSION(knon), INTENT(IN) :: dTs_forcing ! forced temperature difference (w)-(x)
48 REAL, DIMENSION(knon), INTENT(IN) :: dqs_forcing ! forced humidity difference (w)-(x)
49 REAL, DIMENSION(knon,klev), INTENT(IN) :: yt_x, yt_w, yq_x, yq_w
50 REAL, DIMENSION(knon,klev), INTENT(IN) :: yu_x, yu_w, yv_x, yv_w
51 REAL, DIMENSION(knon), INTENT(IN) :: ycdragh_x, ycdragh_w, ycdragq_x, ycdragq_w
52 REAL, DIMENSION(knon), INTENT(IN) :: ycdragm_x, ycdragm_w
53 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT_x, AcoefT_w, AcoefQ_x, AcoefQ_w
54 REAL, DIMENSION(knon), INTENT(IN) :: AcoefU_x, AcoefU_w, AcoefV_x, AcoefV_w
55 REAL, DIMENSION(knon), INTENT(IN) :: BcoefT_x, BcoefT_w, BcoefQ_x, BcoefQ_w
56 REAL, DIMENSION(knon), INTENT(IN) :: BcoefU_x, BcoefU_w, BcoefV_x, BcoefV_w
57 REAL, DIMENSION(knon), INTENT(OUT) :: AcoefT, AcoefQ, AcoefU, AcoefV
58 REAL, DIMENSION(knon), INTENT(OUT) :: BcoefT, BcoefQ, BcoefU, BcoefV
59 REAL, DIMENSION(knon), INTENT(OUT) :: ycdragh, ycdragq, ycdragm
60 REAL, DIMENSION(knon), INTENT(OUT) :: yt1, yq1, yu1, yv1 ! Apparent T, q, u, v at first level, as
61 !seen by surface modules
62 !
63 ! Local variables
64 INTEGER :: j
65 REAL :: dd_Kh
66 REAL :: dd_Kq
67 REAL :: dd_Km
68 REAL :: dd_u
69 REAL :: dd_v
70 REAL :: dd_t
71 REAL :: dd_q
72 !
73 REAL :: LambdaTs, LambdaQs, LambdaUs, LambdaVs
74 !
75 REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region
76 !
77 !
78 sigx(1:knon) = 1.-sigw(1:knon)
79 !
80 !
81 DO j=1,knon
82 !
83 !
84 ! Compute w-x differences
85 dd_t = yt_w(j,1) - yt_x(j,1)
86 dd_q = yq_w(j,1) - yq_x(j,1)
87 dd_u = yu_w(j,1) - yu_x(j,1)
88 dd_v = yv_w(j,1) - yv_x(j,1)
89 !
90 ! Merged exchange coefficients
91 dd_Kh = Kech_h_w(j) - Kech_h_x(j)
92 dd_Kq = Kech_q_w(j) - Kech_q_x(j)
93 dd_Km = Kech_m_w(j) - Kech_m_x(j)
94 !
95 LambdaTs = dd_KTp(j)/Kech_Tp(j)
96 LambdaQs = dd_KQs(j)/Kech_Qs(j)
97 LambdaUs = dd_KUp(j)/Kech_Up(j)
98 LambdaVs = dd_KVp(j)/Kech_Vp(j)
99 !
100 ! Calcul des coef A, B \'equivalents dans la couche 1
101 !
102 ! The dTs_forcing and dqs_forcing terms are added for diagnostic purpose ; they should be zero in normal operation.
103 AcoefT(j) = AcoefT_x(j) + sigw(j)*(1.+sigx(j)*LambdaTs)*(dd_AT(j) - C_p(j)*dTs_forcing(j))
104 AcoefQ(j) = AcoefQ_x(j) + sigw(j)*(1.+sigx(j)*LambdaQs)*(dd_AQ(j) - dqs_forcing(j))
105 AcoefU(j) = AcoefU_x(j) + sigw(j)*(1.+sigx(j)*LambdaUs)*dd_AU(j)
106 AcoefV(j) = AcoefV_x(j) + sigw(j)*(1.+sigx(j)*LambdaVs)*dd_AV(j)
107 !
108 !
109 !! BcoefT(j) = (sigw(j)*Kech_h_w(j)*Kech_T_pw(j)*BcoefT_w(j) + &
110 !! sigx(j)*Kech_h_x(j)*Kech_T_px(j)*BcoefT_x(j) )/(Kech_h(j)*Kech_Tp(j))
111 !! BcoefQ(j) = (sigw(j)*Kech_q_w(j)*Kech_Q_pw(j)*BcoefQ_w(j) + &
112 !! sigx(j)*Kech_q_x(j)*Kech_Q_px(j)*BcoefQ_x(j) )/(Kech_q(j)*Kech_Qp(j))
113 !! BcoefU(j) = (sigw(j)*Kech_m_w(j)*Kech_U_pw(j)*BcoefU_w(j) + &
114 !! sigx(j)*Kech_m_x(j)*Kech_U_px(j)*BcoefU_x(j) )/(Kech_m(j)*Kech_Up(j))
115 !! BcoefV(j) = (sigw(j)*Kech_m_w(j)*Kech_V_pw(j)*BcoefV_w(j) + &
116 !! sigx(j)*Kech_m_x(j)*Kech_V_px(j)*BcoefV_x(j) )/(Kech_m(j)*Kech_Vp(j))
117 !
118 !! Print *,'YYYYpbl0: BcoefT_x, sigw, sigx, dd_Kh, dd_KTp, Kech_h_w ', &
119 !! BcoefT_x, sigw, sigx, dd_Kh, dd_KTp, Kech_h_w
120 !! Print *,'YYYYpbl0: Kech_T_pw, dd_BT, Kech_h, Kech_Tp ', &
121 !! Kech_T_pw, dd_BT, Kech_h, Kech_Tp
122 BcoefT(j) = BcoefT_x(j) + sigw(j)*(sigx(j)*dd_Kh*dd_KTp(j)*BcoefT_x(j) + &
123 Kech_h_w(j)*Kech_T_pw(j)*dd_BT(j))/(Kech_h(j)*Kech_Tp(j))
124 BcoefQ(j) = BcoefQ_x(j) + sigw(j)*(sigx(j)*dd_Kq*dd_KQs(j)*BcoefQ_x(j) + &
125 Kech_q_w(j)*Kech_Q_sw(j)*dd_BQ(j))/(Kech_q(j)*Kech_Qs(j))
126 BcoefU(j) = BcoefU_x(j) + sigw(j)*(sigx(j)*dd_Km*dd_KUp(j)*BcoefU_x(j) + &
127 Kech_m_w(j)*Kech_U_pw(j)*dd_BU(j))/(Kech_m(j)*Kech_Up(j))
128 BcoefV(j) = BcoefV_x(j) + sigw(j)*(sigx(j)*dd_Km*dd_KVp(j)*BcoefV_x(j) + &
129 Kech_m_w(j)*Kech_V_pw(j)*dd_BV(j))/(Kech_m(j)*Kech_Vp(j))
130 !>jyg
131 !
132 !
133 ! Calcul des cdrag \'equivalents dans la couche
134 !
135 ycdragm(j) = ycdragm_x(j) + sigw(j)*dd_Cdragm(j)
136 ycdragh(j) = ycdragh_x(j) + sigw(j)*dd_Cdragh(j)
137 ycdragq(j) = ycdragq_x(j) + sigw(j)*dd_Cdragq(j)
138 !
139 ! Calcul de T, q, u et v \'equivalents dans la couche 1
140 !! yt1(j) = yt_x(j,1) + sigw(j)*dd_t*(1.+sigx(j)*dd_Kh/KCT)
141 !! yq1(j) = yq_x(j,1) + sigw(j)*dd_q*(1.+sigx(j)*dd_Kh/KCQ)
142 !! yu1(j) = yu_x(j,1) + sigw(j)*dd_u*(1.+sigx(j)*dd_Km/KCU)
143 !! yv1(j) = yv_x(j,1) + sigw(j)*dd_v*(1.+sigx(j)*dd_Km/KCV)
144 yt1(j) = yt_x(j,1) + sigw(j)*dd_t
145 yq1(j) = yq_x(j,1) + sigw(j)*dd_q
146 yu1(j) = yu_x(j,1) + sigw(j)*dd_u
147 yv1(j) = yv_x(j,1) + sigw(j)*dd_v
148
149
150 ENDDO
151
152 RETURN
153
154 END SUBROUTINE wx_pbl0_merge
155
156 SUBROUTINE wx_pbl_dts_merge(knon, dtime, ypplay, ypaprs, &
157 sigw, beta, wcstar, wdens, &
158 AT_x, AT_w, &
159 BT_x, BT_w, &
160 AcoefT0, AcoefQ0, BcoefT0, BcoefQ0, &
161 AcoefT, AcoefQ, BcoefT, BcoefQ, &
162 HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, &
163 phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0, &
164 g_T, g_Q, &
165 Gamma_phiT, Gamma_phiQ, &
166 dTs_ins, dqsatsrf_ins &
167 )
168 !
169
170 USE wx_pbl_var_mod
171
172 USE print_control_mod, ONLY: prt_level,lunout
173 !
174 INCLUDE "YOMCST.h"
175 INCLUDE "FCTTRE.h"
176 INCLUDE "YOETHF.h"
177 !
178 INTEGER, INTENT(IN) :: knon ! number of grid cells
179 REAL, INTENT(IN) :: dtime ! time step size (s)
180 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa)
181 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa)
182 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pool fractional area
183 REAL, DIMENSION(knon), INTENT(IN) :: beta ! evaporation by potential evaporation
184 REAL, DIMENSION(knon), INTENT(IN) :: wcstar ! cold pool gust front speed
185 REAL, DIMENSION(knon), INTENT(IN) :: wdens ! cold pool number density
186 REAL, DIMENSION(knon), INTENT(IN) :: AT_x, AT_w
187 REAL, DIMENSION(knon), INTENT(IN) :: BT_x, BT_w
188 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT0, AcoefQ0, BcoefT0, BcoefQ0
189 !
190 REAL, DIMENSION(knon), INTENT(OUT) :: AcoefT, AcoefQ, BcoefT, BcoefQ
191 REAL, DIMENSION(knon), INTENT(OUT) :: HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn
192 REAL, DIMENSION(knon), INTENT(OUT) :: phiT0_b, dphiT0, phiQ0_b, dphiQ0, Rn0_b, dRn0
193 REAL, DIMENSION(knon), INTENT(OUT) :: g_T, g_Q
194 REAL, DIMENSION(knon), INTENT(OUT) :: Gamma_phiT, Gamma_phiQ
195 REAL, DIMENSION(knon), INTENT(OUT) :: dTs_ins, dqsatsrf_ins
196 !
197 ! Local variables
198 REAL, DIMENSION(knon) :: qsat_x
199 REAL, DIMENSION(knon) :: qsat_w
200 REAL, DIMENSION(knon) :: dqsatdT_x
201 REAL, DIMENSION(knon) :: dqsatdT_w
202 !
203 REAL, DIMENSION(knon) :: T10_x
204 REAL, DIMENSION(knon) :: T10_w
205 REAL, DIMENSION(knon) :: phiT0_x
206 REAL, DIMENSION(knon) :: phiT0_w
207 REAL, DIMENSION(knon) :: phiQ0_x
208 REAL, DIMENSION(knon) :: phiQ0_w
209 REAL, DIMENSION(knon) :: Rn0_x
210 REAL, DIMENSION(knon) :: Rn0_w
211 REAL, DIMENSION(knon) :: Rp1_x
212 REAL, DIMENSION(knon) :: Rp1_w
213 REAL, DIMENSION(knon) :: Rps_x
214 REAL, DIMENSION(knon) :: Rps_w
215 !
216 REAL, DIMENSION(knon) :: HTphiT_x
217 REAL, DIMENSION(knon) :: HTphiT_w
218 REAL, DIMENSION(knon) :: HTphiQ_x
219 REAL, DIMENSION(knon) :: HTphiQ_w
220 REAL, DIMENSION(knon) :: HTRn_x
221 REAL, DIMENSION(knon) :: HTRn_w
222 !
223 REAL, DIMENSION(knon) :: HQphiT_x
224 REAL, DIMENSION(knon) :: HQphiT_w
225 REAL, DIMENSION(knon) :: HQphiQ_x
226 REAL, DIMENSION(knon) :: HQphiQ_w
227 REAL, DIMENSION(knon) :: HQRn_x
228 REAL, DIMENSION(knon) :: HQRn_w
229 !
230 REAL, DIMENSION(knon) :: HQphiT_b
231 REAL, DIMENSION(knon) :: dd_HQphiT
232 REAL, DIMENSION(knon) :: HQphiQ_b
233 REAL, DIMENSION(knon) :: dd_HQphiQ
234 REAL, DIMENSION(knon) :: HQRn_b
235 REAL, DIMENSION(knon) :: dd_HQRn
236 !
237
238 REAL, DIMENSION(knon) :: sigx
239 !
240 REAL, DIMENSION(knon) :: Ts, T1
241 !!! REAL, DIMENSION(knon) :: qsat, dqsat_dT
242 !!! REAL, DIMENSION(knon) :: phiT0
243 !
244 !!! REAL, DIMENSION(knon) :: Cp, Lv
245 REAL, DIMENSION(knon) :: tau, Inert
246 !
247 REAL :: dd_Kh
248 REAL :: zdelta, zcvm5, zcor
249 REAL :: qsat
250 !
251 INTEGER :: j
252
253
254 !----------------------------------------------------------------------------
255 ! Reference state
256 ! ---------------
257 ! dqsat_dT_w = dqsat_dT(Ts0_w) dqsat_dT_x = dqsat_dT(Ts0_x)
258 ! T10_w = (AT_w/Cp - Kech_T_w BT_w dtime Ts0_w)/(1 - Kech_T_w BT_w dtime)
259 ! T10_x = (AT_x/Cp - Kech_T_x BT_x dtime Ts0_x)/(1 - Kech_T_x BT_x dtime)
260 ! phiT0_w = Kech_T_pw (AT_w - Cp Ts0_w) phiT0_x = Kech_T_px (AT_x - Cp Ts0_x)
261 ! phiQ0_w = Kech_Q_sw (beta AQ_w - qsatsrf0_w) phiQ0_x = Kech_Q_sx (beta AQ_x - qsatsrf0_x)
262 ! Rn0_w = eps_1 Rsigma T10_w^4 - Rsigma Ts0_w^4 Rn0_x = eps_1 Rsigma T10_x^4 - Rsigma Ts0_x^4
263 ! Rp1_w = 4 eps_1 Rsigma T10_w^3 Rp1_x = 4 eps_1 Rsigma T10_x^3
264 ! Rps_w = 4 Rsigma Ts0_w^3 Rps_x = 4 Rsigma Ts0_x^3
265 !
266 ! phiT0_b = sigw phiT0_w + sigx phiT0_x
267 ! dphiT0 = phiT0_w - phiT0_x
268 ! phiQ0_b = sigw phiQ0_w + sigx phiQ0_x
269 ! dphiQ0 = phiQ0_w - phiQ0_x
270 ! Rn0_b = sigw Rn0_w + sigx Rn0_x
271 dRn0 = Rn0_w - Rn0_x
272 !
273 !
274 !----------------------------------------------------------------------------
275 ! Elementary enthalpy equations
276 ! -----------------------------
277 ! phiT_w = phiT0_w - HTphiT_w (Ts_w-Ts0_w) phiT_x = phiT0_x - HTphiT_x (Ts_x-Ts0_x)
278 ! phiQ_w = phiQ0_w - HTphiQ_w (Ts_w-Ts0_w) phiQ_x = phiQ0_x - HTphiQ_x (Ts_x-Ts0_x)
279 ! Rn_w = Rn0_w - HTRn_w (Ts_w-Ts0_w) Rn_x = Rn0_x - HTRn_x (Ts_x-Ts0_x)
280 ! DFlux_DT coefficients
281 ! ---------------------
282 ! Heat flux equation
283 ! HTphiT_w = Cp Kech_T_pw HTphiT_x = Cp Kech_T_px
284 ! Moisture flux equation
285 ! HTphiQ_w = beta Kech_Q_sw dqsat_dT_w HTphiQ_x = beta Kech_Q_sx dqsat_dT_x
286 ! Radiation equation
287 ! HTRn_w = Rp1_w Kech_T_pw BcoefT_w dtime + Rps_w HTRn_x = Rp1_x Kech_T_px BcoefT_x dtime + Rps_x
288 !
289 !----------------------------------------------------------------------------
290 ! Elementary moisture equations
291 ! -----------------------------
292 ! beta Ts_w = beta Ts0_w + QQ_w (qsatsrf_w-qsatsrf0_w) beta Ts_x = beta Ts0_x + QQ_x (qsatsrf_x-qsatsrf0_x)
293 ! beta phiT_w = beta phiT0_w - HQphiT_w (qsatsrf_w-qsatsrf0_w) beta phiQ_x = beta phiQ0_x - HTphiQ_x (qsatsrf_x-qsatsrf0_x)
294 ! beta phiQ_w = beta phiQ0_w - HQphiQ_w (qsatsrf_w-qsatsrf0_w) beta phiQ_x = beta phiQ0_x - HTphiQ_x (qsatsrf_x-qsatsrf0_x)
295 ! beta Rn_w = beta Rn0_w - HQRn_w (qsatsrf_w-qsatsrf0_w) beta Rn_x = beta Rn0_x - HTRn_x (qsatsrf_x-qsatsrf0_x)
296 ! DFluxDQ coefficients
297 ! ---------------------
298 ! dqsat_dT equation
299 ! QQ_w = 1. / dqsat_dT_w QQ_x = 1. / dqsat_dT_x
300 ! Heat flux equation
301 ! HQphiT_w = Cp Kech_T_pw QQ_w HQphiT_x = Cp Kech_T_px QQ_x
302 ! Moisture flux equation
303 ! HQphiQ_w = beta Kech_Q_sw HQphiQ_x = beta Kech_Q_sx
304 ! Radiation equation
305 ! HQRn_w = (Rp1_w Kech_T_pw BcoefT_w dtime + Rps_w) QQ_w
306 ! HQRn_x = (Rp1_x Kech_T_px BcoefT_x dtime + Rps_x) QQ_x
307 !
308 !----------------------------------------------------------------------------
309 ! Mean values and w-x differences
310 ! -------------------------------
311 ! HTphiT_b = sigw HTphiT_w + sigx HTphiT_x dd_HTphiT = HTphiT_w - HTphiT_x
312 ! HTphiQ_b = sigw HTphiQ_w + sigx HTphiQ_x dd_HTphiQ = HTphiQ_w - HTphiQ_x
313 ! HTRn_b = sigw HTRn_w + sigx HTRn_x dd_HTRn = HTRn_w - HTRn_x
314 !
315 ! QQ_b = sigw QQ_w + sigx QQ_x dd_QQ = QQ_w - QQ_x
316 ! HQphiT_b = sigw HQphiT_w + sigx HQphiT_x dd_HQphiT = HQphiT_w - HQphiT_x
317 ! HQphiQ_b = sigw HQphiQ_w + sigx HQphiQ_x dd_HQphiQ = HQphiQ_w - HQphiQ_x
318 ! HQRn_b = sigw HQRn_w + sigx HQRn_x dd_HQRn = HQRn_w - HQRn_x
319 !
320 !----------------------------------------------------------------------------
321 ! Equations
322 ! ---------
323 ! (1 - g_T) dTs = dTs_ins + Gamma_phiT phiT
324 ! (1 - g_Q) dqsatsrf = dqsatsrf_ins + Gamma_phiQ phiQ
325 !
326 ! Feedback Gains
327 ! --------------
328 ! g_T = - (sqrt(tau)/I) [ HTphiT_b + Lv HTphiQ_b + HTRn_b + &
329 ! (dd_HTphiT + Lv dd_HTphiQ + dd_HTRn) (sigx - sigw - sigw sigx dd_HTphiT/HTphiT_b) ]
330 ! g_Q = - (sqrt(tau)/(I QQ_b)) ( HQphiT_b + Lv HQphiQ_b + HQRn_b ) - &
331 ! (sigx - sigw - sigw sigx dd_HQphiQ/HQphiQ_b) &
332 ! [ dd_QQ/QQ_b + (sqrt(tau)/(I QQ_b))(dd_HQphiT + Lv dd_HQphiQ + dd_HQRn) ]
333 !
334 ! Ts, qs Coupling coefficients /
335 ! ----------------------------
336 ! Gamma_phiT = (sqrt(tau)/(I HTphiT_b)) (dd_HTphiT + Lv dd_HTphiQ + dd_HTRn)
337 ! Gamma_phiQ = (1/(HQphiQ_b QQ_b)) [ dd_QQ + (sqrt(tau)/(I )) (dd_HQphiT + Lv dd_HQphiQ + dd_HQRn) ]
338 !
339 ! Insensitive changes
340 ! -------------------
341 ! dTs_ins = (1 - g_T) dTs0 - Gamma_phiT phiT0_b
342 ! dqsatsrf_ins = (1 - g_Q) dqsatsrf0 - Gamma_phiQ phiQ0_b
343 !
344 !----------------------------------------------------------------------------
345 ! Effective coefficients Acoef and Bcoef
346 ! --------------------------------------
347 ! Equations
348 ! ---------
349 ! Cp Ta = AcoefT + BcoefT phiT dtime
350 ! qa = AcoefQ + BcoefQ phiQ dtime
351 ! Coefficients
352 ! ------------
353 ! AcoefT = AcoefT0 - sigw sigx (dd_KTp/Kech_Tp) Cp dTs_ins/(1 - g_T)
354 ! BcoefT = BcoefT0 - sigw sigx (dd_KTp/Kech_Tp) Cp Gamma_phiT/(1 - g_T)/dtime
355 !
356 ! AcoefQ = AcoefQ0 - sigw sigx (dd_KQp/Kech_Qp) dqs_ins/(1 - g_Q)
357 ! BcoefQ = BcoefQ0 - sigw sigx (dd_KQp/Kech_Qp) Gamma_phiq/(1 - g_Q)/dtime
358 !
359 !==============================================================================
360 !
361 !
362 ! Parameters
363 ! ----------
364 Inert(1:knon) = 2000.
365 tau(1:knon) = sqrt(sigw(1:knon)/max(rpi*wdens(1:knon)*wcstar(1:knon)**2 , &
366 sigw(1:knon)*1.e-12,smallestreal))
367 sigx(1:knon) = 1.-sigw(1:knon)
368 !! Compute Cp, Lv, qsat, dqsat_dT.
369 ! C_p(1:knon) = RCpd
370 ! L_v(1:knon) = RLvtt
371 !
372 ! print *,' AAAA wx_pbl_dTs, C_p(j), qsat0(j), Ts0(j) : ', C_p(:), qsat0(:), Ts0(:)
373 !
374 !
375 T10_x(1:knon) = (AT_x(1:knon)/C_p(1:knon) - Kech_h_x(1:knon)*BT_x(1:knon)*dtime*Ts0_x(1:knon))/ &
376 (1 - Kech_h_x(1:knon)*BT_x(1:knon)*dtime)
377 T10_w(1:knon) = (AT_w(1:knon)/C_p(1:knon) - Kech_h_w(1:knon)*BT_w(1:knon)*dtime*Ts0_w(1:knon))/ &
378 (1 - Kech_h_w(1:knon)*BT_w(1:knon)*dtime)
379 !
380 phiT0_x(1:knon) = Kech_T_px(1:knon)*(AT_x(1:knon) - C_p(1:knon)*Ts0_x(1:knon))
381 phiT0_w(1:knon) = Kech_T_pw(1:knon)*(AT_w(1:knon) - C_p(1:knon)*Ts0_w(1:knon))
382 !
383 phiQ0_x(1:knon) = Kech_Q_sx(1:knon)*(beta(1:knon)*AQ_x(1:knon) - qsatsrf0_x(1:knon))
384 phiQ0_w(1:knon) = Kech_Q_sw(1:knon)*(beta(1:knon)*AQ_w(1:knon) - qsatsrf0_w(1:knon))
385 !
386 Rn0_x(1:knon) = eps_1*Rsigma*T10_x(1:knon)**4 - Rsigma*Ts0_x(1:knon)**4
387 Rn0_w(1:knon) = eps_1*Rsigma*T10_w(1:knon)**4 - Rsigma*Ts0_w(1:knon)**4
388 !
389 Rp1_x(1:knon) = 4*eps_1*Rsigma*T10_x(1:knon)**3
390 Rp1_w(1:knon) = 4*eps_1*Rsigma*T10_w(1:knon)**3
391 !
392 Rps_x(1:knon) = 4*Rsigma*Ts0_x(1:knon)**3
393 Rps_w(1:knon) = 4*Rsigma*Ts0_w(1:knon)**3
394 !
395 ! DFlux_DT coefficients
396 ! ---------------------
397 ! Heat flux equation
398 HTphiT_x(1:knon) = C_p(1:knon)*Kech_T_px(1:knon)
399 HTphiT_w(1:knon) = C_p(1:knon)*Kech_T_pw(1:knon)
400 ! Moisture flux equation
401 HTphiQ_x(1:knon) = beta(1:knon)*Kech_Q_sx(1:knon)*dqsatdT0_x(1:knon)
402 HTphiQ_w(1:knon) = beta(1:knon)*Kech_Q_sw(1:knon)*dqsatdT0_w(1:knon)
403 ! Radiation equation
404 HTRn_x(1:knon) = Rp1_x(1:knon)*Kech_T_px(1:knon)*BT_x(1:knon)*dtime + Rps_x(1:knon)
405 HTRn_w(1:knon) = Rp1_w(1:knon)*Kech_T_pw(1:knon)*BT_w(1:knon)*dtime + Rps_w(1:knon)
406 !
407 ! DFluxDQ coefficients
408 ! ---------------------
409 ! Heat flux equation
410 HQphiT_x(1:knon) = C_p(1:knon)*Kech_T_px(1:knon)*QQ_x(1:knon)
411 HQphiT_w(1:knon) = C_p(1:knon)*Kech_T_pw(1:knon)*QQ_w(1:knon)
412 ! Moisture flux equation
413 HQphiQ_x(1:knon) = beta(1:knon)*Kech_Q_sx(1:knon)
414 HQphiQ_w(1:knon) = beta(1:knon)*Kech_Q_sw(1:knon)
415 ! Radiation equation
416 HQRn_x(1:knon) = (Rp1_x(1:knon)*Kech_T_px(1:knon)*BT_x(1:knon)*dtime + Rps_x(1:knon))*QQ_x(1:knon)
417 HQRn_w(1:knon) = (Rp1_w(1:knon)*Kech_T_pw(1:knon)*BT_w(1:knon)*dtime + Rps_w(1:knon))*QQ_w(1:knon)
418 !
419 ! Mean values and w-x differences
420 ! -------------------------------
421 phiT0_b(1:knon) = sigw(1:knon)*phiT0_w(1:knon) + sigx(1:knon)*phiT0_x(1:knon)
422 phiQ0_b(1:knon) = sigw(1:knon)*phiQ0_w(1:knon) + sigx(1:knon)*phiQ0_x(1:knon)
423 Rn0_b(1:knon) = sigw(1:knon)*Rn0_w(1:knon) + sigx(1:knon)*Rn0_x(1:knon)
424 !
425 dphiT0(1:knon) = phiT0_w(1:knon) - phiT0_x(1:knon)
426 dphiQ0(1:knon) = phiQ0_w(1:knon) - phiQ0_x(1:knon)
427 dRn0(1:knon) = Rn0_w(1:knon) - Rn0_x(1:knon)
428 !
429 HTphiT_b(1:knon) = sigw(1:knon)*HTphiT_w(1:knon) + sigx(1:knon)*HTphiT_x(1:knon)
430 dd_HTphiT(1:knon) = HTphiT_w(1:knon) - HTphiT_x(1:knon)
431 !
432 HTphiQ_b(1:knon) = sigw(1:knon)*HTphiQ_w(1:knon) + sigx(1:knon)*HTphiQ_x(1:knon)
433 dd_HTphiQ(1:knon) = HTphiQ_w(1:knon) - HTphiQ_x(1:knon)
434 !
435 HTRn_b(1:knon) = sigw(1:knon)*HTRn_w(1:knon) + sigx(1:knon)*HTRn_x(1:knon)
436 dd_HTRn(1:knon) = HTRn_w(1:knon) - HTRn_x(1:knon)
437 !
438 HQphiT_b(1:knon) = sigw(1:knon)*HQphiT_w(1:knon) + sigx(1:knon)*HQphiT_x(1:knon)
439 dd_HQphiT(1:knon) = HQphiT_w(1:knon) - HQphiT_x(1:knon)
440 !
441 HQphiQ_b(1:knon) = sigw(1:knon)*HQphiQ_w(1:knon) + sigx(1:knon)*HQphiQ_x(1:knon)
442 dd_HQphiQ(1:knon) = HQphiQ_w - HQphiQ_x(1:knon)
443 !
444 HQRn_b(1:knon) = sigw(1:knon)*HQRn_w(1:knon) + sigx(1:knon)*HQRn_x(1:knon)
445 dd_HQRn(1:knon) = HQRn_w(1:knon) - HQRn_x(1:knon)
446 !
447 ! Feedback Gains
448 ! --------------
449 g_T(1:knon) = - (sqrt(tau(1:knon))/Inert(1:knon)) &
450 * (HTphiT_b(1:knon) + L_v(1:knon)*HTphiQ_b(1:knon) + HTRn_b(1:knon) &
451 + (dd_HTphiT(1:knon) + L_v(1:knon)*dd_HTphiQ(1:knon) + dd_HTRn(1:knon)) &
452 * (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_HTphiT(1:knon)/HTphiT_b(1:knon)) )
453 !
454 !!!! DO j = 1,knon
455 !!!! IF (mod(j,20) .eq.0) THEN
456 !!!! print *, ' j dd_QQ QQ_b dd_HQphiQ dd_HQphiT dd_HQRn HQphiQ_b HQphiT_b HQRn_b '
457 !!!! ENDIF
458 !!!! print 1789, j, dd_QQ(j), QQ_b(j), dd_HQphiQ(j), dd_HQphiT(j), dd_HQRn(j), HQphiQ_b(j), HQphiT_b(j), HQRn_b(j)
459 !!!! 1789 FORMAT( I4, 10(1X,E10.2))
460 !!!! ENDDO
461 g_Q(1:knon) = - (dd_QQ(1:knon)/QQ_b(1:knon)) * &
462 (sigx(1:knon)-sigw(1:knon)-sigw(1:knon)*sigx(1:knon)*dd_KQs(1:knon)/Kech_Qs(1:knon)) &
463 - sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon)) * &
464 ( HQphiT_b(1:knon) + L_v(1:knon)*HQphiQ_b(1:knon) + HQRn_b(1:knon) + &
465 (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_KQs(1:knon)/Kech_Qs(1:knon)) * &
466 (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) )
467
468 !! g_Q(1:knon) = - (dd_QQ(1:knon)/QQ_b(1:knon)) * &
469 !! (sigx(1:knon)-sigw(1:knon)-sigw(1:knon)*sigx(1:knon)*dd_HQphiQ(1:knon)/HQphiQ_b(1:knon)) &
470 !! - sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon)) * &
471 !! ( HQphiT_b(1:knon) + L_v(1:knon)*HQphiQ_b(1:knon) + HQRn_b(1:knon) + &
472 !! (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_HQphiQ(1:knon)/HQphiQ_b(1:knon)) * &
473 !! (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) )
474
475 !! g_Q(1:knon) = - (sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon))) * &
476 !! ( HQphiT_b(1:knon) + L_v(1:knon)*HQphiQ_b(1:knon) + HQRn_b(1:knon) ) &
477 !! - (sigx(1:knon) - sigw(1:knon) - sigw(1:knon)*sigx(1:knon)*dd_HQphiQ(1:knon)/HQphiQ_b(1:knon)) * &
478 !! ( dd_QQ(1:knon)/QQ_b(1:knon) &
479 !! + (sqrt(tau(1:knon))/(Inert(1:knon)*QQ_b(1:knon))) &
480 !! * (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) )
481
482 ! Ts, qs Coupling coefficients /
483 ! ----------------------------
484 Gamma_phiT(1:knon) = (sqrt(tau(1:knon))/(Inert(1:knon)*HTphiT_b(1:knon))) &
485 * (dd_HTphiT(1:knon) + L_v(1:knon)*dd_HTphiQ(1:knon) + dd_HTRn(1:knon))
486 !
487 Gamma_phiQ(1:knon) = (1./(Kech_Qs(1:knon)*QQ_b(1:knon))) * &
488 ( dd_QQ(1:knon) &
489 + (sqrt(tau(1:knon))/(Inert(1:knon))) * &
490 (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) )
491
492 !! Gamma_phiQ(1:knon) = (beta(1:knon)/(HQphiQ_b(1:knon)*QQ_b(1:knon))) * &
493 !! ( dd_QQ(1:knon) &
494 !! + (sqrt(tau(1:knon))/(Inert(1:knon))) * &
495 !! (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) )
496
497 !! Gamma_phiQ(1:knon) = (1/(HQphiQ_b(1:knon)*QQ_b(1:knon))) &
498 !! * ( dd_QQ(1:knon) &
499 !! + (sqrt(tau(1:knon))/(Inert(1:knon))) &
500 !! * (dd_HQphiT(1:knon) + L_v(1:knon)*dd_HQphiQ(1:knon) + dd_HQRn(1:knon)) )
501 !
502 ! Insensitive changes
503 ! -------------------
504 dTs_ins(1:knon) = (sqrt(tau(1:knon))/Inert(1:knon))* &
505 (dphiT0(1:knon) + L_v(1:knon)*dphiQ0(1:knon) + dRn0(1:knon))
506 !
507 dqsatsrf_ins(1:knon) = (beta(1:knon)/QQ_b(1:knon))*dTs_ins(1:knon)
508 !
509 IF (prt_level .Ge. 10) THEN
510 print *,'wx_pbl_merge, tau ', tau
511 print *,'wx_pbl_merge, AcoefT0 ', AcoefT0
512 print *,'wx_pbl_merge, AcoefQ0 ', AcoefQ0
513 print *,'wx_pbl_merge, BcoefT0 ', BcoefT0
514 print *,'wx_pbl_merge, BcoefQ0 ', BcoefQ0
515 print *,'wx_pbl_merge, qsat0_w, qsat0_x ', (qsat0_w(j), qsat0_x(j),j=1,knon)
516 print *,'wx_pbl_merge, dqsatdT0_w, dqsatdT0_x ', (dqsatdT0_w(j), dqsatdT0_x(j),j=1,knon)
517 ENDIF
518 !
519 !----------------------------------------------------------------------------
520 !
521 !------------------------------------------------------------------------------
522 !
523 ! Effective coefficients Acoef and Bcoef
524 ! --------------------------------------
525 DO j = 1,knon
526 AcoefT(j) = AcoefT0(j) - sigw(j)*sigx(j)*(dd_KTp(j)/Kech_Tp(j))*C_p(j)* &
527 (dTs0(j) + (dTs_ins(j)-dTs0(j)-Gamma_phiT(j)*phiT0_b(j))/(1. - g_T(j)))
528 BcoefT(j) = BcoefT0(j) - sigw(j)*sigx(j)*(dd_KTp(j)/Kech_Tp(j))*C_p(j)*Gamma_phiT(j)/(1. - g_T(j))/dtime
529
530 AcoefQ(j) = AcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))* &
531 (dqsatsrf0(j) + (dqsatsrf_ins(j)-(beta(j)/QQ_b(j))*dTs0(j)-Gamma_phiQ(j)*phiQ0_b(j))/(1 - g_Q(j)))/ &
532 max(beta(j),1.e-4)
533 BcoefQ(j) = BcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))*Gamma_phiQ(j)/(1 - g_Q(j))/ &
534 (max(beta(j),1.e-4)*dtime)
535 !! AcoefQ(j) = AcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))* &
536 !! (dqsatsrf0(j) + (dqsatsrf_ins(j)-(beta(j)/QQ_b(j))*dTs0(j)-Gamma_phiQ(j)*phiQ0_b(j))/(1 - g_Q(j)))/ &
537 !! beta(j)
538 !! BcoefQ(j) = BcoefQ0(j) - sigw(j)*sigx(j)*(dd_KQs(j)/Kech_Qs(j))*Gamma_phiQ(j)/(1 - g_Q(j))/(beta(j)*dtime)
539 ENDDO ! j = 1,knon
540
541 IF (prt_level .Ge. 10) THEN
542 print *,'wx_pbl_dts AAAA BcoefQ, BcoefQ0, sigw ', &
543 BcoefQ, BcoefQ0, sigw
544 print *,'wx_pbl_dts_merge, dTs_ins ', dTs_ins
545 print *,'wx_pbl_dts_merge, dqs_ins ', dqsatsrf_ins
546 ENDIF
547
548 RETURN
549
550 END SUBROUTINE wx_pbl_dts_merge
551
552 SUBROUTINE wx_pbl_split(knon, nsrf, dtime, sigw, beta, iflag_split, &
553 g_T, g_Q, &
554 Gamma_phiT, Gamma_phiQ, &
555 dTs_ins, dqsatsrf_ins, &
556 phiT, phiQ, phiU, phiV, &
557 !!!! HTRn_b, dd_HTRn, HTphiT_b, dd_HTphiT, &
558 phiQ0_b, phiT0_b, &
559 phiT_x, phiT_w, &
560 phiQ_x, phiQ_w, &
561 phiU_x, phiU_w, &
562 phiV_x, phiV_w, &
563 philat_x, philat_w, &
564 !!!! Rn_b, dRn, &
565 dqsatsrf, &
566 dTs, delta_qsurf &
567 )
568 !
569
570 USE wx_pbl_var_mod
571
572 USE print_control_mod, ONLY: prt_level,lunout
573 USE indice_sol_mod, ONLY: is_oce
574 !
575 INCLUDE "YOMCST.h"
576 !
577 INTEGER, INTENT(IN) :: knon ! number of grid cells
578 INTEGER, INTENT(IN) :: nsrf ! surface type
579 REAL, INTENT(IN) :: dtime ! time step size (s)
580 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area
581 REAL, DIMENSION(knon), INTENT(IN) :: beta ! aridity factor
582 INTEGER, INTENT(IN) :: iflag_split
583 REAL, DIMENSION(knon), INTENT(IN) :: g_T, g_Q
584 REAL, DIMENSION(knon), INTENT(IN) :: Gamma_phiT, Gamma_phiQ
585 REAL, DIMENSION(knon), INTENT(IN) :: dTs_ins, dqsatsrf_ins
586 REAL, DIMENSION(knon), INTENT(IN) :: phiT, phiQ, phiU, phiV
587 REAL, DIMENSION(knon), INTENT(IN) :: phiQ0_b, phiT0_b
588 !
589 REAL, DIMENSION(knon), INTENT(OUT) :: phiT_x, phiT_w
590 REAL, DIMENSION(knon), INTENT(OUT) :: phiQ_x, phiQ_w
591 REAL, DIMENSION(knon), INTENT(OUT) :: phiU_x, phiU_w
592 REAL, DIMENSION(knon), INTENT(OUT) :: phiV_x, phiV_w
593 REAL, DIMENSION(knon), INTENT(OUT) :: philat_x, philat_w
594 REAL, DIMENSION(knon), INTENT(OUT) :: dqsatsrf ! beta delta(qsat(Ts))
595 REAL, DIMENSION(knon), INTENT(OUT) :: dTs ! Temperature difference at surface
596 REAL, DIMENSION(knon), INTENT(OUT) :: delta_qsurf
597 !
598 !! Local variables
599 INTEGER :: j
600 REAL, DIMENSION(knon) :: dphiT, dphiQ, dphiU, dphiV
601 REAL, DIMENSION(knon) :: q1_x, q1_w
602 !
603 REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region
604
605 !----------------------------------------------------------------------------
606 ! Equations
607 ! ---------
608 !!!!!! (1 - g_T) dTs = dTs_ins + Gamma_phiT phiT
609 !!!!!! (1 - g_Q) dqsatsrf = dqsatsrf_ins + Gamma_phiQ phiQ
610 !!!!!! dphiT = (dd_KTp/KTp) phiT + ( dd_AT - C_p dTs)*KxKwTp/KTp
611 !!!!!! dphiQ = (dd_KQs/KQs) phiQ + (beta dd_AQ - dqsatsrf )*KxKwQs/KQs
612 !!!!!! dphiU = (dd_KUp/KUp) phiU + ( dd_AU )*KxKwUp/KUp
613 !!!!!! dphiV = (dd_KVp/KVp) phiV + ( dd_AV )*KxKwVp/KVp
614 !
615 ! (1 - g_T) (dTs-dTs0) = dTs_ins-dTs0 + Gamma_phiT (phiT-phiT0)
616 ! (1 - g_Q) dqsatsrf = dqsatsrf_ins + Gamma_phiQ phiQ
617 ! dphiT = (dd_KTp/KTp) phiT + ( dd_AT - C_p dTs)*KxKwTp/KTp
618 ! dphiQ = (dd_KQs/KQs) phiQ + (beta dd_AQ - dqsatsrf )*KxKwQs/KQs
619 ! dphiU = (dd_KUp/KUp) phiU + ( dd_AU )*KxKwUp/KUp
620 ! dphiV = (dd_KVp/KVp) phiV + ( dd_AV )*KxKwVp/KVp
621 !
622 !!
623 sigx(:) = 1.-sigw(:)
624 !
625 ! print *,' AAAA wx_pbl_split, C_p(j), qsat0(j), Ts0(j) : ', C_p(:), qsat0(:), Ts0(:)
626 !
627 IF (iflag_split .EQ. 2 .AND. nsrf .NE. is_oce) THEN
628 !
629 ! Delta_tsurf and Delta_qsurf computation
630 ! -----------------------------------------
631 IF (prt_level >=10 ) THEN
632 print *,' wx_pbl_split, dTs_ins, dTs0 , Gamma_phiT, g_T ', dTs_ins, dTs0, Gamma_phiT, g_T
633 print *,' wx_pbl_split, dqsatsrf_ins, Gamma_phiQ, g_q ', dqsatsrf_ins, Gamma_phiQ, g_q
634 ENDIF
635 !
636 DO j = 1,knon
637 dTs(j) = dTs0(j) + (dTs_ins(j) - dTs0(j) + Gamma_phiT(j)*(phiT(j)-phiT0_b(j)) )/(1 - g_T(j))
638 dqsatsrf(j) = dqsatsrf0(j) + (dqsatsrf_ins(j) - (beta(j)/QQ_b(j))*dTs0(j) + &
639 Gamma_phiQ(j)*(phiQ(j)-phiQ0_b(j)) )/(1 - g_Q(j))
640 ENDDO ! j = 1,knon
641 !
642 IF (prt_level >=10 ) THEN
643 print *,' wx_pbl_split, dqsatsrf0, QQ_b ', dqsatsrf0, QQ_b
644 print *,' wx_pbl_split, phiT0_b, phiT, dTs ', phiT0_b, phiT, dTs
645 print *,' wx_pbl_split, phiQ0_b, phiQ, dqsatsrf ', phiQ0_b, phiQ, dqsatsrf
646 ENDIF
647 ELSE
648 dTs(:) = 0.
649 dqsatsrf(:) = 0.
650 ENDIF ! (iflag_split .EQ. 2 .AND. nsrf .NE. is_oce)
651 !
652 DO j = 1,knon
653 dphiT(j) = (phiT(j)*dd_KTp(j) + ( dd_AT(j) - C_p(j)*dTs(j))*KxKwTp(j))/Kech_Tp(j)
654 dphiQ(j) = (phiQ(j)*dd_KQs(j) + (beta(j)*dd_AQ(j) - dqsatsrf(j))*KxKwQs(j))/Kech_Qs(j)
655 dphiU(j) = (phiU(j)*dd_KUp(j) + dd_AU(j) *KxKwUp(j))/Kech_Up(j)
656 dphiV(j) = (phiV(j)*dd_KVp(j) + dd_AV(j) *KxKwVp(j))/Kech_Vp(j)
657 !
658 phiT_x(j)=phiT(j) - sigw(j)*dphiT(j)
659 phiT_w(j)=phiT(j) + sigx(j)*dphiT(j)
660 phiQ_x(j)=phiQ(j) - sigw(j)*dphiQ(j)
661 phiQ_w(j)=phiQ(j) + sigx(j)*dphiQ(j)
662 phiU_x(j)=phiU(j) - sigw(j)*dphiU(j)
663 phiU_w(j)=phiU(j) + sigx(j)*dphiU(j)
664 phiV_x(j)=phiV(j) - sigw(j)*dphiV(j)
665 phiV_w(j)=phiV(j) + sigx(j)*dphiV(j)
666 !
667 philat_x(j)=phiQ_x(j)*RLVTT
668 philat_w(j)=phiQ_w(j)*RLVTT
669 ENDDO ! j = 1,knon
670 !
671 DO j = 1,knon
672 q1_x(j) = AQ_x(j) + BQ_x(j)*phiQ_x(j)*dtime
673 q1_w(j) = AQ_w(j) + BQ_w(j)*phiQ_w(j)*dtime
674 ENDDO ! j = 1,knon
675 DO j = 1,knon
676 delta_qsurf(j) = (1.-beta(j))*(q1_w(j) - q1_x(j)) + dqsatsrf(j)
677 ENDDO ! j = 1,knon
678 !
679 !! Do j = 1,knon
680 !! print *,'XXXsplit : j, q1_x(j), AQ_x(j), BQ_x(j), phiQ_x(j) ', j, q1_x(j), AQ_x(j), BQ_x(j), phiQ_x(j)
681 !! print *,'XXXsplit : j, q1_w(j), AQ_w(j), BQ_w(j), phiQ_w(j) ', j, q1_w(j), AQ_w(j), BQ_w(j), phiQ_w(j)
682 !! ENDDO
683 !
684 IF (prt_level >=10 ) THEN
685 print *,' wx_pbl_split, phiT, dphiT, dTs ', phiT, dphiT, dTs
686 print *,' wx_pbl_split, phiQ, dphiQ, dqsatsrf ', phiQ, dphiQ, dqsatsrf
687 ENDIF
688 !
689 IF (prt_level >=10 ) THEN
690 !! print *,' wx_pbl_split, verif dqsatsrf = beta dqsatdT0 dTs '
691 !! print *,' wx_pbl_split, dqsatsrf, dqsatdT0*dTs ', dqsatsrf, dqsatdT0*dTs
692 ENDIF
693 !
694 !! IF (knon .NE. 0) THEN
695 !! call iophys_ecrit('sigw', 1,'sigw', '.',sigw)
696 !! call iophys_ecrit('phit', 1,'phit', 'W/m2',phit)
697 !! call iophys_ecrit('phit_w', 1,'phit_w', 'W/m2',phit_w)
698 !! call iophys_ecrit('phit_x', 1,'phit_x', 'W/m2',phit_x)
699 !! call iophys_ecrit('phiq', 1,'phiq', 'kg/m2/s',phiq)
700 !! call iophys_ecrit('phiq_w', 1,'phiq_w', 'kg/m2/s',phiq_w)
701 !! call iophys_ecrit('phiq_x', 1,'phiq_x', 'kg/m2/s',phiq_x)
702 !! call iophys_ecrit('q1_w', 1,'q1_w', '.',q1_w)
703 !! call iophys_ecrit('q1_x', 1,'q1_x', '.',q1_x)
704 !! ENDIF ! (knon .NE. 0)
705 !
706 RETURN
707
708 END SUBROUTINE wx_pbl_split
709
710 SUBROUTINE wx_pbl_check( knon, dtime, ypplay, ypaprs, &
711 sigw, beta, iflag_split, &
712 Ts0_b9, dTs09, &
713 qs_b9, Ts_b9, & ! yqsurf, Tsurf_new
714 dTs9, dqsatsrf9, &
715 AcoefT_x, AcoefT_w, &
716 BcoefT_x, BcoefT_w, &
717 AcoefT0, AcoefQ0, BcoefT0, BcoefQ0, &
718 AcoefT, AcoefQ, BcoefT, BcoefQ, &
719 phiT_b9, phiQ_b9, &
720 phiT_x9, phiT_w9, &
721 phiQ_x9, phiQ_w9 &
722 )
723 !
724
725 USE wx_pbl_var_mod
726
727 USE print_control_mod, ONLY: prt_level,lunout
728 !
729 INCLUDE "YOMCST.h"
730 INCLUDE "FCTTRE.h"
731 INCLUDE "YOETHF.h"
732 !
733 INTEGER, INTENT(IN) :: knon ! number of grid cells
734 REAL, INTENT(IN) :: dtime ! time step size (s)
735 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa)
736 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa)
737 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area
738 REAL, DIMENSION(knon), INTENT(IN) :: beta ! aridity factor
739 INTEGER, INTENT(IN) :: iflag_split
740 REAL, DIMENSION(knon), INTENT(IN) :: Ts0_b9, dTs09
741 REAL, DIMENSION(knon), INTENT(IN) :: qs_b9, Ts_b9 ! yqsurf, Tsurf_new
742 REAL, DIMENSION(knon), INTENT(IN) :: dTs9, dqsatsrf9
743 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT_x, AcoefT_w
744 REAL, DIMENSION(knon), INTENT(IN) :: BcoefT_x, BcoefT_w
745 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT0, AcoefQ0, BcoefT0, BcoefQ0
746 !
747 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT, AcoefQ, BcoefT, BcoefQ
748 REAL, DIMENSION(knon), INTENT(IN) :: phiT_b9, phiQ_b9
749 REAL, DIMENSION(knon), INTENT(IN) :: phiT_x9, phiT_w9
750 REAL, DIMENSION(knon), INTENT(IN) :: phiQ_x9, phiQ_w9
751 !
752 !! Local variables
753 INTEGER :: j
754 REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region
755 REAL, DIMENSION(knon) :: AcoefT_b, AcoefQ_b ! mean values of AcoefT and AcoefQ
756 REAL :: zzt, zzq, zzqsat
757 REAL :: zdelta, zcvm5, zcor, qsat
758 REAL, DIMENSION(knon) :: qsat_w, qsat_x
759 REAL, DIMENSION(knon) :: dqsatdT_w, dqsatdT_x
760 REAL, DIMENSION(knon) :: qsat_bs ! qsat(Ts_b)
761 REAL, DIMENSION(knon) :: qsat01, dqsatdT01
762 REAL, DIMENSION(knon) :: Ts_x, Ts_w, qs_x, qs_w
763 REAL, DIMENSION(knon) :: T1_x, T1_w, q1_x, q1_w
764 REAL, DIMENSION(knon) :: Rn_x, Rn_w
765 REAL, DIMENSION(knon) :: phiQ0_x, phiQ0_w
766 REAL, DIMENSION(knon) :: Ta, qa
767 REAL, DIMENSION(knon) :: qsatsrf_w, qsatsrf_x, qsatsrf_b
768 REAL, DIMENSION(knon) :: qsurf_w, qsurf_x
769 REAL :: dphiT, dphiQ
770 REAL :: dqsatsrf1
771 REAL :: phiT_w1, phiT_w2
772 REAL :: phiT_x1, phiT_x2
773 REAL :: phiQ_w1, phiQ_w2, phiQ_w3
774 REAL :: phiQ_x1, phiQ_x2, phiQ_x3
775 REAL :: phiT_b1, phiQ_b1
776 REAL :: Kech_Q_sw1, Kech_Q_sx1
777 REAL :: evap_pot
778
779 !----------------------------------------------------------------------------
780 ! Equations to be checked:
781 ! -----------------------
782 ! Input : Ts0_b, dTs0, Ts_b, dTs, qsatsrf_b, dqsatsrf,
783 ! phiT_b, phiQ_b, phiT_w, phiT_x, phiQ_w, phiQ_x,
784 !
785 ! AcoefT, AcoefQ, AcoefT_w, AcoefQ_w, AcoefT_x, AcoefQ_x,
786 ! BcoefT, BcoefQ, BcoefT_w, BcoefQ_w, BcoefT_x, BcoefQ_x
787 !
788 ! C_p T1_w = AcoefT_w + BcoefT_w phiT_w Delta t C_p T1_x = AcoefT_x + BcoefT_x phiT_x Delta t
789 ! q1_w = AQ_w + BQ_w phiQ_w Delta t q1_x = AQ_x + BQ_x phiQ_x Delta t
790 ! qsatsrf_w = beta qsat(Ts_w) qsatsrf_x = beta qsat(Ts_x)
791 ! qsurf_w = (1-beta) q1_w + qsatsrf_w qsurf_x = (1-beta) q1_x + qsatsrf_x
792 ! phiT_w = Kech_h_w C_p ( T1_w - Ts_w) phiT_x = Kech_h_x C_p ( T1_x - Ts_x)
793 ! phiT_w = Kech_T_pw ( AcoefT_w - C_p Ts_w) phiT_x = Kech_T_px ( AcoefT_x - C_p Ts_x)
794 ! phiq_w = Kech_h_w ( beta q1_w - qsatsrf_w) phiq_x = Kech_h_x ( beta q1_x - qsatsrf_x))
795 ! phiq_w = Kech_Q_sw (beta AQ_w -qsatsrf_w) phiq_x = Kech_Q_sx (beta AQ_x -qsatsrf_x)
796 ! phiq_w = Kech_h_w (q1_w - qsurf_w) phiq_x = Kech_h_x (q1_x - qsurf_x)
797 ! phiT_b = sigw phiT_w + sigx phiT_x dphiT = phiT_w - phiT_x
798 ! phiQ_b = sigw phiQ_w + sigx phiQ_x dphiQ = phiQ_w - phiQ_x
799 ! Ts_b = sigw Ts_w + sigx Ts_x dTs = Ts_w - Ts_x
800 ! qsatsrf_b = sigw qsatsrf_w + sigx qsatsrf_x
801 ! C_p Ta = AcoefT + BcoefT phiT_b Delta t
802 ! qa = AcoefQ + BcoefQ phiQ_b Delta t
803 ! phiT_b = Kech_h C_p (Ta - Ts_b)
804 ! phiQ_b = beta Kech_h (qa - qsatsrf_b)
805 ! dTs = sqrt(tau)/I (dphit + L_v dphiq + dR)
806
807 !----------------------------------------------------------------------------
808 !
809 !!
810 sigx(:) = 1.-sigw(:)
811 AcoefT_b(1:knon) = AcoefT_x(1:knon) + sigw(1:knon)*dd_AT(1:knon)
812 AcoefQ_b(1:knon) = AQ_x(1:knon) + sigw(1:knon)*dd_AQ(1:knon)
813
814 ! Compute the three qsat and dqsatdTs
815 ! ---------------------------------------------
816 !! C_p(1:knon) = RCpd
817 !! L_v(1:knon) = RLvtt
818 IF (prt_level >=10 ) THEN
819 print *,' AAAA wx_pbl_check, C_p(j), qsat0(j), Ts0(j) : ', C_p(:), qsat0(:), Ts0(:)
820 ENDIF ! (prt_level >=10 )
821 !
822 DO j = 1, knon
823 zdelta = MAX(0.,SIGN(1.,RTT-Ts0_b9(j)))
824 zcvm5 = R5LES*(1.-zdelta) + R5IES*zdelta
825 qsat = R2ES*FOEEW(Ts0_b9(j),zdelta)/ypaprs(j,1)
826 qsat = MIN(0.5,qsat)
827 zcor = 1./(1.-RETV*qsat)
828 qsat01(j) = fqsat*qsat*zcor
829 !! dqsatdT0(j) = FOEDE(Ts0_b(j),zdelta,zcvm5,qsat0(j),zcor)/RLVTT ! jyg 20210116
830 !! dqsatdT0(j) = (RLvtt*(1.-zdelta)+RLSTT*zdelta)*qsat0(j)/(Rv*Ts0_b(j)*Ts0_b(j))
831 dqsatdT01(j) = fqsat*FOEDE(Ts0_b9(j),zdelta,zcvm5,qsat01(j),zcor)
832 ENDDO
833 !
834 !--------------------------------------------------------------------------------------------------
835 IF (prt_level >=10 ) THEN
836 !
837 DO j = 1, knon
838 !
839 print *,'wx_pbl_check: Kech_h, Kech_q ', Kech_h(j), Kech_q(j)
840 !
841 Ta(j) = (AcoefT(j) + BcoefT(j)*phiT_b9(j)*dtime)/C_p(j)
842 qa(j) = AcoefQ(j) + BcoefQ(j)*phiQ_b9(j)*dtime
843 print *, 'wx_pbl_check: j, Ta, qa ', Ta(j), qa(j)
844 !
845 qsat_bs(j) = qsat01(j) + dqsatdT01(j)*(Ts_b9(j)-Ts0_b9(j))
846 !
847 print *,'wx_pbl_check: qsat01, qsat_bs ', j,qsat01(j), qsat_bs(j)
848 !
849 Ts_x(j) = Ts_b9(j) - sigw(j)*dTs9(j)
850 Ts_w(j) = Ts_b9(j) + sigx(j)*dTs9(j)
851 print *, 'wx_pbl_check: j, Ts_b9, Ts_w, Ts_x ', j, Ts_b9(j), Ts_w(j), Ts_x(j)
852 !
853 qsat_x(j) = qsat0_x(j) + dqsatdT0_x(j)*(Ts_x(j)-Ts0_x(j))
854 qsat_w(j) = qsat0_w(j) + dqsatdT0_w(j)*(Ts_w(j)-Ts0_w(j))
855 !
856 print *,'wx_pbl_check: qsat0_w, qsat0_x, qsat_w, qsat_x ', qsat0_w(j), qsat0_x(j), qsat_w(j), qsat_x(j)
857 !
858 T1_x(j) = (AcoefT_x(j) + BcoefT_x(j)*phiT_x9(j)*dtime) / C_p(j)
859 T1_w(j) = (AcoefT_w(j) + BcoefT_w(j)*phiT_w9(j)*dtime) / C_p(j)
860 print *, 'wx_pbl_check: j, T1_w, T1_x ', j, T1_w(j), T1_x(j)
861 !
862 q1_x(j) = AQ_x(j) + BQ_x(j)*phiQ_x9(j)*dtime
863 q1_w(j) = AQ_w(j) + BQ_w(j)*phiQ_w9(j)*dtime
864 print *, 'wx_pbl_check: j, q1_w, q1_x ', j, q1_w(j), q1_x(j)
865 !
866 qsatsrf_x(j) = beta(j)*qsat_x(j)
867 qsatsrf_w(j) = beta(j)*qsat_w(j)
868 qsatsrf_b(j) = sigw(j)*qsatsrf_w(j) + sigx(j)*qsatsrf_x(j)
869 !
870 dqsatsrf1 = qsatsrf_w(j) - qsatsrf_x(j)
871 print *, 'wx_pbl_check: j, qsatsrf_w, qsatsrf_x, dqsatsrf1, dqsatsrf9 ', &
872 qsatsrf_w(j), qsatsrf_x(j), dqsatsrf1, dqsatsrf9(j)
873 !
874 qsurf_x(j) = (1-beta(j))*q1_x(j) + qsatsrf_x(j)
875 qsurf_w(j) = (1-beta(j))*q1_w(j) + qsatsrf_w(j)
876 print *, 'wx_pbl_check: j, qsurf_w, qsurf_x ', j, qsurf_w(j), qsurf_x(j)
877 !
878 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
879 ! Test qsat01 = qsat0 et dqsatdT01 = dqsatdT0
880 !------------------------------------------------------------------------------------------------------
881 print *, 'wx_pbl_check: j, qsat01(j), qsat0(j) ', j, qsat01(j), qsat0(j)
882 print *, 'wx_pbl_check: j, dqsatdT01(j), dqsatdT0(j) ', j, dqsatdT01(j), dqsatdT0(j)
883 !
884 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
885 ! Test Kexh_Q_sw = Kech_q_w/(1.-beta*Kech_q_w*BcoefQ) Kexh_Q_sx = Kech_q_x/(1.-beta*Kech_q_x*BcoefQ)
886 !------------------------------------------------------------------------------------------------------
887 Kech_Q_sx1 = Kech_q_x(j)/(1.-beta(j)*Kech_q_x(j)*BQ_x(j)*dtime)
888 Kech_Q_sw1 = Kech_q_w(j)/(1.-beta(j)*Kech_q_w(j)*BQ_w(j)*dtime)
889 print *, 'wx_pbl_check: j, Kech_Q_sx1, Kech_Q_sx(j)', j, Kech_Q_sx1, Kech_Q_sx(j)
890 print *, 'wx_pbl_check: j, Kech_Q_sw1, Kech_Q_sw(j)', j, Kech_Q_sw1, Kech_Q_sw(j)
891 !
892 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
893 ! Test phiT_w = Kech_h_w*C_p(j)*(T1_w(j)-Ts_w(j)) phiT_x = Kech_h_x*C_p(j)*(T1_x(j)-Ts_x(j))
894 !-----------------------------------------------------
895 phiT_x1 = Kech_h_x(j)*C_p(j)*(T1_x(j)-Ts_x(j))
896 phiT_w1 = Kech_h_w(j)*C_p(j)*(T1_w(j)-Ts_w(j))
897 !
898 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
899 ! Test phiT_w = Kech_T_pw*(AcoefT_w(j)-C_p(j)*Ts_w(j)) phiT_x = Kech_T_px*(AcoefT_x(j)-C_p(j)*Ts_x(j))
900 !-----------------------------------------------------
901 phiT_x2 = Kech_T_px(j)*(AcoefT_x(j)-C_p(j)*Ts_x(j))
902 phiT_w2 = Kech_T_pw(j)*(AcoefT_w(j)-C_p(j)*Ts_w(j))
903 print *, 'wx_pbl_check: j, phiT_w1, phiT_w2, phiT_w9 ', j, phiT_w1, phiT_w2, phiT_w9(j)
904 print *, 'wx_pbl_check: j, phiT_x1, phiT_x2, phiT_x9 ', j, phiT_x1, phiT_x2, phiT_x9(j)
905 !
906 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
907 ! Test phiq_w = Kech_q_w ( beta q1_w - qsatsrf_w) phiq_x = Kech_q_x ( beta q1_x - qsatsrf_x))
908 !--------------------------------------------------------------
909 phiq_x1 = Kech_q_x(j)*( beta(j)*q1_x(j) - qsatsrf_x(j))
910 phiq_w1 = Kech_q_w(j)*( beta(j)*q1_w(j) - qsatsrf_w(j))
911 !
912 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
913 ! Test phiq_w = Kech_Q_sw (beta AQ_w -qsatsrf_w) phiq_x = Kech_Q_sx (beta AQ_x -qsatsrf_x)
914 !--------------------------------------------------------------
915 phiq_x2 = Kech_Q_sx(j)*(beta(j)*AQ_x(j) -qsatsrf_x(j))
916 phiq_w2 = Kech_Q_sw(j)*(beta(j)*AQ_w(j) -qsatsrf_w(j))
917 !
918 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
919 ! Test phiq_w = Kech_q_w ( q1_w - qsurf_w) phiq_x = Kech_q_x ( q1_x - qsurf_x))
920 !--------------------------------------------------------------
921 phiq_x3 = Kech_q_x(j)*( q1_x(j) - qsurf_x(j))
922 phiq_w3 = Kech_q_w(j)*( q1_w(j) - qsurf_w(j))
923 print *, 'wx_pbl_check: j, phiQ_w1, phiQ_w2, phiQ_w3, phiQ_w9 ', j, phiQ_w1, phiQ_w2, phiQ_w3, phiQ_w9(j)
924 print *, 'wx_pbl_check: j, phiQ_x1, phiQ_x2, phiQ_x3, phiQ_x9 ', j, phiQ_x1, phiQ_x2, phiQ_x3, phiQ_x9(j)
925 !
926 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
927 ! Test phiT_b = Kech_h C_p (Ta - Ts_b)
928 !--------------------------------------------------------------
929 phiT_b1 = Kech_h(j)*C_p(j)*(Ta(j) - Ts_b9(j))
930 print *, 'wx_pbl_check: j, phiT_b1, PhiT_b9 ', j, phiT_b1, PhiT_b9(j)
931 !
932 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
933 ! Test phiQ_b = beta Kech_q (qa - qsat_bs)
934 !--------------------------------------------------------------
935 evap_pot = Kech_q(j)*(qa(j) - qsat_bs(j))
936 phiQ_b1 = beta(j)*Kech_q(j)*(qa(j) - qsat_bs(j))
937 print *, 'wx_pbl_check: j, beta, evap_pot, phiQ_b1, PhiQ_b9 ', j, beta(j), evap_pot, phiQ_b1, PhiQ_b9(j)
938 !
939 !
940 ENDDO ! j = 1, knon
941
942 ENDIF ! (prt_level >=10 )
943 !--------------------------------------------------------------------------------------------------
944
945 RETURN
946
947 END SUBROUTINE wx_pbl_check
948
949 SUBROUTINE wx_pbl_dts_check( knon, dtime, ypplay, ypaprs, &
950 sigw, beta, iflag_split, &
951 Ts0_b9, dTs09, &
952 qs_b9, Ts_b9, & ! yqsurf, Tsurf_new
953 dqsatsrf9, dTs9, delta_qsurf9, &
954 AcoefT_x, AcoefT_w, &
955 BcoefT_x, BcoefT_w, &
956 AcoefT0, AcoefQ0, BcoefT0, BcoefQ0, &
957 AcoefT, AcoefQ, BcoefT, BcoefQ, &
958 HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn, &
959 phiT0_b9, dphiT09, phiQ0_b9, dphiQ09, Rn0_b9, dRn09, &
960 g_T, g_Q, &
961 Gamma_phiT, Gamma_phiQ, &
962 dTs_ins, dqsatsrf_ins, &
963 phiT_b9, phiQ_b9, &
964 phiT_x9, phiT_w9, &
965 phiQ_x9, phiQ_w9 &
966 )
967 !
968
969 USE wx_pbl_var_mod
970
971 USE print_control_mod, ONLY: prt_level,lunout
972 !
973 INCLUDE "YOMCST.h"
974 INCLUDE "FCTTRE.h"
975 INCLUDE "YOETHF.h"
976 !
977 INTEGER, INTENT(IN) :: knon ! number of grid cells
978 REAL, INTENT(IN) :: dtime ! time step size (s)
979 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypplay ! mid-layer pressure (Pa)
980 REAL, DIMENSION(knon,klev), INTENT(IN) :: ypaprs ! pressure at layer interfaces (pa)
981 REAL, DIMENSION(knon), INTENT(IN) :: sigw ! cold pools fractional area
982 REAL, DIMENSION(knon), INTENT(IN) :: beta ! aridity factor
983 INTEGER, INTENT(IN) :: iflag_split
984 REAL, DIMENSION(knon), INTENT(IN) :: Ts0_b9, dTs09
985 REAL, DIMENSION(knon), INTENT(IN) :: qs_b9, Ts_b9 ! yqsurf, Tsurf_new
986 REAL, DIMENSION(knon), INTENT(IN) :: dTs9, dqsatsrf9
987 REAL, DIMENSION(knon), INTENT(IN) :: delta_qsurf9
988 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT_x, AcoefT_w
989 REAL, DIMENSION(knon), INTENT(IN) :: BcoefT_x, BcoefT_w
990 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT0, AcoefQ0, BcoefT0, BcoefQ0
991 !
992 REAL, DIMENSION(knon), INTENT(IN) :: AcoefT, AcoefQ, BcoefT, BcoefQ
993 REAL, DIMENSION(knon), INTENT(IN) :: HTphiT_b, dd_HTphiT, HTphiQ_b, dd_HTphiQ, HTRn_b, dd_HTRn
994 REAL, DIMENSION(knon), INTENT(IN) :: phiT0_b9, dphiT09, phiQ0_b9, dphiQ09, Rn0_b9, dRn09
995 REAL, DIMENSION(knon), INTENT(IN) :: g_T, g_Q
996 REAL, DIMENSION(knon), INTENT(IN) :: Gamma_phiT, Gamma_phiQ
997 REAL, DIMENSION(knon), INTENT(IN) :: dTs_ins, dqsatsrf_ins
998 REAL, DIMENSION(knon), INTENT(IN) :: phiT_b9, phiQ_b9
999 REAL, DIMENSION(knon), INTENT(IN) :: phiT_x9, phiT_w9
1000 REAL, DIMENSION(knon), INTENT(IN) :: phiQ_x9, phiQ_w9
1001 !
1002 !! Local variables
1003 INTEGER :: j
1004 REAL, DIMENSION(knon) :: sigx ! fractional area of (x) region
1005 REAL, DIMENSION(knon) :: AcoefT_b, AcoefQ_b ! mean values of AcoefT and AcoefQ
1006 REAL :: zzt, zzq, zzqsat
1007 REAL :: zdelta, zcvm5, zcor, qsat
1008 REAL, DIMENSION(knon) :: qsat_w, qsat_x
1009 REAL, DIMENSION(knon) :: Ts_x, Ts_w, qs_x, qs_w
1010 REAL, DIMENSION(knon) :: T1_x, T1_w, q1_x, q1_w
1011 REAL, DIMENSION(knon) :: Rn_x, Rn_w
1012 REAL, DIMENSION(knon) :: Rn_b, dRn
1013 REAL, DIMENSION(knon) :: phiQ0_x, phiQ0_w
1014 REAL, DIMENSION(knon) :: Ta, qa
1015 REAL, DIMENSION(knon) :: err_phiT_w, err_phiT_x
1016 REAL, DIMENSION(knon) :: err_phiq_w, err_phiq_x
1017 REAL, DIMENSION(knon) :: err_phiT_b
1018 REAL, DIMENSION(knon) :: err_phiQ_b
1019 REAL, DIMENSION(knon) :: err2_phiT_b
1020 REAL :: T1A_x, T1A_w, q1A_x, q1A_w
1021 REAL :: qsatsrf_w, qsatsrf_x, qsatsrfb, qsbA
1022 REAL :: dphiT, dphiQ
1023 REAL :: dphiT_H, dphiQ_H
1024 REAL :: phiQ_pot
1025 REAL :: phiQ_w_m_phiQ0_w
1026 REAL :: phiQ_x_m_phiQ0_x
1027 REAL :: dphiQ_m_dphiQ0
1028 REAL :: dphiT_m_dphiT0
1029 REAL :: dRN_m_dRn0
1030 REAL :: phiTb_m_phiT0b
1031
1032 !----------------------------------------------------------------------------
1033 ! Equations to be checked:
1034 ! -----------------------
1035 ! Input : Ts0_b, dTs0, Ts_b, dTs, qsatsrf_b, dqsatsrf,
1036 ! phiT_b, phiQ_b, phiT_w, phiT_x, phiQ_w, phiQ_x,
1037 !
1038 ! AcoefT, AcoefQ, AcoefT_w, AcoefQ_w, AcoefT_x, AcoefQ_x,
1039 ! BcoefT, BcoefQ, BcoefT_w, BcoefQ_w, BcoefT_x, BcoefQ_x
1040 !
1041 ! Ts_w = Ts_b + sigx dTs Ts_x = Ts_b - sigw dTs
1042 ! T1_w = AcoefT_w + BcoefT_w phiT_w Delta t T1_x = AcoefT_x + BcoefT_x phiT_x Delta t
1043 ! q1_w = AcoefQ_w + BcoefQ_w phiQ_w Delta t q1_x = AcoefQ_x + BcoefQ_x phiQ_x Delta t
1044 ! phiT_w = Kech_h_w ( T1_w - Ts_w) phiT_x = Kech_h_x ( T1_x - Ts_x)
1045 ! phiq_w = beta Kech_h_w ( q1_w - qsat(Ts_w)) phiq_x = beta Kech_h_x ( q1_x - qsat(Ts_x))
1046 ! phiT_b = sigw phiT_w + sigx phiT_x dphiT = phiT_w - phiT_x
1047 ! phiQ_b = sigw phiQ_w + sigx phiQ_x dphiQ = phiQ_w - phiQ_x
1048 ! Ts_b = sigw Ts_w + sigx Ts_x dTs = Ts_w - Ts_x
1049 ! Ta = AcoefT + BcoefT phiT_b Delta t
1050 ! qa = AcoefQ + BcoefQ phiQ_b Delta t
1051 ! phiT_b = Kech_h (Ta - Ts_b)
1052 ! phiQ_b = beta Kech_h (qa - qsat(Ts_b))
1053 ! dTs = sqrt(tau)/I (dphit + L_v dphiq + dR)
1054
1055 !----------------------------------------------------------------------------
1056 !
1057 !!
1058 sigx(:) = 1.-sigw(:)
1059 AcoefT_b(1:knon) = AcoefT_x(1:knon) + sigw(1:knon)*dd_AT(1:knon)
1060 AcoefQ_b(1:knon) = AQ_x(1:knon) + sigw(1:knon)*dd_AQ(1:knon)
1061
1062 IF (prt_level >=10 ) THEN
1063 print *,'->wx_pbl_dts_check, HTphiT_b, HTphiQ_b, HTRn_b ', &
1064 HTphiT_b, HTphiQ_b, HTRn_b
1065 print *,'->wx_pbl_dts_check, dd_HTphiT, dd_HTphiQ, dd_HTRn ', &
1066 dd_HTphiT, dd_HTphiQ, dd_HTRn
1067 ENDIF ! (prt_level >=10 )
1068 !
1069 ! Compute the three qsat and dqsatdTs
1070 ! ---------------------------------------------
1071 !! print *,' AAAA wx_pbl_dts_check, C_p(j), qsat0(j), Ts0(j) : ', &
1072 !! (C_p(j), qsat0(j), Ts0(j), j = 1,knon)
1073 !
1074 !
1075 !--------------------------------------------------------------------------------------------------
1076 IF (prt_level >=10 ) THEN
1077 !
1078 DO j = 1, knon
1079 Ts_x(j) = Ts_b9(j) - sigw(j)*dTs9(j)
1080 Ts_w(j) = Ts_b9(j) + sigx(j)*dTs9(j)
1081 print *, 'wx_pbl_dts_check: j, Ts_b9, Ts_w, Ts_x ', j, Ts_b9(j), Ts_w(j), Ts_x(j)
1082 !
1083 qsat_x(j) = qsat0_x(j) + dqsatdT0_x(j)*(Ts_x(j)-Ts0_x(j))
1084 qsat_w(j) = qsat0_w(j) + dqsatdT0_w(j)*(Ts_w(j)-Ts0_w(j))
1085 !
1086 T1_x(j) = (AcoefT_x(j) + BcoefT_x(j)*phiT_x9(j)*dtime) / C_p(j)
1087 T1_w(j) = (AcoefT_w(j) + BcoefT_w(j)*phiT_w9(j)*dtime) / C_p(j)
1088 print *, 'wx_pbl_dts_check: j, T1_w, T1_x ', j, T1_w(j), T1_x(j)
1089 !
1090 q1_x(j) = AQ_x(j) + BQ_x(j)*phiQ_x9(j)*dtime
1091 q1_w(j) = AQ_w(j) + BQ_w(j)*phiQ_w9(j)*dtime
1092 print *, 'wx_pbl_dts_check: j, q1_w, q1_x ', j, q1_w(j), q1_x(j)
1093 !
1094 Rn_x(j) = eps_1*Rsigma*T1_x(j)**4 - Rsigma*Ts_x(j)**4
1095 Rn_w(j) = eps_1*Rsigma*T1_w(j)**4 - Rsigma*Ts_w(j)**4
1096 Rn_b(j) = sigw(j)*Rn_w(j) + sigx(j)*Rn_x(j)
1097 dRn(j) = dRn09(j) - ( HTRn_b(j) &
1098 +(sigx(j)-sigw(j))*dd_HTRn(j) &
1099 -sigw(j)*sigx(j)*dd_HTRn(j)*dd_HTphiT(j)/HTphiT_b(j) &
1100 )*(dTs9(j)-dTs09(j)) &
1101 + dd_HTRn(j)/HTphiT_b(j)*(phiT_b9(j)-phiT0_b9(j))
1102 !
1103 print *,'wx_pbl_dts_check, dphiT, L_v*dphiQ, dRn, dTs ', &
1104 phiT_w9(j)-phiT_x9(j), L_v(j)*(phiQ_w9(j)-phiQ_x9(j)), dRn(j), dTs9(j)
1105 !
1106 phiQ0_x(j) = PhiQ0_b9(j) - sigw(j)*dphiQ09(j)
1107 phiQ0_w(j) = PhiQ0_b9(j) + sigx(j)*dphiQ09(j)
1108 !
1109 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1110 ! Test phiQ_w-phiQ0_w = -beta*Kech_Q_sw*dqsatdT_w*(Ts_w-Ts0_w)
1111 !--------------------------------------------------------------
1112 print *,'wx_pbl_dts_check: beta(j), Kech_Q_sw(j), dqsatdT0_w(j), Ts_w(j), Ts0_w(j) ', &
1113 beta(j), Kech_Q_sw(j), dqsatdT0_w(j), Ts_w(j), Ts0_w(j)
1114 phiQ_w_m_phiQ0_w = -beta(j)*Kech_Q_sw(j)*dqsatdT0_w(j)*(Ts_w(j)-Ts0_w(j))
1115 print *,'wx_pbl_dts_check: j, phiQ_w9-phiQ0_w, phiQ_w_m_phiQ0_w ', &
1116 j, phiQ_w9(j)-phiQ0_w(j), phiQ_w_m_phiQ0_w
1117 !
1118 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1119 ! Test phiQ_x-phiQ0_x = -beta*Kech_Q_sx*dqsatdT_x*(Ts_x-Ts0_x)
1120 !--------------------------------------------------------------
1121 phiQ_x_m_phiQ0_x = -beta(j)*Kech_Q_sx(j)*dqsatdT0_x(j)*(Ts_x(j)-Ts0_x(j))
1122 print *,'wx_pbl_dts_check: j, phiQ_x9-phiQ0_x, phiQ_x_m_phiQ0_x ', &
1123 j, phiQ_x9(j)-phiQ0_x(j), phiQ_x_m_phiQ0_x
1124 !
1125 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1126 ! Test dphiT-dphiT0 = -(HTphiT_b+(sigx-sigw)*dd_HTphiT)*(dTs-dTs0) - dd_HTphiT*(Ts_b-Ts0_b)
1127 !-------------------------------------------------------------------------------------------
1128 dphiT = phiT_w9(j) - phiT_x9(j)
1129 dphiT_m_dphiT0 = -(HTphiT_b(j)+(sigx(j)-sigw(j))*dd_HTphiT(j))*(dTs9(j)-dTs09(j)) &
1130 - dd_HTphiT(j)*(Ts_b9(j)-Ts0_b9(j))
1131 print *,'wx_pbl_dts_check: j, dphiT-dphiT09, dphiT_m_dphiT0 ',j, dphiT-dphiT09(j), dphiT_m_dphiT0
1132 !
1133 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1134 ! Test dphiQ-dphiQ0 = -(HTphiQ_b+(sigx-sigw)*dd_HTphiQ)*(dTs-dTs0) - dd_HTphiQ*(Ts_b-Ts0_b)
1135 !-------------------------------------------------------------------------------------------
1136 dphiQ = phiQ_w9(j) - phiQ_x9(j)
1137 dphiQ_m_dphiQ0 = -(HTphiQ_b(j)+(sigx(j)-sigw(j))*dd_HTphiQ(j))*(dTs9(j)-dTs09(j)) &
1138 - dd_HTphiQ(j)*(Ts_b9(j)-Ts0_b9(j))
1139 print *,'wx_pbl_dts_check: j, dphiQ-dphiQ09, dphiQ_m_dphiQ0 ',j, dphiQ-dphiQ09(j), dphiQ_m_dphiQ0
1140 !
1141 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1142 ! Test dRn-dRn0 = -(HTRn_b+(sigx-sigw)*dd_HTRn)*(dTs-dTs0) - dd_HTRn*(Ts_b-Ts0_b)
1143 !-------------------------------------------------------------------------------------------
1144 dRn_m_dRn0 = -(HTRn_b(j)+(sigx(j)-sigw(j))*dd_HTRn(j))*(dTs9(j)-dTs09(j)) &
1145 - dd_HTRn(j)*(Ts_b9(j)-Ts0_b9(j))
1146 print *,'wx_pbl_dts_check: j, dRn-dRn09, dRn_m_dRn0 ',j, dRn-dRn09(j), dRn_m_dRn0
1147 !
1148 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1149 ! Test phiT_b-phiT0_b = -sigx*sigw*dd_HTphiT*(dTs-dTs0) - HTphiT_b*(Ts_b-Ts0_b)
1150 !-------------------------------------------------------------------------------
1151 phiTb_m_phiT0b = -sigx(j)*sigw(j)*dd_HTphiT(j)*(dTs9(j)-dTs09(j)) - HTphiT_b(j)*(Ts_b9(j)-Ts0_b9(j))
1152 print *,'wx_pbl_dts_check: j, phiT_b9-phiT0_b9, phiTb_m_phiT0b ',j ,phiT_b9(j)-phiT0_b9(j), phiTb_m_phiT0b
1153 !
1154 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1155 ! Test phiT_w, phiT_x, dphiT from HTphiT
1156 !------------------------------------------
1157 ! phiT_w = Kech_h_w C_p ( T1_w - Ts_w) phiT_x = Kech_h_x C_p ( T1_x - Ts_x)
1158 err_phiT_x(j) = Kech_h_x(j)*C_p(j)*(T1_x(j) - Ts_x(j)) - phiT_x9(j)
1159 err_phiT_w(j) = Kech_h_w(j)*C_p(j)*(T1_w(j) - Ts_w(j)) - phiT_w9(j)
1160 print *, 'wx_pbl_dts_check: j, phiT_w9, phiT_x9, err_phiT_w, err_phiT_x ', &
1161 j, phiT_w9(j), phiT_x9(j), err_phiT_w(j), err_phiT_x(j)
1162 dphiT = phiT_w9(j) - phiT_x9(j)
1163 dphiT_H = dphiT09(j) - ( HTphiT_b(j) &
1164 +(sigx(j)-sigw(j))*dd_HTphiT(j) &
1165 -sigw(j)*sigx(j)*dd_HTphiT(j)*dd_HTphiT(j)/HTphiT_b(j) &
1166 )*(dTs9(j)-dTs09(j)) &
1167 + dd_HTphiT(j)/HTphiT_b(j)*(phiT_b9(j)-phiT0_b9(j))
1168 print *,'wx_pbl_dts_check: j, dphiT, dphiT_H ', j, dphiT, dphiT_H
1169
1170 !
1171 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1172 ! Test phiq_w, phiq_x, dphiq from HTphiq
1173 !------------------------------------------
1174 !
1175 ! phiq_w = beta Kech_q_w ( q1_w - qsat(Ts_w)) phiq_x = beta Kech_q_x ( q1_x - qsat(Ts_x))
1176 err_phiq_x(j) = beta(j)*Kech_q_x(j)*( q1_x(j) - qsat_x(j)) - phiq_x9(j)
1177 err_phiq_w(j) = beta(j)*Kech_q_w(j)*( q1_w(j) - qsat_w(j)) - phiq_w9(j)
1178 dphiQ = phiQ_w9(j) - phiQ_x9(j)
1179 dphiQ_H = dphiQ09(j) - ( HTphiQ_b(j) &
1180 +(sigx(j)-sigw(j))*dd_HTphiQ(j) &
1181 -sigw(j)*sigx(j)*dd_HTphiQ(j)*dd_HTphiT(j)/HTphiT_b(j) &
1182 )*(dTs9(j)-dTs09(j)) &
1183 + dd_HTphiQ(j)/HTphiT_b(j)*(phiT_b9(j)-phiT0_b9(j))
1184 print *,'wx_pbl_dts_check: j, dphiQ, dphiQ_H ', j, dphiQ, dphiQ_H
1185 !
1186 ! phiT_b = sigw phiT_w + sigx phiT_x dphiT = phiT_w - phiT_x
1187 err_phiT_b(j) = sigw(j)*phiT_w9(j) + sigx(j)*phiT_x9(j) - phiT_b9(j)
1188 !
1189 ! phiQ_b = sigw phiQ_w + sigx phiQ_x dphiQ = phiQ_w - phiQ_x
1190 err_phiQ_b(j) = sigw(j)*phiQ_w9(j) + sigx(j)*phiQ_x9(j) - phiQ_b9(j)
1191 !
1192 ! Ta = AcoefT + BcoefT phiT_b Delta t
1193 ! phiT_b = Kech_h C_p (Ta - Ts_b)
1194 Ta(j) = (AcoefT(j) + BcoefT(j)*phiT_b9(j)*dtime) / C_p(j)
1195 err2_phiT_b(j) = Kech_h(j)*C_p(j)*(Ta(j) - Ts_b9(j)) - phiT_b9(j)
1196 print *, 'wx_pbl_dts_check: j, Ta, phiT_b9, err2_phiT_b ', &
1197 j, Ta(j), phiT_b9(j), err2_phiT_b(j)
1198 !
1199 ENDDO ! j = 1, knon
1200
1201 ENDIF ! (prt_level >=10 )
1202 !--------------------------------------------------------------------------------------------------
1203 RETURN
1204
1205 END SUBROUTINE wx_pbl_dts_check
1206
1207 SUBROUTINE wx_evappot(knon, q1, Ts, evap_pot)
1208
1209 USE wx_pbl_var_mod
1210
1211 INTEGER, INTENT(IN) :: knon ! number of grid cells
1212 REAL, DIMENSION(knon), INTENT(IN) :: q1 ! specific humidity in layer 1
1213 REAL, DIMENSION(knon), INTENT(IN) :: Ts ! surface temperature
1214 !
1215 REAL, DIMENSION(knon), INTENT(OUT) :: evap_pot ! potential evaporation
1216 !
1217 INTEGER :: j
1218 REAL :: qsat_bs
1219 !
1220 DO j = 1,knon
1221 evap_pot(j) = Kech_q(j)*(qsat0(j)+dqsatdT0(j)*(Ts(j)-Ts0(j))-q1(j))
1222 !
1223 qsat_bs = qsat0(j)+dqsatdT0(j)*(Ts(j)-Ts0(j))
1224 !! print *,'wx_evappot : Kech_q, qsat_bs, qa, evap_pot ', Kech_q(j), qsat_bs, q1(j), evap_pot(j)
1225 ENDDO
1226 !
1227 RETURN
1228 END SUBROUTINE wx_evappot
1229
1230 END MODULE wx_pbl_mod
1231
1232