GCC Code Coverage Report


Directory: ./
File: phys/wake.f90
Date: 2022-01-11 19:19:34
Exec Total Coverage
Lines: 858 965 88.9%
Branches: 671 890 75.4%

Line Branch Exec Source
1
2 ! $Id: wake.F90 3648 2020-03-16 15:36:59Z jghattas $
3
4 768437836 SUBROUTINE wake(znatsurf, p, ph, pi, dtime, &
5 te0, qe0, omgb, &
6 dtdwn, dqdwn, amdwn, amup, dta, dqa, wgen, &
7 sigd_con, Cin, &
8 480 deltatw, deltaqw, sigmaw, awdens, wdens, & ! state variables
9 dth, hw, wape, fip, gfl, &
10 dtls, dqls, ktopw, omgbdth, dp_omgb, tu, qu, &
11 dtke, dqke, omg, dp_deltomg, spread, cstar, &
12 d_deltat_gw, &
13 d_deltatw2, d_deltaqw2, d_sigmaw2, d_awdens2, d_wdens2) ! tendencies
14
15
16 ! **************************************************************
17 ! *
18 ! WAKE *
19 ! retour a un Pupper fixe *
20 ! *
21 ! written by : GRANDPEIX Jean-Yves 09/03/2000 *
22 ! modified by : ROEHRIG Romain 01/29/2007 *
23 ! **************************************************************
24
25 USE ioipsl_getin_p_mod, ONLY : getin_p
26 USE dimphy
27 use mod_phys_lmdz_para
28 USE print_control_mod, ONLY: prt_level
29 IMPLICIT NONE
30 ! ============================================================================
31
32
33 ! But : Decrire le comportement des poches froides apparaissant dans les
34 ! grands systemes convectifs, et fournir l'energie disponible pour
35 ! le declenchement de nouvelles colonnes convectives.
36
37 ! State variables :
38 ! deltatw : temperature difference between wake and off-wake regions
39 ! deltaqw : specific humidity difference between wake and off-wake regions
40 ! sigmaw : fractional area covered by wakes.
41 ! wdens : number of wakes per unit area
42
43 ! Variable de sortie :
44
45 ! wape : WAke Potential Energy
46 ! fip : Front Incident Power (W/m2) - ALP
47 ! gfl : Gust Front Length per unit area (m-1)
48 ! dtls : large scale temperature tendency due to wake
49 ! dqls : large scale humidity tendency due to wake
50 ! hw : wake top hight (given by hw*deltatw(1)/2=wape)
51 ! dp_omgb : vertical gradient of large scale omega
52 ! awdens : densite de poches actives
53 ! wdens : densite de poches
54 ! omgbdth: flux of Delta_Theta transported by LS omega
55 ! dtKE : differential heating (wake - unpertubed)
56 ! dqKE : differential moistening (wake - unpertubed)
57 ! omg : Delta_omg =vertical velocity diff. wake-undist. (Pa/s)
58 ! dp_deltomg : vertical gradient of omg (s-1)
59 ! spread : spreading term in d_t_wake and d_q_wake
60 ! deltatw : updated temperature difference (T_w-T_u).
61 ! deltaqw : updated humidity difference (q_w-q_u).
62 ! sigmaw : updated wake fractional area.
63 ! d_deltat_gw : delta T tendency due to GW
64
65 ! Variables d'entree :
66
67 ! aire : aire de la maille
68 ! te0 : temperature dans l'environnement (K)
69 ! qe0 : humidite dans l'environnement (kg/kg)
70 ! omgb : vitesse verticale moyenne sur la maille (Pa/s)
71 ! dtdwn: source de chaleur due aux descentes (K/s)
72 ! dqdwn: source d'humidite due aux descentes (kg/kg/s)
73 ! dta : source de chaleur due courants satures et detrain (K/s)
74 ! dqa : source d'humidite due aux courants satures et detra (kg/kg/s)
75 ! wgen : number of wakes generated per unit area and per sec (/m^2/s)
76 ! amdwn: flux de masse total des descentes, par unite de
77 ! surface de la maille (kg/m2/s)
78 ! amup : flux de masse total des ascendances, par unite de
79 ! surface de la maille (kg/m2/s)
80 ! sigd_con:
81 ! Cin : convective inhibition
82 ! p : pressions aux milieux des couches (Pa)
83 ! ph : pressions aux interfaces (Pa)
84 ! pi : (p/p_0)**kapa (adim)
85 ! dtime: increment temporel (s)
86
87 ! Variables internes :
88
89 ! rhow : masse volumique de la poche froide
90 ! rho : environment density at P levels
91 ! rhoh : environment density at Ph levels
92 ! te : environment temperature | may change within
93 ! qe : environment humidity | sub-time-stepping
94 ! the : environment potential temperature
95 ! thu : potential temperature in undisturbed area
96 ! tu : temperature in undisturbed area
97 ! qu : humidity in undisturbed area
98 ! dp_omgb: vertical gradient og LS omega
99 ! omgbw : wake average vertical omega
100 ! dp_omgbw: vertical gradient of omgbw
101 ! omgbdq : flux of Delta_q transported by LS omega
102 ! dth : potential temperature diff. wake-undist.
103 ! th1 : first pot. temp. for vertical advection (=thu)
104 ! th2 : second pot. temp. for vertical advection (=thw)
105 ! q1 : first humidity for vertical advection
106 ! q2 : second humidity for vertical advection
107 ! d_deltatw : terme de redistribution pour deltatw
108 ! d_deltaqw : terme de redistribution pour deltaqw
109 ! deltatw0 : deltatw initial
110 ! deltaqw0 : deltaqw initial
111 ! hw0 : wake top hight (defined as the altitude at which deltatw=0)
112 ! amflux : horizontal mass flux through wake boundary
113 ! wdens_ref: initial number of wakes per unit area (3D) or per
114 ! unit length (2D), at the beginning of each time step
115 ! Tgw : 1 sur la p�riode de onde de gravit�
116 ! Cgw : vitesse de propagation de onde de gravit�
117 ! LL : distance entre 2 poches
118
119 ! -------------------------------------------------------------------------
120 ! D�claration de variables
121 ! -------------------------------------------------------------------------
122
123 include "YOMCST.h"
124 include "cvthermo.h"
125
126 ! Arguments en entree
127 ! --------------------
128
129 INTEGER, DIMENSION (klon), INTENT(IN) :: znatsurf
130 REAL, DIMENSION (klon, klev), INTENT(IN) :: p, pi
131 REAL, DIMENSION (klon, klev+1), INTENT(IN) :: ph
132 REAL, DIMENSION (klon, klev), INTENT(IN) :: omgb
133 REAL, INTENT(IN) :: dtime
134 REAL, DIMENSION (klon, klev), INTENT(IN) :: te0, qe0
135 REAL, DIMENSION (klon, klev), INTENT(IN) :: dtdwn, dqdwn
136 REAL, DIMENSION (klon, klev), INTENT(IN) :: amdwn, amup
137 REAL, DIMENSION (klon, klev), INTENT(IN) :: dta, dqa
138 REAL, DIMENSION (klon), INTENT(IN) :: wgen
139 REAL, DIMENSION (klon), INTENT(IN) :: sigd_con
140 REAL, DIMENSION (klon), INTENT(IN) :: Cin
141
142 !
143 ! Input/Output
144 ! State variables
145 REAL, DIMENSION (klon, klev), INTENT(INOUT) :: deltatw, deltaqw
146 REAL, DIMENSION (klon), INTENT(INOUT) :: sigmaw
147 REAL, DIMENSION (klon), INTENT(INOUT) :: awdens
148 REAL, DIMENSION (klon), INTENT(INOUT) :: wdens
149
150 ! Sorties
151 ! --------
152
153 REAL, DIMENSION (klon, klev), INTENT(OUT) :: dth
154 REAL, DIMENSION (klon, klev), INTENT(OUT) :: tu, qu
155 REAL, DIMENSION (klon, klev), INTENT(OUT) :: dtls, dqls
156 REAL, DIMENSION (klon, klev), INTENT(OUT) :: dtke, dqke
157 REAL, DIMENSION (klon, klev), INTENT(OUT) :: spread ! unused (jyg)
158 REAL, DIMENSION (klon, klev), INTENT(OUT) :: omgbdth, omg
159 REAL, DIMENSION (klon, klev), INTENT(OUT) :: dp_omgb, dp_deltomg
160 REAL, DIMENSION (klon, klev), INTENT(OUT) :: d_deltat_gw
161 REAL, DIMENSION (klon), INTENT(OUT) :: hw, wape, fip, gfl, cstar
162 INTEGER, DIMENSION (klon), INTENT(OUT) :: ktopw
163 ! Tendencies of state variables
164 REAL, DIMENSION (klon, klev), INTENT(OUT) :: d_deltatw2, d_deltaqw2
165 REAL, DIMENSION (klon), INTENT(OUT) :: d_sigmaw2, d_awdens2, d_wdens2
166
167 ! Variables internes
168 ! -------------------
169
170 ! Variables � fixer
171 INTEGER, SAVE :: igout
172 !$OMP THREADPRIVATE(igout)
173 LOGICAL, SAVE :: first = .TRUE.
174 !$OMP THREADPRIVATE(first)
175 !jyg<
176 !! REAL, SAVE :: stark, wdens_ref, coefgw, alpk
177 REAL, SAVE, DIMENSION(2) :: wdens_ref
178 REAL, SAVE :: stark, coefgw, alpk
179 !>jyg
180 REAL, SAVE :: crep_upper, crep_sol
181 !$OMP THREADPRIVATE(stark, wdens_ref, coefgw, alpk, crep_upper, crep_sol)
182
183 REAL, SAVE :: tau_cv
184 !$OMP THREADPRIVATE(tau_cv)
185 REAL, SAVE :: rzero, aa0 ! minimal wake radius and area
186 !$OMP THREADPRIVATE(rzero, aa0)
187
188 LOGICAL, SAVE :: flag_wk_check_trgl
189 !$OMP THREADPRIVATE(flag_wk_check_trgl)
190 INTEGER, SAVE :: iflag_wk_check_trgl
191 !$OMP THREADPRIVATE(iflag_wk_check_trgl)
192 INTEGER, SAVE :: iflag_wk_pop_dyn
193 !$OMP THREADPRIVATE(iflag_wk_pop_dyn)
194
195 REAL :: delta_t_min
196 INTEGER :: nsub
197 REAL :: dtimesub
198 REAL, SAVE :: wdensmin
199 !$OMP THREADPRIVATE(wdensmin)
200 REAL, SAVE :: sigmad, hwmin, wapecut, cstart
201 !$OMP THREADPRIVATE(sigmad, hwmin, wapecut, cstart)
202 REAL, SAVE :: sigmaw_max
203 !$OMP THREADPRIVATE(sigmaw_max)
204 REAL, SAVE :: dens_rate
205 !$OMP THREADPRIVATE(dens_rate)
206 REAL :: wdens0
207 ! IM 080208
208 960 LOGICAL, DIMENSION (klon) :: gwake
209
210 ! Variables de sauvegarde
211 960 REAL, DIMENSION (klon, klev) :: deltatw0
212 960 REAL, DIMENSION (klon, klev) :: deltaqw0
213 960 REAL, DIMENSION (klon, klev) :: te, qe
214 !! REAL, DIMENSION (klon) :: sigmaw1
215
216 ! Variables liees a la dynamique de population
217 960 REAL, DIMENSION(klon) :: act
218 960 REAL, DIMENSION(klon) :: rad_wk, tau_wk_inv
219 960 REAL, DIMENSION(klon) :: f_shear
220 960 REAL, DIMENSION(klon) :: drdt
221 960 REAL, DIMENSION(klon) :: d_sig_gen, d_sig_death, d_sig_col
222 960 REAL, DIMENSION(klon) :: wape1_act, wape2_act
223 960 LOGICAL, DIMENSION (klon) :: kill_wake
224 INTEGER, SAVE :: iflag_wk_act
225 !$OMP THREADPRIVATE(iflag_wk_act)
226 REAL :: drdt_pos
227 REAL :: tau_wk_inv_min
228
229 ! Variables pour les GW
230 960 REAL, DIMENSION (klon) :: ll
231 960 REAL, DIMENSION (klon, klev) :: n2
232 960 REAL, DIMENSION (klon, klev) :: cgw
233 960 REAL, DIMENSION (klon, klev) :: tgw
234
235 ! Variables liees au calcul de hw
236 960 REAL, DIMENSION (klon) :: ptop_provis, ptop, ptop_new
237 960 REAL, DIMENSION (klon) :: sum_dth
238 960 REAL, DIMENSION (klon) :: dthmin
239 960 REAL, DIMENSION (klon) :: z, dz, hw0
240 960 INTEGER, DIMENSION (klon) :: ktop, kupper
241
242 ! Variables liees au test de la forme triangulaire du profil de Delta_theta
243 960 REAL, DIMENSION (klon) :: sum_half_dth
244 960 REAL, DIMENSION (klon) :: dz_half
245
246 ! Sub-timestep tendencies and related variables
247 960 REAL, DIMENSION (klon, klev) :: d_deltatw, d_deltaqw
248 960 REAL, DIMENSION (klon, klev) :: d_te, d_qe
249 960 REAL, DIMENSION (klon) :: d_awdens, d_wdens
250 960 REAL, DIMENSION (klon) :: d_sigmaw, alpha
251 960 REAL, DIMENSION (klon) :: q0_min, q1_min
252 960 LOGICAL, DIMENSION (klon) :: wk_adv, ok_qx_qw
253 REAL, SAVE :: epsilon
254 !$OMP THREADPRIVATE(epsilon)
255 DATA epsilon/1.E-15/
256
257 ! Autres variables internes
258 INTEGER ::isubstep, k, i
259
260 REAL :: wdens_targ
261 REAL :: sigmaw_targ
262
263 960 REAL, DIMENSION (klon) :: sum_thu, sum_tu, sum_qu, sum_thvu
264 960 REAL, DIMENSION (klon) :: sum_dq, sum_rho
265 960 REAL, DIMENSION (klon) :: sum_dtdwn, sum_dqdwn
266 960 REAL, DIMENSION (klon) :: av_thu, av_tu, av_qu, av_thvu
267 960 REAL, DIMENSION (klon) :: av_dth, av_dq, av_rho
268 960 REAL, DIMENSION (klon) :: av_dtdwn, av_dqdwn
269
270 960 REAL, DIMENSION (klon, klev) :: rho, rhow
271 960 REAL, DIMENSION (klon, klev+1) :: rhoh
272 960 REAL, DIMENSION (klon, klev) :: rhow_moyen
273 960 REAL, DIMENSION (klon, klev) :: zh
274 960 REAL, DIMENSION (klon, klev+1) :: zhh
275 960 REAL, DIMENSION (klon, klev) :: epaisseur1, epaisseur2
276
277 960 REAL, DIMENSION (klon, klev) :: the, thu
278
279 960 REAL, DIMENSION (klon, klev) :: omgbw
280 960 REAL, DIMENSION (klon) :: pupper
281 960 REAL, DIMENSION (klon) :: omgtop
282 960 REAL, DIMENSION (klon, klev) :: dp_omgbw
283 960 REAL, DIMENSION (klon) :: ztop, dztop
284 960 REAL, DIMENSION (klon, klev) :: alpha_up
285
286 960 REAL, DIMENSION (klon) :: rre1, rre2
287 REAL :: rrd1, rrd2
288 960 REAL, DIMENSION (klon, klev) :: th1, th2, q1, q2
289 960 REAL, DIMENSION (klon, klev) :: d_th1, d_th2, d_dth
290 960 REAL, DIMENSION (klon, klev) :: d_q1, d_q2, d_dq
291 960 REAL, DIMENSION (klon, klev) :: omgbdq
292
293 960 REAL, DIMENSION (klon) :: ff, gg
294 960 REAL, DIMENSION (klon) :: wape2, cstar2, heff
295
296 960 REAL, DIMENSION (klon, klev) :: crep
297
298 960 REAL, DIMENSION (klon, klev) :: ppi
299
300 ! cc nrlmd
301 960 REAL, DIMENSION (klon) :: death_rate
302 !! REAL, DIMENSION (klon) :: nat_rate
303 960 REAL, DIMENSION (klon, klev) :: entr
304 960 REAL, DIMENSION (klon, klev) :: detr
305
306 960 REAL, DIMENSION(klon) :: sigmaw_in ! pour les prints
307 480 REAL, DIMENSION(klon) :: awdens_in, wdens_in ! pour les prints
308
309 ! -------------------------------------------------------------------------
310 ! Initialisations
311 ! -------------------------------------------------------------------------
312
313 ! print*, 'wake initialisations'
314
315 ! Essais d'initialisation avec sigmaw = 0.02 et hw = 10.
316 ! -------------------------------------------------------------------------
317
318 !! DATA wapecut, sigmad, hwmin/5., .02, 10./
319 !! DATA wapecut, sigmad, hwmin/1., .02, 10./
320 DATA sigmad, hwmin/.02, 10./
321 !! DATA wdensmin/1.e-12/
322 DATA wdensmin/1.e-14/
323 ! cc nrlmd
324 DATA sigmaw_max/0.4/
325 DATA dens_rate/0.1/
326 ! cc
327 ! Longueur de maille (en m)
328 ! -------------------------------------------------------------------------
329
330 ! ALON = 3.e5
331 ! alon = 1.E6
332
333 ! Provisionnal; to be suppressed when f_shear is parameterized
334
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 f_shear(:) = 1. ! 0. for strong shear, 1. for weak shear
335
336
337 ! Configuration de coefgw,stark,wdens (22/02/06 by YU Jingmei)
338
339 ! coefgw : Coefficient pour les ondes de gravit�
340 ! stark : Coefficient k dans Cstar=k*sqrt(2*WAPE)
341 ! wdens : Densit� surfacique de poche froide
342 ! -------------------------------------------------------------------------
343
344 ! cc nrlmd coefgw=10
345 ! coefgw=1
346 ! wdens0 = 1.0/(alon**2)
347 ! cc nrlmd wdens = 1.0/(alon**2)
348 ! cc nrlmd stark = 0.50
349 ! CRtest
350 ! cc nrlmd alpk=0.1
351 ! alpk = 1.0
352 ! alpk = 0.5
353 ! alpk = 0.05
354
355
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 479 times.
480 if (first) then
356
357 1 igout = klon/2+1/klon
358
359 1 crep_upper = 0.9
360 1 crep_sol = 1.0
361
362 ! Get wapecut from parameter file
363 1 wapecut = 1.
364 1 CALL getin_p('wapecut', wapecut)
365
366 ! cc nrlmd Lecture du fichier wake_param.data
367 1 stark=0.33
368 1 CALL getin_p('stark',stark)
369 1 cstart = stark*sqrt(2.*wapecut)
370
371 1 alpk=0.25
372 1 CALL getin_p('alpk',alpk)
373 !jyg<
374 !! wdens_ref=8.E-12
375 !! CALL getin_p('wdens_ref',wdens_ref)
376 1 wdens_ref(1)=8.E-12
377 1 wdens_ref(2)=8.E-12
378 1 CALL getin_p('wdens_ref_o',wdens_ref(1)) !wake number per unit area ; ocean
379 1 CALL getin_p('wdens_ref_l',wdens_ref(2)) !wake number per unit area ; land
380 !>jyg
381 !
382 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
383 !!!!!!!!! Population dynamics parameters !!!!!!!!!!!!!!!!!!!!!!!!!!!!
384 !------------------------------------------------------------------------
385
386 1 iflag_wk_pop_dyn = 0
387 1 CALL getin_p('iflag_wk_pop_dyn',iflag_wk_pop_dyn) ! switch between wdens prescribed
388 ! and wdens prognostic
389 1 iflag_wk_act = 0
390 1 CALL getin_p('iflag_wk_act',iflag_wk_act) ! 0: act(:)=0.
391 ! 1: act(:)=1.
392 ! 2: act(:)=f(Wape)
393
394 1 rzero = 5000.
395 1 CALL getin_p('rzero_wk', rzero)
396 1 aa0 = 3.14*rzero*rzero
397 !
398 1 tau_cv = 4000.
399 1 CALL getin_p('tau_cv', tau_cv)
400
401 !------------------------------------------------------------------------
402
403 1 coefgw=4.
404 1 CALL getin_p('coefgw',coefgw)
405
406 1 WRITE(*,*) 'stark=', stark
407 1 WRITE(*,*) 'alpk=', alpk
408 !jyg<
409 !! WRITE(*,*) 'wdens_ref=', wdens_ref
410 1 WRITE(*,*) 'wdens_ref_o=', wdens_ref(1)
411 1 WRITE(*,*) 'wdens_ref_l=', wdens_ref(2)
412 !>jyg
413 1 WRITE(*,*) 'iflag_wk_pop_dyn=',iflag_wk_pop_dyn
414 1 WRITE(*,*) 'iflag_wk_act',iflag_wk_act
415 1 WRITE(*,*) 'coefgw=', coefgw
416
417 1 flag_wk_check_trgl=.false.
418 1 CALL getin_p('flag_wk_check_trgl ', flag_wk_check_trgl)
419 1 WRITE(*,*) 'flag_wk_check_trgl=', flag_wk_check_trgl
420 1 WRITE(*,*) 'flag_wk_check_trgl OBSOLETE. Utilisr iflag_wk_check_trgl plutot'
421
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 1 times.
1 iflag_wk_check_trgl=0 ; IF (flag_wk_check_trgl) iflag_wk_check_trgl=1
422 1 CALL getin_p('iflag_wk_check_trgl ', iflag_wk_check_trgl)
423 1 WRITE(*,*) 'iflag_wk_check_trgl=', iflag_wk_check_trgl
424
425 1 first=.false.
426 endif
427
428
1/2
✓ Branch 0 taken 480 times.
✗ Branch 1 not taken.
480 IF (iflag_wk_pop_dyn == 0) THEN
429 ! Initialisation de toutes des densites a wdens_ref.
430 ! Les densites peuvent evoluer si les poches debordent
431 ! (voir au tout debut de la boucle sur les substeps)
432 !jyg<
433 !! wdens(:) = wdens_ref
434
2/2
✓ Branch 0 taken 480 times.
✓ Branch 1 taken 477120 times.
477600 DO i = 1,klon
435 477600 wdens(i) = wdens_ref(znatsurf(i)+1)
436 ENDDO
437 !>jyg
438 ENDIF ! (iflag_wk_pop_dyn == 0)
439
440 ! print*,'stark',stark
441 ! print*,'alpk',alpk
442 ! print*,'wdens',wdens
443 ! print*,'coefgw',coefgw
444 ! cc
445 ! Minimum value for |T_wake - T_undist|. Used for wake top definition
446 ! -------------------------------------------------------------------------
447
448 delta_t_min = 0.2
449
450 ! 1. - Save initial values, initialize tendencies, initialize output fields
451 ! ------------------------------------------------------------------------
452
453 !jyg<
454 !! DO k = 1, klev
455 !! DO i = 1, klon
456 !! ppi(i, k) = pi(i, k)
457 !! deltatw0(i, k) = deltatw(i, k)
458 !! deltaqw0(i, k) = deltaqw(i, k)
459 !! te(i, k) = te0(i, k)
460 !! qe(i, k) = qe0(i, k)
461 !! dtls(i, k) = 0.
462 !! dqls(i, k) = 0.
463 !! d_deltat_gw(i, k) = 0.
464 !! d_te(i, k) = 0.
465 !! d_qe(i, k) = 0.
466 !! d_deltatw(i, k) = 0.
467 !! d_deltaqw(i, k) = 0.
468 !! ! IM 060508 beg
469 !! d_deltatw2(i, k) = 0.
470 !! d_deltaqw2(i, k) = 0.
471 !! ! IM 060508 end
472 !! END DO
473 !! END DO
474
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 ppi(:,:) = pi(:,:)
475
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 deltatw0(:,:) = deltatw(:,:)
476
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 deltaqw0(:,:) = deltaqw(:,:)
477
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 te(:,:) = te0(:,:)
478
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 qe(:,:) = qe0(:,:)
479
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 dtls(:,:) = 0.
480
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 dqls(:,:) = 0.
481
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 d_deltat_gw(:,:) = 0.
482
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 d_te(:,:) = 0.
483
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 d_qe(:,:) = 0.
484
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 d_deltatw(:,:) = 0.
485
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 d_deltaqw(:,:) = 0.
486
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 d_deltatw2(:,:) = 0.
487
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 d_deltaqw2(:,:) = 0.
488
489
1/2
✓ Branch 0 taken 480 times.
✗ Branch 1 not taken.
480 IF (iflag_wk_act == 0) THEN
490
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 act(:) = 0.
491 ELSEIF (iflag_wk_act == 1) THEN
492 act(:) = 1.
493 ENDIF
494
495 !! DO i = 1, klon
496 !! sigmaw_in(i) = sigmaw(i)
497 !! END DO
498
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 sigmaw_in(:) = sigmaw(:)
499 !>jyg
500
501 ! sigmaw1=sigmaw
502 ! IF (sigd_con.GT.sigmaw1) THEN
503 ! print*, 'sigmaw,sigd_con', sigmaw, sigd_con
504 ! ENDIF
505
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
480 IF (iflag_wk_pop_dyn >=1) THEN
506 DO i = 1, klon
507 wdens_targ = max(wdens(i),wdensmin)
508 d_wdens2(i) = wdens_targ - wdens(i)
509 wdens(i) = wdens_targ
510 END DO
511 ELSE
512
2/2
✓ Branch 0 taken 480 times.
✓ Branch 1 taken 477120 times.
477600 DO i = 1, klon
513 477120 d_awdens2(i) = 0.
514 477600 d_wdens2(i) = 0.
515 END DO
516 ENDIF ! (iflag_wk_pop_dyn >=1)
517 !
518
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
519 ! c sigmaw(i) = amax1(sigmaw(i),sigd_con(i))
520 !jyg<
521 !! sigmaw(i) = amax1(sigmaw(i), sigmad)
522 !! sigmaw(i) = amin1(sigmaw(i), 0.99)
523 477120 sigmaw_targ = min(max(sigmaw(i), sigmad),0.99)
524 477120 d_sigmaw2(i) = sigmaw_targ - sigmaw(i)
525 477600 sigmaw(i) = sigmaw_targ
526 !>jyg
527 END DO
528
529 !
530
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
480 IF (iflag_wk_pop_dyn >= 1) THEN
531 awdens_in(:) = awdens(:)
532 wdens_in(:) = wdens(:)
533 !! wdens(:) = wdens(:) + wgen(:)*dtime
534 !! d_wdens2(:) = wgen(:)*dtime
535 !! ELSE
536 ENDIF ! (iflag_wk_pop_dyn >= 1)
537
538
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 wape(:) = 0.
539
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 wape2(:) = 0.
540
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 d_sigmaw(:) = 0.
541
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 ktopw(:) = 0
542 !
543 !<jyg
544
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 dth(:,:) = 0.
545
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 tu(:,:) = 0.
546
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 qu(:,:) = 0.
547
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 dtke(:,:) = 0.
548
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 dqke(:,:) = 0.
549
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 spread(:,:) = 0.
550
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 omgbdth(:,:) = 0.
551
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 omg(:,:) = 0.
552
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 dp_omgb(:,:) = 0.
553
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 dp_deltomg(:,:) = 0.
554
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 hw(:) = 0.
555
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 wape(:) = 0.
556
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 fip(:) = 0.
557
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 gfl(:) = 0.
558
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 cstar(:) = 0.
559
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 ktopw(:) = 0
560 !
561 ! Vertical advection local variables
562
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 omgbw(:,:) = 0.
563
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 omgtop(:) = 0
564
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 dp_omgbw(:,:) = 0.
565
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 18607680 times.
✓ Branch 3 taken 18720 times.
18626880 omgbdq(:,:) = 0.
566 !>jyg
567 !
568
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
480 IF (prt_level>=10) THEN
569 PRINT *, 'wake-1, sigmaw(igout) ', sigmaw(igout)
570 PRINT *, 'wake-1, deltatw(igout,k) ', (k,deltatw(igout,k), k=1,klev)
571 PRINT *, 'wake-1, deltaqw(igout,k) ', (k,deltaqw(igout,k), k=1,klev)
572 PRINT *, 'wake-1, dowwdraughts, amdwn(igout,k) ', (k,amdwn(igout,k), k=1,klev)
573 PRINT *, 'wake-1, dowwdraughts, dtdwn(igout,k) ', (k,dtdwn(igout,k), k=1,klev)
574 PRINT *, 'wake-1, dowwdraughts, dqdwn(igout,k) ', (k,dqdwn(igout,k), k=1,klev)
575 PRINT *, 'wake-1, updraughts, amup(igout,k) ', (k,amup(igout,k), k=1,klev)
576 PRINT *, 'wake-1, updraughts, dta(igout,k) ', (k,dta(igout,k), k=1,klev)
577 PRINT *, 'wake-1, updraughts, dqa(igout,k) ', (k,dqa(igout,k), k=1,klev)
578 ENDIF
579
580 ! 2. - Prognostic part
581 ! --------------------
582
583
584 ! 2.1 - Undisturbed area and Wake integrals
585 ! ---------------------------------------------------------
586
587
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
588 477120 z(i) = 0.
589 477120 ktop(i) = 0
590 477120 kupper(i) = 0
591 477120 sum_thu(i) = 0.
592 477120 sum_tu(i) = 0.
593 477120 sum_qu(i) = 0.
594 477120 sum_thvu(i) = 0.
595 477120 sum_dth(i) = 0.
596 477120 sum_dq(i) = 0.
597 477120 sum_rho(i) = 0.
598 477120 sum_dtdwn(i) = 0.
599 477120 sum_dqdwn(i) = 0.
600
601 477120 av_thu(i) = 0.
602 477120 av_tu(i) = 0.
603 477120 av_qu(i) = 0.
604 477120 av_thvu(i) = 0.
605 477120 av_dth(i) = 0.
606 477120 av_dq(i) = 0.
607 477120 av_rho(i) = 0.
608 477120 av_dtdwn(i) = 0.
609 477600 av_dqdwn(i) = 0.
610 END DO
611
612 ! Distance between wakes
613
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
614 477600 ll(i) = (1-sqrt(sigmaw(i)))/sqrt(wdens(i))
615 END DO
616 ! Potential temperatures and humidity
617 ! ----------------------------------------------------------
618
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = 1, klev
619
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
620 ! write(*,*)'wake 1',i,k,rd,te(i,k)
621 18607680 rho(i, k) = p(i, k)/(rd*te(i,k))
622 ! write(*,*)'wake 2',rho(i,k)
623
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 18130560 times.
18607680 IF (k==1) THEN
624 ! write(*,*)'wake 3',i,k,rd,te(i,k)
625 477120 rhoh(i, k) = ph(i, k)/(rd*te(i,k))
626 ! write(*,*)'wake 4',i,k,rd,te(i,k)
627 477120 zhh(i, k) = 0
628 ELSE
629 ! write(*,*)'wake 5',rd,(te(i,k)+te(i,k-1))
630 18130560 rhoh(i, k) = ph(i, k)*2./(rd*(te(i,k)+te(i,k-1)))
631 ! write(*,*)'wake 6',(-rhoh(i,k)*RG)+zhh(i,k-1)
632 18130560 zhh(i, k) = (ph(i,k)-ph(i,k-1))/(-rhoh(i,k)*rg) + zhh(i, k-1)
633 END IF
634 ! write(*,*)'wake 7',ppi(i,k)
635 18607680 the(i, k) = te(i, k)/ppi(i, k)
636 18607680 thu(i, k) = (te(i,k)-deltatw(i,k)*sigmaw(i))/ppi(i, k)
637 18607680 tu(i, k) = te(i, k) - deltatw(i, k)*sigmaw(i)
638 18607680 qu(i, k) = qe(i, k) - deltaqw(i, k)*sigmaw(i)
639 ! write(*,*)'wake 8',(rd*(te(i,k)+deltatw(i,k)))
640 18607680 rhow(i, k) = p(i, k)/(rd*(te(i,k)+deltatw(i,k)))
641 18626400 dth(i, k) = deltatw(i, k)/ppi(i, k)
642 END DO
643 END DO
644
645
2/2
✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 480 times.
18720 DO k = 1, klev - 1
646
2/2
✓ Branch 0 taken 18130560 times.
✓ Branch 1 taken 18240 times.
18149280 DO i = 1, klon
647
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 17653440 times.
18130560 IF (k==1) THEN
648 477120 n2(i, k) = 0
649 ELSE
650 n2(i, k) = amax1(0., -rg**2/the(i,k)*rho(i,k)*(the(i,k+1)-the(i,k-1))/ &
651 17653440 (p(i,k+1)-p(i,k-1)))
652 END IF
653 18130560 zh(i, k) = (zhh(i,k)+zhh(i,k+1))/2
654
655 18130560 cgw(i, k) = sqrt(n2(i,k))*zh(i, k)
656 18148800 tgw(i, k) = coefgw*cgw(i, k)/ll(i)
657 END DO
658 END DO
659
660
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
661 477120 n2(i, klev) = 0
662 477120 zh(i, klev) = 0
663 477120 cgw(i, klev) = 0
664 477600 tgw(i, klev) = 0
665 END DO
666
667 ! Calcul de la masse volumique moyenne de la colonne (bdlmd)
668 ! -----------------------------------------------------------------
669
670
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = 1, klev
671
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
672 18607680 epaisseur1(i, k) = 0.
673 18626400 epaisseur2(i, k) = 0.
674 END DO
675 END DO
676
677
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
678 477120 epaisseur1(i, 1) = -(ph(i,2)-ph(i,1))/(rho(i,1)*rg) + 1.
679 477120 epaisseur2(i, 1) = -(ph(i,2)-ph(i,1))/(rho(i,1)*rg) + 1.
680 477600 rhow_moyen(i, 1) = rhow(i, 1)
681 END DO
682
683
2/2
✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 480 times.
18720 DO k = 2, klev
684
2/2
✓ Branch 0 taken 18130560 times.
✓ Branch 1 taken 18240 times.
18149280 DO i = 1, klon
685 18130560 epaisseur1(i, k) = -(ph(i,k+1)-ph(i,k))/(rho(i,k)*rg) + 1.
686 18130560 epaisseur2(i, k) = epaisseur2(i, k-1) + epaisseur1(i, k)
687 rhow_moyen(i, k) = (rhow_moyen(i,k-1)*epaisseur2(i,k-1)+rhow(i,k)* &
688 18148800 epaisseur1(i,k))/epaisseur2(i, k)
689 END DO
690 END DO
691
692
693 ! Choose an integration bound well above wake top
694 ! -----------------------------------------------------------------
695
696 ! Pupper = 50000. ! melting level
697 ! Pupper = 60000.
698 ! Pupper = 80000. ! essais pour case_e
699
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
700 477120 pupper(i) = 0.6*ph(i, 1)
701 477600 pupper(i) = max(pupper(i), 45000.)
702 ! cc Pupper(i) = 60000.
703 END DO
704
705
706 ! Determine Wake top pressure (Ptop) from buoyancy integral
707 ! --------------------------------------------------------
708
709 ! -1/ Pressure of the level where dth becomes less than delta_t_min.
710
711
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
712 477600 ptop_provis(i) = ph(i, 1)
713 END DO
714
2/2
✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 480 times.
18720 DO k = 2, klev
715
2/2
✓ Branch 0 taken 18130560 times.
✓ Branch 1 taken 18240 times.
18149280 DO i = 1, klon
716
717 ! IM v3JYG; ptop_provis(i).LT. ph(i,1)
718
719
6/6
✓ Branch 0 taken 17951405 times.
✓ Branch 1 taken 179155 times.
✓ Branch 2 taken 90867 times.
✓ Branch 3 taken 17860538 times.
✓ Branch 4 taken 89579 times.
✓ Branch 5 taken 1288 times.
18130560 IF (dth(i,k)>-delta_t_min .AND. dth(i,k-1)<-delta_t_min .AND. &
720 18240 ptop_provis(i)==ph(i,1)) THEN
721 ptop_provis(i) = ((dth(i,k)+delta_t_min)*p(i,k-1)- &
722 89579 (dth(i,k-1)+delta_t_min)*p(i,k))/(dth(i,k)-dth(i,k-1))
723 END IF
724 END DO
725 END DO
726
727 ! -2/ dth integral
728
729
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
730 477120 sum_dth(i) = 0.
731 477120 dthmin(i) = -delta_t_min
732 477600 z(i) = 0.
733 END DO
734
735
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = 1, klev
736
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
737 18607680 dz(i) = -(amax1(ph(i,k+1),ptop_provis(i))-ph(i,k))/(rho(i,k)*rg)
738
2/2
✓ Branch 0 taken 327175 times.
✓ Branch 1 taken 18280505 times.
18626400 IF (dz(i)>0) THEN
739 327175 z(i) = z(i) + dz(i)
740 327175 sum_dth(i) = sum_dth(i) + dth(i, k)*dz(i)
741 327175 dthmin(i) = amin1(dthmin(i), dth(i,k))
742 END IF
743 END DO
744 END DO
745
746 ! -3/ height of triangle with area= sum_dth and base = dthmin
747
748
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
749 477120 hw0(i) = 2.*sum_dth(i)/amin1(dthmin(i), -0.5)
750 477600 hw0(i) = amax1(hwmin, hw0(i))
751 END DO
752
753 ! -4/ now, get Ptop
754
755
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
756 477120 z(i) = 0.
757 477600 ptop(i) = ph(i, 1)
758 END DO
759
760
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = 1, klev
761
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
762 18607680 dz(i) = amin1(-(ph(i,k+1)-ph(i,k))/(rho(i,k)*rg), hw0(i)-z(i))
763
2/2
✓ Branch 0 taken 731172 times.
✓ Branch 1 taken 17876508 times.
18626400 IF (dz(i)>0) THEN
764 731172 z(i) = z(i) + dz(i)
765 731172 ptop(i) = ph(i, k) - rho(i, k)*rg*dz(i)
766 END IF
767 END DO
768 END DO
769
770
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
480 IF (prt_level>=10) THEN
771 PRINT *, 'wake-2, ptop_provis(igout), ptop(igout) ', ptop_provis(igout), ptop(igout)
772 ENDIF
773
774
775 ! -5/ Determination de ktop et kupper
776
777
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = klev, 1, -1
778
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
779
2/2
✓ Branch 0 taken 18353628 times.
✓ Branch 1 taken 254052 times.
18607680 IF (ph(i,k+1)<ptop(i)) ktop(i) = k
780
2/2
✓ Branch 0 taken 13829755 times.
✓ Branch 1 taken 4777925 times.
18626400 IF (ph(i,k+1)<pupper(i)) kupper(i) = k
781 END DO
782 END DO
783
784 ! On evite kupper = 1 et kupper = klev
785
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
786 477120 kupper(i) = max(kupper(i), 2)
787 477600 kupper(i) = min(kupper(i), klev-1)
788 END DO
789
790
791 ! -6/ Correct ktop and ptop
792
793
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
794 477600 ptop_new(i) = ptop(i)
795 END DO
796
2/2
✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 480 times.
18720 DO k = klev, 2, -1
797
2/2
✓ Branch 0 taken 18130560 times.
✓ Branch 1 taken 18240 times.
18149280 DO i = 1, klon
798 IF (k<=ktop(i) .AND. ptop_new(i)==ptop(i) .AND. &
799
8/8
✓ Branch 0 taken 254052 times.
✓ Branch 1 taken 17876508 times.
✓ Branch 2 taken 130169 times.
✓ Branch 3 taken 123883 times.
✓ Branch 4 taken 87426 times.
✓ Branch 5 taken 42743 times.
✓ Branch 6 taken 69598 times.
✓ Branch 7 taken 17828 times.
18148800 dth(i,k)>-delta_t_min .AND. dth(i,k-1)<-delta_t_min) THEN
800 ptop_new(i) = ((dth(i,k)+delta_t_min)*p(i,k-1)-(dth(i, &
801 69598 k-1)+delta_t_min)*p(i,k))/(dth(i,k)-dth(i,k-1))
802 END IF
803 END DO
804 END DO
805
806
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
807 477600 ptop(i) = ptop_new(i)
808 END DO
809
810
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = klev, 1, -1
811
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
812
2/2
✓ Branch 0 taken 18392088 times.
✓ Branch 1 taken 215592 times.
18626400 IF (ph(i,k+1)<ptop(i)) ktop(i) = k
813 END DO
814 END DO
815
816
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
480 IF (prt_level>=10) THEN
817 PRINT *, 'wake-3, ktop(igout), kupper(igout) ', ktop(igout), kupper(igout)
818 ENDIF
819
820 ! -5/ Set deltatw & deltaqw to 0 above kupper
821
822
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = 1, klev
823
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
824
2/2
✓ Branch 0 taken 13829755 times.
✓ Branch 1 taken 4777925 times.
18626400 IF (k>=kupper(i)) THEN
825 13829755 deltatw(i, k) = 0.
826 13829755 deltaqw(i, k) = 0.
827 13829755 d_deltatw2(i,k) = -deltatw0(i,k)
828 13829755 d_deltaqw2(i,k) = -deltaqw0(i,k)
829 END IF
830 END DO
831 END DO
832
833
834 ! Vertical gradient of LS omega
835
836
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = 1, klev
837
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
838
2/2
✓ Branch 0 taken 5255045 times.
✓ Branch 1 taken 13352635 times.
18626400 IF (k<=kupper(i)) THEN
839 5255045 dp_omgb(i, k) = (omgb(i,k+1)-omgb(i,k))/(ph(i,k+1)-ph(i,k))
840 END IF
841 END DO
842 END DO
843
844 ! Integrals (and wake top level number)
845 ! --------------------------------------
846
847 ! Initialize sum_thvu to 1st level virt. pot. temp.
848
849
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
850 477120 z(i) = 1.
851 477120 dz(i) = 1.
852 477120 sum_thvu(i) = thu(i, 1)*(1.+epsim1*qu(i,1))*dz(i)
853 477600 sum_dth(i) = 0.
854 END DO
855
856
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = 1, klev
857
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
858 18607680 dz(i) = -(amax1(ph(i,k+1),ptop(i))-ph(i,k))/(rho(i,k)*rg)
859
2/2
✓ Branch 0 taken 692712 times.
✓ Branch 1 taken 17914968 times.
18626400 IF (dz(i)>0) THEN
860 692712 z(i) = z(i) + dz(i)
861 692712 sum_thu(i) = sum_thu(i) + thu(i, k)*dz(i)
862 692712 sum_tu(i) = sum_tu(i) + tu(i, k)*dz(i)
863 692712 sum_qu(i) = sum_qu(i) + qu(i, k)*dz(i)
864 692712 sum_thvu(i) = sum_thvu(i) + thu(i, k)*(1.+epsim1*qu(i,k))*dz(i)
865 692712 sum_dth(i) = sum_dth(i) + dth(i, k)*dz(i)
866 692712 sum_dq(i) = sum_dq(i) + deltaqw(i, k)*dz(i)
867 692712 sum_rho(i) = sum_rho(i) + rhow(i, k)*dz(i)
868 692712 sum_dtdwn(i) = sum_dtdwn(i) + dtdwn(i, k)*dz(i)
869 692712 sum_dqdwn(i) = sum_dqdwn(i) + dqdwn(i, k)*dz(i)
870 END IF
871 END DO
872 END DO
873
874
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
875 477600 hw0(i) = z(i)
876 END DO
877
878
879 ! 2.1 - WAPE and mean forcing computation
880 ! ---------------------------------------
881
882 ! ---------------------------------------
883
884 ! Means
885
886
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
887 477120 av_thu(i) = sum_thu(i)/hw0(i)
888 477120 av_tu(i) = sum_tu(i)/hw0(i)
889 477120 av_qu(i) = sum_qu(i)/hw0(i)
890 477120 av_thvu(i) = sum_thvu(i)/hw0(i)
891 ! av_thve = sum_thve/hw0
892 477120 av_dth(i) = sum_dth(i)/hw0(i)
893 477120 av_dq(i) = sum_dq(i)/hw0(i)
894 477120 av_rho(i) = sum_rho(i)/hw0(i)
895 477120 av_dtdwn(i) = sum_dtdwn(i)/hw0(i)
896 477120 av_dqdwn(i) = sum_dqdwn(i)/hw0(i)
897
898 wape(i) = -rg*hw0(i)*(av_dth(i)+ &
899 477600 epsim1*(av_thu(i)*av_dq(i)+av_dth(i)*av_qu(i)+av_dth(i)*av_dq(i)))/av_thvu(i)
900
901 END DO
902
903 ! 2.2 Prognostic variable update
904 ! ------------------------------
905
906 ! Filter out bad wakes
907
908
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = 1, klev
909
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
910
2/2
✓ Branch 0 taken 6106503 times.
✓ Branch 1 taken 12501177 times.
18626400 IF (wape(i)<0.) THEN
911 6106503 deltatw(i, k) = 0.
912 6106503 deltaqw(i, k) = 0.
913 6106503 dth(i, k) = 0.
914 6106503 d_deltatw2(i,k) = -deltatw0(i,k)
915 6106503 d_deltaqw2(i,k) = -deltaqw0(i,k)
916 END IF
917 END DO
918 END DO
919
920
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
921
2/2
✓ Branch 0 taken 156577 times.
✓ Branch 1 taken 320543 times.
477600 IF (wape(i)<0.) THEN
922 156577 wape(i) = 0.
923 156577 cstar(i) = 0.
924 156577 hw(i) = hwmin
925 !jyg<
926 !! sigmaw(i) = amax1(sigmad, sigd_con(i))
927 156577 sigmaw_targ = max(sigmad, sigd_con(i))
928 156577 d_sigmaw2(i) = d_sigmaw2(i) + sigmaw_targ - sigmaw(i)
929 156577 sigmaw(i) = sigmaw_targ
930 !>jyg
931 156577 fip(i) = 0.
932 156577 gwake(i) = .FALSE.
933 ELSE
934 320543 hw(i) = hw0(i)
935 320543 cstar(i) = stark*sqrt(2.*wape(i))
936 320543 gwake(i) = .TRUE.
937 END IF
938 END DO
939
940
941 ! Check qx and qw positivity
942 ! --------------------------
943
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
944 q0_min(i) = min((qe(i,1)-sigmaw(i)*deltaqw(i,1)), &
945 477600 (qe(i,1)+(1.-sigmaw(i))*deltaqw(i,1)))
946 END DO
947
2/2
✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 480 times.
18720 DO k = 2, klev
948
2/2
✓ Branch 0 taken 18130560 times.
✓ Branch 1 taken 18240 times.
18149280 DO i = 1, klon
949 q1_min(i) = min((qe(i,k)-sigmaw(i)*deltaqw(i,k)), &
950 18130560 (qe(i,k)+(1.-sigmaw(i))*deltaqw(i,k)))
951
2/2
✓ Branch 0 taken 10639991 times.
✓ Branch 1 taken 7490569 times.
18148800 IF (q1_min(i)<=q0_min(i)) THEN
952 10639991 q0_min(i) = q1_min(i)
953 END IF
954 END DO
955 END DO
956
957
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
958 477120 ok_qx_qw(i) = q0_min(i) >= 0.
959 477600 alpha(i) = 1.
960 END DO
961
962
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
480 IF (prt_level>=10) THEN
963 PRINT *, 'wake-4, sigmaw(igout), cstar(igout), wape(igout), ktop(igout) ', &
964 sigmaw(igout), cstar(igout), wape(igout), ktop(igout)
965 ENDIF
966
967
968 ! C -----------------------------------------------------------------
969 ! Sub-time-stepping
970 ! -----------------
971
972 nsub = 10
973 480 dtimesub = dtime/nsub
974
975 ! ------------------------------------------------------------
976
2/2
✓ Branch 0 taken 4800 times.
✓ Branch 1 taken 480 times.
5280 DO isubstep = 1, nsub
977 ! ------------------------------------------------------------
978
979 ! wk_adv is the logical flag enabling wake evolution in the time advance
980 ! loop
981
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
982
2/4
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 4771200 times.
4776000 wk_adv(i) = ok_qx_qw(i) .AND. alpha(i) >= 1.
983 END DO
984
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 4800 times.
4800 IF (prt_level>=10) THEN
985 PRINT *, 'wake-4.1, isubstep,wk_adv(igout),cstar(igout),wape(igout), ptop(igout) ', &
986 isubstep,wk_adv(igout),cstar(igout),wape(igout), ptop(igout)
987 ENDIF
988
989 ! cc nrlmd Ajout d'un recalcul de wdens dans le cas d'un entrainement
990 ! n�gatif de ktop � kupper --------
991 ! cc On calcule pour cela une densit� wdens0 pour laquelle on
992 ! aurait un entrainement nul ---
993 !jyg<
994 ! Dans la configuration avec wdens prognostique, il s'agit d'un cas ou
995 ! les poches sont insuffisantes pour accueillir tout le flux de masse
996 ! des descentes unsaturees. Nous faisons alors l'hypothese que la
997 ! convection profonde cree directement de nouvelles poches, sans passer
998 ! par les thermiques. La nouvelle valeur de wdens est alors impos�e.
999
1000
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1001 ! c print *,' isubstep,wk_adv(i),cstar(i),wape(i) ',
1002 ! c $ isubstep,wk_adv(i),cstar(i),wape(i)
1003
3/4
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 973300 times.
✓ Branch 3 taken 3797900 times.
4776000 IF (wk_adv(i) .AND. cstar(i)>0.01) THEN
1004 omg(i, kupper(i)+1) = -rg*amdwn(i, kupper(i)+1)/sigmaw(i) + &
1005 973300 rg*amup(i, kupper(i)+1)/(1.-sigmaw(i))
1006 wdens0 = (sigmaw(i)/(4.*3.14))* &
1007 973300 ((1.-sigmaw(i))*omg(i,kupper(i)+1)/((ph(i,1)-pupper(i))*cstar(i)))**(2)
1008
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 973300 times.
973300 IF (prt_level >= 10) THEN
1009 print*,'omg(i,kupper(i)+1),wdens0,wdens(i),cstar(i), ph(i,1)-pupper(i)', &
1010 omg(i,kupper(i)+1),wdens0,wdens(i),cstar(i), ph(i,1)-pupper(i)
1011 ENDIF
1012
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 973300 times.
973300 IF (wdens(i)<=wdens0*1.1) THEN
1013 IF (iflag_wk_pop_dyn >= 1) THEN
1014 d_wdens2(i) = d_wdens2(i) + wdens0 - wdens(i)
1015 ENDIF
1016 wdens(i) = wdens0
1017 END IF
1018 END IF
1019 END DO
1020
1021
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1022
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1023 4771200 gfl(i) = 2.*sqrt(3.14*wdens(i)*sigmaw(i))
1024 4771200 rad_wk(i) = sqrt(sigmaw(i)/(3.14*wdens(i)))
1025 !jyg<
1026 !! sigmaw(i) = amin1(sigmaw(i), sigmaw_max)
1027 4771200 sigmaw_targ = min(sigmaw(i), sigmaw_max)
1028 4771200 d_sigmaw2(i) = d_sigmaw2(i) + sigmaw_targ - sigmaw(i)
1029 4771200 sigmaw(i) = sigmaw_targ
1030 !>jyg
1031 END IF
1032 END DO
1033
1034
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 4800 times.
4800 IF (iflag_wk_pop_dyn >= 1) THEN
1035 ! The variable "death_rate" is significant only when iflag_wk_pop_dyn = 0.
1036 ! Here, it has to be set to zero.
1037 death_rate(:) = 0.
1038
1039 IF (iflag_wk_act ==2) THEN
1040 DO i = 1, klon
1041 IF (wk_adv(i)) THEN
1042 wape1_act(i) = abs(cin(i))
1043 wape2_act(i) = 2.*wape1_act(i) + 1.
1044 act(i) = min(1., max(0., (wape(i)-wape1_act(i)) / (wape2_act(i)-wape1_act(i)) ))
1045 ENDIF ! (wk_adv(i))
1046 ENDDO
1047 ENDIF ! (iflag_wk_act ==2)
1048
1049
1050 DO i = 1, klon
1051 IF (wk_adv(i)) THEN
1052 !! tau_wk(i) = max(rad_wk(i)/(3.*cstar(i))*((cstar(i)/cstart)**1.5 - 1), 100.)
1053 tau_wk_inv(i) = max( (3.*cstar(i))/(rad_wk(i)*((cstar(i)/cstart)**1.5 - 1)), 0.)
1054 tau_wk_inv_min = min(tau_wk_inv(i), 1./dtimesub)
1055 drdt(i) = (cstar(i) - wgen(i)*(sigmaw(i)/wdens(i)-aa0)/gfl(i)) / &
1056 (1 + 2*f_shear(i)*(2.*sigmaw(i)-aa0*wdens(i)) - 2.*sigmaw(i))
1057 !! (1 - 2*sigmaw(i)*(1.-f_shear(i)))
1058 drdt_pos=max(drdt(i),0.)
1059
1060 !! d_wdens(i) = ( wgen(i)*(1.+2.*(sigmaw(i)-sigmad)) &
1061 !! - wdens(i)*tau_wk_inv_min &
1062 !! - 2.*gfl(i)*wdens(i)*Cstar(i) )*dtimesub
1063 d_awdens(i) = ( wgen(i) - (1./tau_cv)*(awdens(i) - act(i)*wdens(i)) )*dtimesub
1064 d_wdens(i) = ( wgen(i) - (wdens(i)-awdens(i))*tau_wk_inv_min - &
1065 2.*wdens(i)*gfl(i)*drdt_pos )*dtimesub
1066 d_wdens(i) = max(d_wdens(i), wdensmin-wdens(i))
1067
1068 !! d_sigmaw(i) = ( (1.-2*f_shear(i)*sigmaw(i))*(gfl(i)*Cstar(i)+wgen(i)*sigmad/wdens(i)) &
1069 !! + 2.*f_shear(i)*wgen(i)*sigmaw(i)**2/wdens(i) &
1070 !! - sigmaw(i)*tau_wk_inv_min )*dtimesub
1071 d_sig_gen(i) = wgen(i)*aa0
1072 d_sig_death(i) = - sigmaw(i)*(1.-awdens(i)/wdens(i))*tau_wk_inv_min
1073 !! d_sig_col(i) = - 2*f_shear(i)*sigmaw(i)*gfl(i)*drdt_pos
1074 d_sig_col(i) = - 2*f_shear(i)*(2.*sigmaw(i)-wdens(i)*aa0)*gfl(i)*drdt_pos
1075 d_sigmaw(i) = ( d_sig_gen(i) + d_sig_death(i) + d_sig_col(i) + gfl(i)*cstar(i) )*dtimesub
1076 d_sigmaw(i) = max(d_sigmaw(i), sigmad-sigmaw(i))
1077 ENDIF
1078 ENDDO
1079
1080 IF (prt_level >= 10) THEN
1081 print *,'wake, cstar(1), cstar(1)/cstart, rad_wk(1), tau_wk_inv(1), drdt(1) ', &
1082 cstar(1), cstar(1)/cstart, rad_wk(1), tau_wk_inv(1), drdt(1)
1083 print *,'wake, wdens(1), awdens(1), act(1), d_awdens(1) ', &
1084 wdens(1), awdens(1), act(1), d_awdens(1)
1085 print *,'wake, wgen, -(wdens-awdens)*tau_wk_inv, -2.*wdens*gfl*drdt_pos, d_wdens ', &
1086 wgen(1), -(wdens(1)-awdens(1))*tau_wk_inv(1), -2.*wdens(1)*gfl(1)*drdt_pos, d_wdens(1)
1087 print *,'wake, d_sig_gen(1), d_sig_death(1), d_sig_col(1), d_sigmaw(1) ', &
1088 d_sig_gen(1), d_sig_death(1), d_sig_col(1), d_sigmaw(1)
1089 ENDIF
1090
1091 ELSE ! (iflag_wk_pop_dyn >= 1)
1092
1093 ! cc nrlmd
1094
1095
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1096
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1097 ! cc nrlmd Introduction du taux de mortalit� des poches et
1098 ! test sur sigmaw_max=0.4
1099 ! cc d_sigmaw(i) = gfl(i)*Cstar(i)*dtimesub
1100
2/2
✓ Branch 0 taken 644397 times.
✓ Branch 1 taken 4126803 times.
4771200 IF (sigmaw(i)>=sigmaw_max) THEN
1101 644397 death_rate(i) = gfl(i)*cstar(i)/sigmaw(i)
1102 ELSE
1103 4126803 death_rate(i) = 0.
1104 END IF
1105
1106 d_sigmaw(i) = gfl(i)*cstar(i)*dtimesub - death_rate(i)*sigmaw(i)* &
1107 4771200 dtimesub
1108 ! $ - nat_rate(i)*sigmaw(i)*dtimesub
1109 ! c print*, 'd_sigmaw(i),sigmaw(i),gfl(i),Cstar(i),wape(i),
1110 ! c $ death_rate(i),ktop(i),kupper(i)',
1111 ! c $ d_sigmaw(i),sigmaw(i),gfl(i),Cstar(i),wape(i),
1112 ! c $ death_rate(i),ktop(i),kupper(i)
1113
1114 ! sigmaw(i) =sigmaw(i) + gfl(i)*Cstar(i)*dtimesub
1115 ! sigmaw(i) =min(sigmaw(i),0.99) !!!!!!!!
1116 ! wdens = wdens0/(10.*sigmaw)
1117 ! sigmaw =max(sigmaw,sigd_con)
1118 ! sigmaw =max(sigmaw,sigmad)
1119 END IF
1120 END DO
1121
1122 ENDIF ! (iflag_wk_pop_dyn >= 1)
1123
1124
1125 ! calcul de la difference de vitesse verticale poche - zone non perturbee
1126 ! IM 060208 differences par rapport au code initial; init. a 0 dp_deltomg
1127 ! IM 060208 et omg sur les niveaux de 1 a klev+1, alors que avant l'on definit
1128 ! IM 060208 au niveau k=1..?
1129 !JYG 161013 Correction : maintenant omg est dimensionne a klev.
1130
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1131
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1132
1/2
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
186264000 IF (wk_adv(i)) THEN !!! nrlmd
1133 186076800 dp_deltomg(i, k) = 0.
1134 END IF
1135 END DO
1136 END DO
1137
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1138
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1139
1/2
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
186264000 IF (wk_adv(i)) THEN !!! nrlmd
1140 186076800 omg(i, k) = 0.
1141 END IF
1142 END DO
1143 END DO
1144
1145
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1146
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1147 4771200 z(i) = 0.
1148 4771200 omg(i, 1) = 0.
1149 4771200 dp_deltomg(i, 1) = -(gfl(i)*cstar(i))/(sigmaw(i)*(1-sigmaw(i)))
1150 END IF
1151 END DO
1152
1153
2/2
✓ Branch 0 taken 182400 times.
✓ Branch 1 taken 4800 times.
187200 DO k = 2, klev
1154
2/2
✓ Branch 0 taken 181305600 times.
✓ Branch 1 taken 182400 times.
181492800 DO i = 1, klon
1155
3/4
✓ Branch 0 taken 181305600 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2021926 times.
✓ Branch 3 taken 179283674 times.
181488000 IF (wk_adv(i) .AND. k<=ktop(i)) THEN
1156 2021926 dz(i) = -(ph(i,k)-ph(i,k-1))/(rho(i,k-1)*rg)
1157 2021926 z(i) = z(i) + dz(i)
1158 2021926 dp_deltomg(i, k) = dp_deltomg(i, 1)
1159 2021926 omg(i, k) = dp_deltomg(i, 1)*z(i)
1160 END IF
1161 END DO
1162 END DO
1163
1164
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1165
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1166 4771200 dztop(i) = -(ptop(i)-ph(i,ktop(i)))/(rho(i,ktop(i))*rg)
1167 4771200 ztop(i) = z(i) + dztop(i)
1168 4771200 omgtop(i) = dp_deltomg(i, 1)*ztop(i)
1169 END IF
1170 END DO
1171
1172
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 4800 times.
4800 IF (prt_level>=10) THEN
1173 PRINT *, 'wake-4.2, omg(igout,k) ', (k,omg(igout,k), k=1,klev)
1174 PRINT *, 'wake-4.2, omgtop(igout), ptop(igout), ktop(igout) ', &
1175 omgtop(igout), ptop(igout), ktop(igout)
1176 ENDIF
1177
1178 ! -----------------
1179 ! From m/s to Pa/s
1180 ! -----------------
1181
1182
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1183
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1184 4771200 omgtop(i) = -rho(i, ktop(i))*rg*omgtop(i)
1185 4771200 dp_deltomg(i, 1) = omgtop(i)/(ptop(i)-ph(i,1))
1186 END IF
1187 END DO
1188
1189
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1190
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1191
3/4
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 6793126 times.
✓ Branch 3 taken 179283674 times.
186264000 IF (wk_adv(i) .AND. k<=ktop(i)) THEN
1192 6793126 omg(i, k) = -rho(i, k)*rg*omg(i, k)
1193 6793126 dp_deltomg(i, k) = dp_deltomg(i, 1)
1194 END IF
1195 END DO
1196 END DO
1197
1198 ! raccordement lineaire de omg de ptop a pupper
1199
1200
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1201
3/4
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 4771197 times.
✓ Branch 3 taken 3 times.
4776000 IF (wk_adv(i) .AND. kupper(i)>ktop(i)) THEN
1202 omg(i, kupper(i)+1) = -rg*amdwn(i, kupper(i)+1)/sigmaw(i) + &
1203 4771197 rg*amup(i, kupper(i)+1)/(1.-sigmaw(i))
1204 dp_deltomg(i, kupper(i)) = (omgtop(i)-omg(i,kupper(i)+1))/ &
1205 4771197 (ptop(i)-pupper(i))
1206 END IF
1207 END DO
1208
1209 ! c DO i=1,klon
1210 ! c print*,'Pente entre 0 et kupper (r�f�rence)'
1211 ! c $ ,omg(i,kupper(i)+1)/(pupper(i)-ph(i,1))
1212 ! c print*,'Pente entre ktop et kupper'
1213 ! c $ ,(omg(i,kupper(i)+1)-omgtop(i))/(pupper(i)-ptop(i))
1214 ! c ENDDO
1215 ! c
1216
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1217
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1218
5/6
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 179283674 times.
✓ Branch 3 taken 6793126 times.
✓ Branch 4 taken 45757324 times.
✓ Branch 5 taken 133526350 times.
186264000 IF (wk_adv(i) .AND. k>ktop(i) .AND. k<=kupper(i)) THEN
1219 45757324 dp_deltomg(i, k) = dp_deltomg(i, kupper(i))
1220 45757324 omg(i, k) = omgtop(i) + (ph(i,k)-ptop(i))*dp_deltomg(i, kupper(i))
1221 END IF
1222 END DO
1223 END DO
1224 !! print *,'omg(igout,k) ', (k,omg(igout,k),k=1,klev)
1225 ! cc nrlmd
1226 ! c DO i=1,klon
1227 ! c print*,'deltaw_ktop,deltaw_conv',omgtop(i),omg(i,kupper(i)+1)
1228 ! c END DO
1229 ! cc
1230
1231
1232 ! -- Compute wake average vertical velocity omgbw
1233
1234
1235
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1236
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1237
1/2
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
186264000 IF (wk_adv(i)) THEN
1238 186076800 omgbw(i, k) = omgb(i, k) + (1.-sigmaw(i))*omg(i, k)
1239 END IF
1240 END DO
1241 END DO
1242 ! -- and its vertical gradient dp_omgbw
1243
1244
2/2
✓ Branch 0 taken 182400 times.
✓ Branch 1 taken 4800 times.
187200 DO k = 1, klev-1
1245
2/2
✓ Branch 0 taken 181305600 times.
✓ Branch 1 taken 182400 times.
181492800 DO i = 1, klon
1246
1/2
✓ Branch 0 taken 181305600 times.
✗ Branch 1 not taken.
181488000 IF (wk_adv(i)) THEN
1247 181305600 dp_omgbw(i, k) = (omgbw(i,k+1)-omgbw(i,k))/(ph(i,k+1)-ph(i,k))
1248 END IF
1249 END DO
1250 END DO
1251
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1252
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1253 4771200 dp_omgbw(i, klev) = 0.
1254 END IF
1255 END DO
1256
1257 ! -- Upstream coefficients for omgb velocity
1258 ! -- (alpha_up(k) is the coefficient of the value at level k)
1259 ! -- (1-alpha_up(k) is the coefficient of the value at level k-1)
1260
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1261
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1262
1/2
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
186264000 IF (wk_adv(i)) THEN
1263 186076800 alpha_up(i, k) = 0.
1264
2/2
✓ Branch 0 taken 90720030 times.
✓ Branch 1 taken 95356770 times.
186076800 IF (omgb(i,k)>0.) alpha_up(i, k) = 1.
1265 END IF
1266 END DO
1267 END DO
1268
1269 ! Matrix expressing [The,deltatw] from [Th1,Th2]
1270
1271
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1272
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1273 4771200 rre1(i) = 1. - sigmaw(i)
1274 4771200 rre2(i) = sigmaw(i)
1275 END IF
1276 END DO
1277 rrd1 = -1.
1278 rrd2 = 1.
1279
1280 ! -- Get [Th1,Th2], dth and [q1,q2]
1281
1282
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1283
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1284
3/4
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 57321650 times.
✓ Branch 3 taken 128755150 times.
186264000 IF (wk_adv(i) .AND. k<=kupper(i)+1) THEN
1285 57321650 dth(i, k) = deltatw(i, k)/ppi(i, k)
1286 57321650 th1(i, k) = the(i, k) - sigmaw(i)*dth(i, k) ! undisturbed area
1287 57321650 th2(i, k) = the(i, k) + (1.-sigmaw(i))*dth(i, k) ! wake
1288 57321650 q1(i, k) = qe(i, k) - sigmaw(i)*deltaqw(i, k) ! undisturbed area
1289 57321650 q2(i, k) = qe(i, k) + (1.-sigmaw(i))*deltaqw(i, k) ! wake
1290 END IF
1291 END DO
1292 END DO
1293
1294
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1295
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN !!! nrlmd
1296 4771200 d_th1(i, 1) = 0.
1297 4771200 d_th2(i, 1) = 0.
1298 4771200 d_dth(i, 1) = 0.
1299 4771200 d_q1(i, 1) = 0.
1300 4771200 d_q2(i, 1) = 0.
1301 4771200 d_dq(i, 1) = 0.
1302 END IF
1303 END DO
1304
1305
2/2
✓ Branch 0 taken 182400 times.
✓ Branch 1 taken 4800 times.
187200 DO k = 2, klev
1306
2/2
✓ Branch 0 taken 181305600 times.
✓ Branch 1 taken 182400 times.
181492800 DO i = 1, klon
1307
3/4
✓ Branch 0 taken 181305600 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 52550450 times.
✓ Branch 3 taken 128755150 times.
181488000 IF (wk_adv(i) .AND. k<=kupper(i)+1) THEN
1308 52550450 d_th1(i, k) = th1(i, k-1) - th1(i, k)
1309 52550450 d_th2(i, k) = th2(i, k-1) - th2(i, k)
1310 52550450 d_dth(i, k) = dth(i, k-1) - dth(i, k)
1311 52550450 d_q1(i, k) = q1(i, k-1) - q1(i, k)
1312 52550450 d_q2(i, k) = q2(i, k-1) - q2(i, k)
1313 52550450 d_dq(i, k) = deltaqw(i, k-1) - deltaqw(i, k)
1314 END IF
1315 END DO
1316 END DO
1317
1318
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1319
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1320 4771200 omgbdth(i, 1) = 0.
1321 4771200 omgbdq(i, 1) = 0.
1322 END IF
1323 END DO
1324
1325
2/2
✓ Branch 0 taken 182400 times.
✓ Branch 1 taken 4800 times.
187200 DO k = 2, klev
1326
2/2
✓ Branch 0 taken 181305600 times.
✓ Branch 1 taken 182400 times.
181492800 DO i = 1, klon
1327
3/4
✓ Branch 0 taken 181305600 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 52550450 times.
✓ Branch 3 taken 128755150 times.
181488000 IF (wk_adv(i) .AND. k<=kupper(i)+1) THEN ! loop on interfaces
1328 52550450 omgbdth(i, k) = omgb(i, k)*(dth(i,k-1)-dth(i,k))
1329 52550450 omgbdq(i, k) = omgb(i, k)*(deltaqw(i,k-1)-deltaqw(i,k))
1330 END IF
1331 END DO
1332 END DO
1333
1334
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 4800 times.
4800 IF (prt_level>=10) THEN
1335 PRINT *, 'wake-4.3, th1(igout,k) ', (k,th1(igout,k), k=1,klev)
1336 PRINT *, 'wake-4.3, th2(igout,k) ', (k,th2(igout,k), k=1,klev)
1337 PRINT *, 'wake-4.3, dth(igout,k) ', (k,dth(igout,k), k=1,klev)
1338 PRINT *, 'wake-4.3, omgbdth(igout,k) ', (k,omgbdth(igout,k), k=1,klev)
1339 ENDIF
1340
1341 ! -----------------------------------------------------------------
1342
2/2
✓ Branch 0 taken 182400 times.
✓ Branch 1 taken 4800 times.
187200 DO k = 1, klev-1
1343
2/2
✓ Branch 0 taken 181305600 times.
✓ Branch 1 taken 182400 times.
181492800 DO i = 1, klon
1344
3/4
✓ Branch 0 taken 181305600 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 47779250 times.
✓ Branch 3 taken 133526350 times.
181488000 IF (wk_adv(i) .AND. k<=kupper(i)-1) THEN
1345 ! -----------------------------------------------------------------
1346
1347 ! Compute redistribution (advective) term
1348
1349 d_deltatw(i, k) = dtimesub/(ph(i,k)-ph(i,k+1))* &
1350 (rrd1*omg(i,k)*sigmaw(i)*d_th1(i,k) - &
1351 rrd2*omg(i,k+1)*(1.-sigmaw(i))*d_th2(i,k+1)- &
1352 (1.-alpha_up(i,k))*omgbdth(i,k)- &
1353 47779250 alpha_up(i,k+1)*omgbdth(i,k+1))*ppi(i, k)
1354 ! print*,'d_deltatw=', k, d_deltatw(i,k)
1355
1356 d_deltaqw(i, k) = dtimesub/(ph(i,k)-ph(i,k+1))* &
1357 (rrd1*omg(i,k)*sigmaw(i)*d_q1(i,k)- &
1358 rrd2*omg(i,k+1)*(1.-sigmaw(i))*d_q2(i,k+1)- &
1359 (1.-alpha_up(i,k))*omgbdq(i,k)- &
1360 47779250 alpha_up(i,k+1)*omgbdq(i,k+1))
1361 ! print*,'d_deltaqw=', k, d_deltaqw(i,k)
1362
1363 ! and increment large scale tendencies
1364
1365
1366
1367
1368 ! C
1369 ! -----------------------------------------------------------------
1370 d_te(i, k) = dtimesub*((rre1(i)*omg(i,k)*sigmaw(i)*d_th1(i,k)- &
1371 rre2(i)*omg(i,k+1)*(1.-sigmaw(i))*d_th2(i,k+1))/ &
1372 (ph(i,k)-ph(i,k+1)) &
1373 -sigmaw(i)*(1.-sigmaw(i))*dth(i,k)*(omg(i,k)-omg(i,k+1))/ &
1374 47779250 (ph(i,k)-ph(i,k+1)) )*ppi(i, k)
1375
1376 d_qe(i, k) = dtimesub*((rre1(i)*omg(i,k)*sigmaw(i)*d_q1(i,k)- &
1377 rre2(i)*omg(i,k+1)*(1.-sigmaw(i))*d_q2(i,k+1))/ &
1378 (ph(i,k)-ph(i,k+1)) &
1379 -sigmaw(i)*(1.-sigmaw(i))*deltaqw(i,k)*(omg(i,k)-omg(i,k+1))/ &
1380 47779250 (ph(i,k)-ph(i,k+1)) )
1381
3/4
✓ Branch 0 taken 133526350 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 4771200 times.
✓ Branch 3 taken 128755150 times.
133526350 ELSE IF (wk_adv(i) .AND. k==kupper(i)) THEN
1382 4771200 d_te(i, k) = dtimesub*(rre1(i)*omg(i,k)*sigmaw(i)*d_th1(i,k)/(ph(i,k)-ph(i,k+1)))*ppi(i, k)
1383
1384 4771200 d_qe(i, k) = dtimesub*(rre1(i)*omg(i,k)*sigmaw(i)*d_q1(i,k)/(ph(i,k)-ph(i,k+1)))
1385
1386 END IF
1387 ! cc
1388 END DO
1389 END DO
1390 ! ------------------------------------------------------------------
1391
1392
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 4800 times.
4800 IF (prt_level>=10) THEN
1393 PRINT *, 'wake-4.3, d_deltatw(igout,k) ', (k,d_deltatw(igout,k), k=1,klev)
1394 PRINT *, 'wake-4.3, d_deltaqw(igout,k) ', (k,d_deltaqw(igout,k), k=1,klev)
1395 ENDIF
1396
1397 ! Increment state variables
1398 !jyg<
1399
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 4800 times.
4800 IF (iflag_wk_pop_dyn >= 1) THEN
1400 DO k = 1, klev
1401 DO i = 1, klon
1402 IF (wk_adv(i) .AND. k<=kupper(i)) THEN
1403 detr(i,k) = - d_sig_death(i) - d_sig_col(i)
1404 entr(i,k) = d_sig_gen(i)
1405 ENDIF
1406 ENDDO
1407 ENDDO
1408 ELSE ! (iflag_wk_pop_dyn >= 1)
1409
2/2
✓ Branch 0 taken 4800 times.
✓ Branch 1 taken 187200 times.
192000 DO k = 1, klev
1410
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1411
3/4
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 52550450 times.
✓ Branch 3 taken 133526350 times.
186264000 IF (wk_adv(i) .AND. k<=kupper(i)) THEN
1412 52550450 detr(i, k) = 0.
1413
1414 52550450 entr(i, k) = 0.
1415 ENDIF
1416 ENDDO
1417 ENDDO
1418 ENDIF ! (iflag_wk_pop_dyn >= 1)
1419
1420
1421
1422
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1423
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1424 ! cc nrlmd IF( wk_adv(i) .AND. k .LE. kupper(i)-1) THEN
1425
3/4
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 52550450 times.
✓ Branch 3 taken 133526350 times.
186264000 IF (wk_adv(i) .AND. k<=kupper(i)) THEN
1426 ! cc
1427
1428
1429
1430 ! Coefficient de r�partition
1431
1432 crep(i, k) = crep_sol*(ph(i,kupper(i))-ph(i,k))/ &
1433 52550450 (ph(i,kupper(i))-ph(i,1))
1434 crep(i, k) = crep(i, k) + crep_upper*(ph(i,1)-ph(i,k))/ &
1435 52550450 (p(i,1)-ph(i,kupper(i)))
1436
1437
1438 ! Reintroduce compensating subsidence term.
1439
1440 ! dtKE(k)=(dtdwn(k)*Crep(k))/sigmaw
1441 ! dtKE(k)=dtKE(k)-(dtdwn(k)*(1-Crep(k))+dta(k))
1442 ! . /(1-sigmaw)
1443 ! dqKE(k)=(dqdwn(k)*Crep(k))/sigmaw
1444 ! dqKE(k)=dqKE(k)-(dqdwn(k)*(1-Crep(k))+dqa(k))
1445 ! . /(1-sigmaw)
1446
1447 ! dtKE(k)=(dtdwn(k)*Crep(k)+(1-Crep(k))*dta(k))/sigmaw
1448 ! dtKE(k)=dtKE(k)-(dtdwn(k)*(1-Crep(k))+dta(k)*Crep(k))
1449 ! . /(1-sigmaw)
1450 ! dqKE(k)=(dqdwn(k)*Crep(k)+(1-Crep(k))*dqa(k))/sigmaw
1451 ! dqKE(k)=dqKE(k)-(dqdwn(k)*(1-Crep(k))+dqa(k)*Crep(k))
1452 ! . /(1-sigmaw)
1453
1454 52550450 dtke(i, k) = (dtdwn(i,k)/sigmaw(i)-dta(i,k)/(1.-sigmaw(i)))
1455 52550450 dqke(i, k) = (dqdwn(i,k)/sigmaw(i)-dqa(i,k)/(1.-sigmaw(i)))
1456 ! print*,'dtKE= ',dtKE(i,k),' dqKE= ',dqKE(i,k)
1457
1458 !
1459
1460 ! cc nrlmd Prise en compte du taux de mortalit�
1461 ! cc D�finitions de entr, detr
1462 !jyg<
1463 !! detr(i, k) = 0.
1464 !!
1465 !! entr(i, k) = detr(i, k) + gfl(i)*cstar(i) + &
1466 !! sigmaw(i)*(1.-sigmaw(i))*dp_deltomg(i, k)
1467 !!
1468 entr(i, k) = entr(i,k) + gfl(i)*cstar(i) + &
1469 52550450 sigmaw(i)*(1.-sigmaw(i))*dp_deltomg(i, k)
1470 !>jyg
1471 52550450 spread(i, k) = (entr(i,k)-detr(i,k))/sigmaw(i)
1472
1473 ! cc spread(i,k) =
1474 ! (1.-sigmaw(i))*dp_deltomg(i,k)+gfl(i)*Cstar(i)/
1475 ! cc $ sigmaw(i)
1476
1477
1478 ! ajout d'un effet onde de gravit� -Tgw(k)*deltatw(k) 03/02/06 YU
1479 ! Jingmei
1480
1481 ! write(lunout,*)'wake.F ',i,k, dtimesub,d_deltat_gw(i,k),
1482 ! & Tgw(i,k),deltatw(i,k)
1483 d_deltat_gw(i, k) = d_deltat_gw(i, k) - tgw(i, k)*deltatw(i, k)* &
1484 52550450 dtimesub
1485 ! write(lunout,*)'wake.F ',i,k, dtimesub,d_deltatw(i,k)
1486 52550450 ff(i) = d_deltatw(i, k)/dtimesub
1487
1488 ! Sans GW
1489
1490 ! deltatw(k)=deltatw(k)+dtimesub*(ff+dtKE(k)-spread(k)*deltatw(k))
1491
1492 ! GW formule 1
1493
1494 ! deltatw(k) = deltatw(k)+dtimesub*
1495 ! $ (ff+dtKE(k) - spread(k)*deltatw(k)-Tgw(k)*deltatw(k))
1496
1497 ! GW formule 2
1498
1499
2/2
✓ Branch 0 taken 8072060 times.
✓ Branch 1 taken 44478390 times.
52550450 IF (dtimesub*tgw(i,k)<1.E-10) THEN
1500 d_deltatw(i, k) = dtimesub*(ff(i)+dtke(i,k) - &
1501 entr(i,k)*deltatw(i,k)/sigmaw(i) - &
1502 (death_rate(i)*sigmaw(i)+detr(i,k))*deltatw(i,k)/(1.-sigmaw(i)) - & ! cc
1503 8072060 tgw(i,k)*deltatw(i,k) )
1504 ELSE
1505 d_deltatw(i, k) = 1/tgw(i, k)*(1-exp(-dtimesub*tgw(i,k)))* &
1506 (ff(i)+dtke(i,k) - &
1507 entr(i,k)*deltatw(i,k)/sigmaw(i) - &
1508 (death_rate(i)*sigmaw(i)+detr(i,k))*deltatw(i,k)/(1.-sigmaw(i)) - &
1509 44478390 tgw(i,k)*deltatw(i,k) )
1510 END IF
1511
1512 52550450 dth(i, k) = deltatw(i, k)/ppi(i, k)
1513
1514 52550450 gg(i) = d_deltaqw(i, k)/dtimesub
1515
1516 d_deltaqw(i, k) = dtimesub*(gg(i)+dqke(i,k) - &
1517 entr(i,k)*deltaqw(i,k)/sigmaw(i) - &
1518 52550450 (death_rate(i)*sigmaw(i)+detr(i,k))*deltaqw(i,k)/(1.-sigmaw(i)))
1519 ! cc
1520
1521 ! cc nrlmd
1522 ! cc d_deltatw2(i,k)=d_deltatw2(i,k)+d_deltatw(i,k)
1523 ! cc d_deltaqw2(i,k)=d_deltaqw2(i,k)+d_deltaqw(i,k)
1524 ! cc
1525 END IF
1526 END DO
1527 END DO
1528
1529
1530 ! Scale tendencies so that water vapour remains positive in w and x.
1531
1532 CALL wake_vec_modulation(klon, klev, wk_adv, epsilon, qe, d_qe, deltaqw, &
1533 4800 d_deltaqw, sigmaw, d_sigmaw, alpha)
1534
1535 ! cc nrlmd
1536 ! c print*,'alpha'
1537 ! c do i=1,klon
1538 ! c print*,alpha(i)
1539 ! c end do
1540 ! cc
1541
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1542
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1543
3/4
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 52550450 times.
✓ Branch 3 taken 133526350 times.
186264000 IF (wk_adv(i) .AND. k<=kupper(i)) THEN
1544 52550450 d_te(i, k) = alpha(i)*d_te(i, k)
1545 52550450 d_qe(i, k) = alpha(i)*d_qe(i, k)
1546 52550450 d_deltatw(i, k) = alpha(i)*d_deltatw(i, k)
1547 52550450 d_deltaqw(i, k) = alpha(i)*d_deltaqw(i, k)
1548 52550450 d_deltat_gw(i, k) = alpha(i)*d_deltat_gw(i, k)
1549 END IF
1550 END DO
1551 END DO
1552
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1553
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1554 4771200 d_sigmaw(i) = alpha(i)*d_sigmaw(i)
1555 END IF
1556 END DO
1557
1558 ! Update large scale variables and wake variables
1559 ! IM 060208 manque DO i + remplace DO k=1,kupper(i)
1560 ! IM 060208 DO k = 1,kupper(i)
1561
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1562
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1563
3/4
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 52550450 times.
✓ Branch 3 taken 133526350 times.
186264000 IF (wk_adv(i) .AND. k<=kupper(i)) THEN
1564 52550450 dtls(i, k) = dtls(i, k) + d_te(i, k)
1565 52550450 dqls(i, k) = dqls(i, k) + d_qe(i, k)
1566 ! cc nrlmd
1567 52550450 d_deltatw2(i, k) = d_deltatw2(i, k) + d_deltatw(i, k)
1568 52550450 d_deltaqw2(i, k) = d_deltaqw2(i, k) + d_deltaqw(i, k)
1569 ! cc
1570 END IF
1571 END DO
1572 END DO
1573
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1574
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1575
3/4
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 52550450 times.
✓ Branch 3 taken 133526350 times.
186264000 IF (wk_adv(i) .AND. k<=kupper(i)) THEN
1576 52550450 te(i, k) = te0(i, k) + dtls(i, k)
1577 52550450 qe(i, k) = qe0(i, k) + dqls(i, k)
1578 52550450 the(i, k) = te(i, k)/ppi(i, k)
1579 52550450 deltatw(i, k) = deltatw(i, k) + d_deltatw(i, k)
1580 52550450 deltaqw(i, k) = deltaqw(i, k) + d_deltaqw(i, k)
1581 52550450 dth(i, k) = deltatw(i, k)/ppi(i, k)
1582 ! c print*,'k,qx,qw',k,qe(i,k)-sigmaw(i)*deltaqw(i,k)
1583 ! c $ ,qe(i,k)+(1-sigmaw(i))*deltaqw(i,k)
1584 END IF
1585 END DO
1586 END DO
1587 !
1588
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1589
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1590 4771200 sigmaw(i) = sigmaw(i) + d_sigmaw(i)
1591 4771200 d_sigmaw2(i) = d_sigmaw2(i) + d_sigmaw(i)
1592 END IF
1593 END DO
1594 !jyg<
1595
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 4800 times.
4800 IF (iflag_wk_pop_dyn >= 1) THEN
1596 DO i = 1, klon
1597 IF (wk_adv(i)) THEN
1598 awdens(i) = awdens(i) + d_awdens(i)
1599 wdens(i) = wdens(i) + d_wdens(i)
1600 d_awdens2(i) = d_awdens2(i) + d_awdens(i)
1601 d_wdens2(i) = d_wdens2(i) + d_wdens(i)
1602 END IF
1603 END DO
1604 DO i = 1, klon
1605 IF (wk_adv(i)) THEN
1606 wdens_targ = max(wdens(i),wdensmin)
1607 d_wdens2(i) = d_wdens2(i) + wdens_targ - wdens(i)
1608 wdens(i) = wdens_targ
1609 !
1610 wdens_targ = min( max(awdens(i),0.), wdens(i) )
1611 d_awdens2(i) = d_awdens2(i) + wdens_targ - awdens(i)
1612 awdens(i) = wdens_targ
1613 END IF
1614 END DO
1615 DO i = 1, klon
1616 IF (wk_adv(i)) THEN
1617 sigmaw_targ = max(sigmaw(i),sigmad)
1618 d_sigmaw2(i) = d_sigmaw2(i) + sigmaw_targ - sigmaw(i)
1619 sigmaw(i) = sigmaw_targ
1620 END IF
1621 END DO
1622 ENDIF ! (iflag_wk_pop_dyn >= 1)
1623 !>jyg
1624
1625
1626 ! Determine Ptop from buoyancy integral
1627 ! ---------------------------------------
1628
1629 ! - 1/ Pressure of the level where dth changes sign.
1630
1631
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1632
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1633 4771200 ptop_provis(i) = ph(i, 1)
1634 END IF
1635 END DO
1636
1637
2/2
✓ Branch 0 taken 182400 times.
✓ Branch 1 taken 4800 times.
187200 DO k = 2, klev
1638
2/2
✓ Branch 0 taken 181305600 times.
✓ Branch 1 taken 182400 times.
181492800 DO i = 1, klon
1639 IF (wk_adv(i) .AND. ptop_provis(i)==ph(i,1) .AND. &
1640
7/8
✓ Branch 0 taken 181305600 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 150956965 times.
✓ Branch 3 taken 30348635 times.
✓ Branch 4 taken 149358023 times.
✓ Branch 5 taken 1598942 times.
✓ Branch 6 taken 863829 times.
✓ Branch 7 taken 148494194 times.
181488000 dth(i,k)>-delta_t_min .AND. dth(i,k-1)<-delta_t_min) THEN
1641 ptop_provis(i) = ((dth(i,k)+delta_t_min)*p(i,k-1) - &
1642 863829 (dth(i,k-1)+delta_t_min)*p(i,k))/(dth(i,k)-dth(i,k-1))
1643 END IF
1644 END DO
1645 END DO
1646
1647 ! - 2/ dth integral
1648
1649
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1650
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN !!! nrlmd
1651 4771200 sum_dth(i) = 0.
1652 4771200 dthmin(i) = -delta_t_min
1653 4771200 z(i) = 0.
1654 END IF
1655 END DO
1656
1657
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1658
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1659
1/2
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
186264000 IF (wk_adv(i)) THEN
1660 186076800 dz(i) = -(amax1(ph(i,k+1),ptop_provis(i))-ph(i,k))/(rho(i,k)*rg)
1661
2/2
✓ Branch 0 taken 2948975 times.
✓ Branch 1 taken 183127825 times.
186076800 IF (dz(i)>0) THEN
1662 2948975 z(i) = z(i) + dz(i)
1663 2948975 sum_dth(i) = sum_dth(i) + dth(i, k)*dz(i)
1664 2948975 dthmin(i) = amin1(dthmin(i), dth(i,k))
1665 END IF
1666 END IF
1667 END DO
1668 END DO
1669
1670 ! - 3/ height of triangle with area= sum_dth and base = dthmin
1671
1672
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1673
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1674 4771200 hw(i) = 2.*sum_dth(i)/amin1(dthmin(i), -0.5)
1675 4771200 hw(i) = amax1(hwmin, hw(i))
1676 END IF
1677 END DO
1678
1679 ! - 4/ now, get Ptop
1680
1681
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1682
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN !!! nrlmd
1683 4771200 ktop(i) = 0
1684 4771200 z(i) = 0.
1685 END IF
1686 END DO
1687
1688
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1689
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1690
1/2
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
186264000 IF (wk_adv(i)) THEN
1691 186076800 dz(i) = amin1(-(ph(i,k+1)-ph(i,k))/(rho(i,k)*rg), hw(i)-z(i))
1692
2/2
✓ Branch 0 taken 7162340 times.
✓ Branch 1 taken 178914460 times.
186076800 IF (dz(i)>0) THEN
1693 7162340 z(i) = z(i) + dz(i)
1694 7162340 ptop(i) = ph(i, k) - rho(i, k)*rg*dz(i)
1695 7162340 ktop(i) = k
1696 END IF
1697 END IF
1698 END DO
1699 END DO
1700
1701 ! 4.5/Correct ktop and ptop
1702
1703
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1704
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1705 4771200 ptop_new(i) = ptop(i)
1706 END IF
1707 END DO
1708
1709
2/2
✓ Branch 0 taken 182400 times.
✓ Branch 1 taken 4800 times.
187200 DO k = klev, 2, -1
1710
2/2
✓ Branch 0 taken 181305600 times.
✓ Branch 1 taken 182400 times.
181492800 DO i = 1, klon
1711 ! IM v3JYG; IF (k .GE. ktop(i)
1712 IF (wk_adv(i) .AND. k<=ktop(i) .AND. ptop_new(i)==ptop(i) .AND. &
1713
9/10
✓ Branch 0 taken 181305600 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2391140 times.
✓ Branch 3 taken 178914460 times.
✓ Branch 4 taken 1199289 times.
✓ Branch 5 taken 1191851 times.
✓ Branch 6 taken 839441 times.
✓ Branch 7 taken 359848 times.
✓ Branch 8 taken 740622 times.
✓ Branch 9 taken 98819 times.
181488000 dth(i,k)>-delta_t_min .AND. dth(i,k-1)<-delta_t_min) THEN
1714 ptop_new(i) = ((dth(i,k)+delta_t_min)*p(i,k-1) - &
1715 740622 (dth(i,k-1)+delta_t_min)*p(i,k))/(dth(i,k)-dth(i,k-1))
1716 END IF
1717 END DO
1718 END DO
1719
1720
1721
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1722
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN
1723 4771200 ptop(i) = ptop_new(i)
1724 END IF
1725 END DO
1726
1727
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = klev, 1, -1
1728
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1729
1/2
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
186264000 IF (wk_adv(i)) THEN !!! nrlmd
1730
2/2
✓ Branch 0 taken 184074525 times.
✓ Branch 1 taken 2002275 times.
186076800 IF (ph(i,k+1)<ptop(i)) ktop(i) = k
1731 END IF
1732 END DO
1733 END DO
1734
1735 ! 5/ Set deltatw & deltaqw to 0 above kupper
1736
1737
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1738
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1739
3/4
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 138297550 times.
✓ Branch 3 taken 47779250 times.
186264000 IF (wk_adv(i) .AND. k>=kupper(i)) THEN
1740 138297550 deltatw(i, k) = 0.
1741 138297550 deltaqw(i, k) = 0.
1742 138297550 d_deltatw2(i,k) = -deltatw0(i,k)
1743 138297550 d_deltaqw2(i,k) = -deltaqw0(i,k)
1744 END IF
1745 END DO
1746 END DO
1747
1748
1749 ! -------------Cstar computation---------------------------------
1750
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1751
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN !!! nrlmd
1752 4771200 sum_thu(i) = 0.
1753 4771200 sum_tu(i) = 0.
1754 4771200 sum_qu(i) = 0.
1755 4771200 sum_thvu(i) = 0.
1756 4771200 sum_dth(i) = 0.
1757 4771200 sum_dq(i) = 0.
1758 4771200 sum_rho(i) = 0.
1759 4771200 sum_dtdwn(i) = 0.
1760 4771200 sum_dqdwn(i) = 0.
1761
1762 4771200 av_thu(i) = 0.
1763 4771200 av_tu(i) = 0.
1764 4771200 av_qu(i) = 0.
1765 4771200 av_thvu(i) = 0.
1766 4771200 av_dth(i) = 0.
1767 4771200 av_dq(i) = 0.
1768 4771200 av_rho(i) = 0.
1769 4771200 av_dtdwn(i) = 0.
1770 4771200 av_dqdwn(i) = 0.
1771 END IF
1772 END DO
1773
1774 ! Integrals (and wake top level number)
1775 ! --------------------------------------
1776
1777 ! Initialize sum_thvu to 1st level virt. pot. temp.
1778
1779
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1780
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN !!! nrlmd
1781 4771200 z(i) = 1.
1782 4771200 dz(i) = 1.
1783 4771200 sum_thvu(i) = thu(i, 1)*(1.+epsim1*qu(i,1))*dz(i)
1784 4771200 sum_dth(i) = 0.
1785 END IF
1786 END DO
1787
1788
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1789
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1790
1/2
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
186264000 IF (wk_adv(i)) THEN !!! nrlmd
1791 186076800 dz(i) = -(max(ph(i,k+1),ptop(i))-ph(i,k))/(rho(i,k)*rg)
1792
2/2
✓ Branch 0 taken 6773474 times.
✓ Branch 1 taken 179303326 times.
186076800 IF (dz(i)>0) THEN
1793 6773474 z(i) = z(i) + dz(i)
1794 6773474 sum_thu(i) = sum_thu(i) + thu(i, k)*dz(i)
1795 6773474 sum_tu(i) = sum_tu(i) + tu(i, k)*dz(i)
1796 6773474 sum_qu(i) = sum_qu(i) + qu(i, k)*dz(i)
1797 6773474 sum_thvu(i) = sum_thvu(i) + thu(i, k)*(1.+epsim1*qu(i,k))*dz(i)
1798 6773474 sum_dth(i) = sum_dth(i) + dth(i, k)*dz(i)
1799 6773474 sum_dq(i) = sum_dq(i) + deltaqw(i, k)*dz(i)
1800 6773474 sum_rho(i) = sum_rho(i) + rhow(i, k)*dz(i)
1801 6773474 sum_dtdwn(i) = sum_dtdwn(i) + dtdwn(i, k)*dz(i)
1802 6773474 sum_dqdwn(i) = sum_dqdwn(i) + dqdwn(i, k)*dz(i)
1803 END IF
1804 END IF
1805 END DO
1806 END DO
1807
1808
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1809
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN !!! nrlmd
1810 4771200 hw0(i) = z(i)
1811 END IF
1812 END DO
1813
1814
1815 ! - WAPE and mean forcing computation
1816 ! ---------------------------------------
1817
1818 ! ---------------------------------------
1819
1820 ! Means
1821
1822
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776000 DO i = 1, klon
1823
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN !!! nrlmd
1824 4771200 av_thu(i) = sum_thu(i)/hw0(i)
1825 4771200 av_tu(i) = sum_tu(i)/hw0(i)
1826 4771200 av_qu(i) = sum_qu(i)/hw0(i)
1827 4771200 av_thvu(i) = sum_thvu(i)/hw0(i)
1828 4771200 av_dth(i) = sum_dth(i)/hw0(i)
1829 4771200 av_dq(i) = sum_dq(i)/hw0(i)
1830 4771200 av_rho(i) = sum_rho(i)/hw0(i)
1831 4771200 av_dtdwn(i) = sum_dtdwn(i)/hw0(i)
1832 4771200 av_dqdwn(i) = sum_dqdwn(i)/hw0(i)
1833
1834 wape(i) = -rg*hw0(i)*(av_dth(i)+epsim1*(av_thu(i)*av_dq(i) + &
1835 4771200 av_dth(i)*av_qu(i)+av_dth(i)*av_dq(i)))/av_thvu(i)
1836 END IF
1837 END DO
1838
1839 ! Filter out bad wakes
1840
1841
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, klev
1842
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, klon
1843
1/2
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
186264000 IF (wk_adv(i)) THEN !!! nrlmd
1844
2/2
✓ Branch 0 taken 150384 times.
✓ Branch 1 taken 185926416 times.
186076800 IF (wape(i)<0.) THEN
1845 150384 deltatw(i, k) = 0.
1846 150384 deltaqw(i, k) = 0.
1847 150384 dth(i, k) = 0.
1848 150384 d_deltatw2(i,k) = -deltatw0(i,k)
1849 150384 d_deltaqw2(i,k) = -deltaqw0(i,k)
1850 END IF
1851 END IF
1852 END DO
1853 END DO
1854
1855
2/2
✓ Branch 0 taken 4771200 times.
✓ Branch 1 taken 4800 times.
4776480 DO i = 1, klon
1856
1/2
✓ Branch 0 taken 4771200 times.
✗ Branch 1 not taken.
4776000 IF (wk_adv(i)) THEN !!! nrlmd
1857
2/2
✓ Branch 0 taken 3856 times.
✓ Branch 1 taken 4767344 times.
4771200 IF (wape(i)<0.) THEN
1858 3856 wape(i) = 0.
1859 3856 cstar(i) = 0.
1860 3856 hw(i) = hwmin
1861 !jyg<
1862 !! sigmaw(i) = max(sigmad, sigd_con(i))
1863 3856 sigmaw_targ = max(sigmad, sigd_con(i))
1864 3856 d_sigmaw2(i) = d_sigmaw2(i) + sigmaw_targ - sigmaw(i)
1865 3856 sigmaw(i) = sigmaw_targ
1866 !>jyg
1867 3856 fip(i) = 0.
1868 3856 gwake(i) = .FALSE.
1869 ELSE
1870 4767344 cstar(i) = stark*sqrt(2.*wape(i))
1871 4767344 gwake(i) = .TRUE.
1872 END IF
1873 END IF
1874 END DO
1875
1876 END DO ! end sub-timestep loop
1877
1878
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
480 IF (prt_level>=10) THEN
1879 PRINT *, 'wake-5, sigmaw(igout), cstar(igout), wape(igout), ptop(igout) ', &
1880 sigmaw(igout), cstar(igout), wape(igout), ptop(igout)
1881 ENDIF
1882
1883
1884 ! ----------------------------------------------------------
1885 ! Determine wake final state; recompute wape, cstar, ktop;
1886 ! filter out bad wakes.
1887 ! ----------------------------------------------------------
1888
1889 ! 2.1 - Undisturbed area and Wake integrals
1890 ! ---------------------------------------------------------
1891
1892
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
1893 ! cc nrlmd if (wk_adv(i)) then !!! nrlmd
1894
1/2
✓ Branch 0 taken 477120 times.
✗ Branch 1 not taken.
477600 IF (ok_qx_qw(i)) THEN
1895 ! cc
1896 477120 z(i) = 0.
1897 477120 sum_thu(i) = 0.
1898 477120 sum_tu(i) = 0.
1899 477120 sum_qu(i) = 0.
1900 477120 sum_thvu(i) = 0.
1901 477120 sum_dth(i) = 0.
1902 477120 sum_half_dth(i) = 0.
1903 477120 sum_dq(i) = 0.
1904 477120 sum_rho(i) = 0.
1905 477120 sum_dtdwn(i) = 0.
1906 477120 sum_dqdwn(i) = 0.
1907
1908 477120 av_thu(i) = 0.
1909 477120 av_tu(i) = 0.
1910 477120 av_qu(i) = 0.
1911 477120 av_thvu(i) = 0.
1912 477120 av_dth(i) = 0.
1913 477120 av_dq(i) = 0.
1914 477120 av_rho(i) = 0.
1915 477120 av_dtdwn(i) = 0.
1916 477120 av_dqdwn(i) = 0.
1917
1918 477120 dthmin(i) = -delta_t_min
1919 END IF
1920 END DO
1921 ! Potential temperatures and humidity
1922 ! ----------------------------------------------------------
1923
1924
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = 1, klev
1925
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
1926 ! cc nrlmd IF ( wk_adv(i)) THEN
1927
1/2
✓ Branch 0 taken 18607680 times.
✗ Branch 1 not taken.
18626400 IF (ok_qx_qw(i)) THEN
1928 ! cc
1929 18607680 rho(i, k) = p(i, k)/(rd*te(i,k))
1930
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 18130560 times.
18607680 IF (k==1) THEN
1931 477120 rhoh(i, k) = ph(i, k)/(rd*te(i,k))
1932 477120 zhh(i, k) = 0
1933 ELSE
1934 18130560 rhoh(i, k) = ph(i, k)*2./(rd*(te(i,k)+te(i,k-1)))
1935 18130560 zhh(i, k) = (ph(i,k)-ph(i,k-1))/(-rhoh(i,k)*rg) + zhh(i, k-1)
1936 END IF
1937 18607680 the(i, k) = te(i, k)/ppi(i, k)
1938 18607680 thu(i, k) = (te(i,k)-deltatw(i,k)*sigmaw(i))/ppi(i, k)
1939 18607680 tu(i, k) = te(i, k) - deltatw(i, k)*sigmaw(i)
1940 18607680 qu(i, k) = qe(i, k) - deltaqw(i, k)*sigmaw(i)
1941 18607680 rhow(i, k) = p(i, k)/(rd*(te(i,k)+deltatw(i,k)))
1942 18607680 dth(i, k) = deltatw(i, k)/ppi(i, k)
1943 END IF
1944 END DO
1945 END DO
1946
1947 ! Integrals (and wake top level number)
1948 ! -----------------------------------------------------------
1949
1950 ! Initialize sum_thvu to 1st level virt. pot. temp.
1951
1952
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
1953 ! cc nrlmd IF ( wk_adv(i)) THEN
1954
1/2
✓ Branch 0 taken 477120 times.
✗ Branch 1 not taken.
477600 IF (ok_qx_qw(i)) THEN
1955 ! cc
1956 477120 z(i) = 1.
1957 477120 dz(i) = 1.
1958 477120 dz_half(i) = 1.
1959 477120 sum_thvu(i) = thu(i, 1)*(1.+epsim1*qu(i,1))*dz(i)
1960 477120 sum_dth(i) = 0.
1961 END IF
1962 END DO
1963
1964
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = 1, klev
1965
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
1966 ! cc nrlmd IF ( wk_adv(i)) THEN
1967
1/2
✓ Branch 0 taken 18607680 times.
✗ Branch 1 not taken.
18626400 IF (ok_qx_qw(i)) THEN
1968 ! cc
1969 18607680 dz(i) = -(amax1(ph(i,k+1),ptop(i))-ph(i,k))/(rho(i,k)*rg)
1970 18607680 dz_half(i) = -(amax1(ph(i,k+1),0.5*(ptop(i)+ph(i,1)))-ph(i,k))/(rho(i,k)*rg)
1971
2/2
✓ Branch 0 taken 673061 times.
✓ Branch 1 taken 17934619 times.
18607680 IF (dz(i)>0) THEN
1972 673061 z(i) = z(i) + dz(i)
1973 673061 sum_thu(i) = sum_thu(i) + thu(i, k)*dz(i)
1974 673061 sum_tu(i) = sum_tu(i) + tu(i, k)*dz(i)
1975 673061 sum_qu(i) = sum_qu(i) + qu(i, k)*dz(i)
1976 673061 sum_thvu(i) = sum_thvu(i) + thu(i, k)*(1.+epsim1*qu(i,k))*dz(i)
1977 673061 sum_dth(i) = sum_dth(i) + dth(i, k)*dz(i)
1978 673061 sum_dq(i) = sum_dq(i) + deltaqw(i, k)*dz(i)
1979 673061 sum_rho(i) = sum_rho(i) + rhow(i, k)*dz(i)
1980 673061 sum_dtdwn(i) = sum_dtdwn(i) + dtdwn(i, k)*dz(i)
1981 673061 sum_dqdwn(i) = sum_dqdwn(i) + dqdwn(i, k)*dz(i)
1982 !
1983 673061 dthmin(i) = min(dthmin(i), dth(i,k))
1984 END IF
1985
2/2
✓ Branch 0 taken 580251 times.
✓ Branch 1 taken 18027429 times.
18607680 IF (dz_half(i)>0) THEN
1986 580251 sum_half_dth(i) = sum_half_dth(i) + dth(i, k)*dz_half(i)
1987 END IF
1988 END IF
1989 END DO
1990 END DO
1991
1992
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
1993 ! cc nrlmd IF ( wk_adv(i)) THEN
1994
1/2
✓ Branch 0 taken 477120 times.
✗ Branch 1 not taken.
477600 IF (ok_qx_qw(i)) THEN
1995 ! cc
1996 477120 hw0(i) = z(i)
1997 END IF
1998 END DO
1999
2000 ! - WAPE and mean forcing computation
2001 ! -------------------------------------------------------------
2002
2003 ! Means
2004
2005
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
2006 ! cc nrlmd IF ( wk_adv(i)) THEN
2007
1/2
✓ Branch 0 taken 477120 times.
✗ Branch 1 not taken.
477600 IF (ok_qx_qw(i)) THEN
2008 ! cc
2009 477120 av_thu(i) = sum_thu(i)/hw0(i)
2010 477120 av_tu(i) = sum_tu(i)/hw0(i)
2011 477120 av_qu(i) = sum_qu(i)/hw0(i)
2012 477120 av_thvu(i) = sum_thvu(i)/hw0(i)
2013 477120 av_dth(i) = sum_dth(i)/hw0(i)
2014 477120 av_dq(i) = sum_dq(i)/hw0(i)
2015 477120 av_rho(i) = sum_rho(i)/hw0(i)
2016 477120 av_dtdwn(i) = sum_dtdwn(i)/hw0(i)
2017 477120 av_dqdwn(i) = sum_dqdwn(i)/hw0(i)
2018
2019 wape2(i) = -rg*hw0(i)*(av_dth(i)+epsim1*(av_thu(i)*av_dq(i) + &
2020 477120 av_dth(i)*av_qu(i)+av_dth(i)*av_dq(i)))/av_thvu(i)
2021 END IF
2022 END DO
2023
2024
2025
2026 ! Prognostic variable update
2027 ! ------------------------------------------------------------
2028
2029 ! Filter out bad wakes
2030
2031
1/2
✓ Branch 0 taken 480 times.
✗ Branch 1 not taken.
480 IF (iflag_wk_check_trgl>=1) THEN
2032 ! Check triangular shape of dth profile
2033
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
2034
1/2
✓ Branch 0 taken 477120 times.
✗ Branch 1 not taken.
477600 IF (ok_qx_qw(i)) THEN
2035 !! print *,'wake, hw0(i), dthmin(i) ', hw0(i), dthmin(i)
2036 !! print *,'wake, 2.*sum_dth(i)/(hw0(i)*dthmin(i)) ', &
2037 !! 2.*sum_dth(i)/(hw0(i)*dthmin(i))
2038 !! print *,'wake, sum_half_dth(i), sum_dth(i) ', &
2039 !! sum_half_dth(i), sum_dth(i)
2040
3/4
✓ Branch 0 taken 477120 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 390227 times.
✓ Branch 3 taken 86893 times.
477120 IF ((hw0(i) < 1.) .or. (dthmin(i) >= -delta_t_min) ) THEN
2041 390227 wape2(i) = -1.
2042 !! print *,'wake, rej 1'
2043
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 86893 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
86893 ELSE IF (iflag_wk_check_trgl==1.AND.abs(2.*sum_dth(i)/(hw0(i)*dthmin(i)) - 1.) > 0.5) THEN
2044 wape2(i) = -1.
2045 !! print *,'wake, rej 2'
2046
2/2
✓ Branch 0 taken 1239 times.
✓ Branch 1 taken 85654 times.
86893 ELSE IF (abs(sum_half_dth(i)) < 0.5*abs(sum_dth(i)) ) THEN
2047 1239 wape2(i) = -1.
2048 !! print *,'wake, rej 3'
2049 END IF
2050 END IF
2051 END DO
2052 END IF
2053
2054
2055
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = 1, klev
2056
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
2057 ! cc nrlmd IF ( wk_adv(i) .AND. wape2(i) .LT. 0.) THEN
2058
3/4
✓ Branch 0 taken 18607680 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 15267174 times.
✓ Branch 3 taken 3340506 times.
18626400 IF (ok_qx_qw(i) .AND. wape2(i)<0.) THEN
2059 ! cc
2060 15267174 deltatw(i, k) = 0.
2061 15267174 deltaqw(i, k) = 0.
2062 15267174 dth(i, k) = 0.
2063 15267174 d_deltatw2(i,k) = -deltatw0(i,k)
2064 15267174 d_deltaqw2(i,k) = -deltaqw0(i,k)
2065 END IF
2066 END DO
2067 END DO
2068
2069
2070
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
2071 ! cc nrlmd IF ( wk_adv(i)) THEN
2072
1/2
✓ Branch 0 taken 477120 times.
✗ Branch 1 not taken.
477600 IF (ok_qx_qw(i)) THEN
2073 ! cc
2074
2/2
✓ Branch 0 taken 391466 times.
✓ Branch 1 taken 85654 times.
477120 IF (wape2(i)<0.) THEN
2075 391466 wape2(i) = 0.
2076 391466 cstar2(i) = 0.
2077 391466 hw(i) = hwmin
2078 !jyg<
2079 !! sigmaw(i) = amax1(sigmad, sigd_con(i))
2080 391466 sigmaw_targ = max(sigmad, sigd_con(i))
2081 391466 d_sigmaw2(i) = d_sigmaw2(i) + sigmaw_targ - sigmaw(i)
2082 391466 sigmaw(i) = sigmaw_targ
2083 !>jyg
2084 391466 fip(i) = 0.
2085 391466 gwake(i) = .FALSE.
2086 ELSE
2087
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 85654 times.
85654 IF (prt_level>=10) PRINT *, 'wape2>0'
2088 85654 cstar2(i) = stark*sqrt(2.*wape2(i))
2089 85654 gwake(i) = .TRUE.
2090 END IF
2091 END IF
2092 END DO
2093
2094
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
2095 ! cc nrlmd IF ( wk_adv(i)) THEN
2096
1/2
✓ Branch 0 taken 477120 times.
✗ Branch 1 not taken.
477600 IF (ok_qx_qw(i)) THEN
2097 ! cc
2098 477120 ktopw(i) = ktop(i)
2099 END IF
2100 END DO
2101
2102
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
2103 ! cc nrlmd IF ( wk_adv(i)) THEN
2104
1/2
✓ Branch 0 taken 477120 times.
✗ Branch 1 not taken.
477600 IF (ok_qx_qw(i)) THEN
2105 ! cc
2106
3/4
✓ Branch 0 taken 477120 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 85654 times.
✓ Branch 3 taken 391466 times.
477120 IF (ktopw(i)>0 .AND. gwake(i)) THEN
2107
2108 ! jyg1 Utilisation d'un h_efficace constant ( ~ feeding layer)
2109 ! cc heff = 600.
2110 ! Utilisation de la hauteur hw
2111 ! c heff = 0.7*hw
2112 85654 heff(i) = hw(i)
2113
2114 fip(i) = 0.5*rho(i, ktopw(i))*cstar2(i)**3*heff(i)*2* &
2115 85654 sqrt(sigmaw(i)*wdens(i)*3.14)
2116 85654 fip(i) = alpk*fip(i)
2117 ! jyg2
2118 ELSE
2119 391466 fip(i) = 0.
2120 END IF
2121 END IF
2122 END DO
2123
2124 ! Limitation de sigmaw
2125
2126 ! cc nrlmd
2127 ! DO i=1,klon
2128 ! IF (OK_qx_qw(i)) THEN
2129 ! IF (sigmaw(i).GE.sigmaw_max) sigmaw(i)=sigmaw_max
2130 ! ENDIF
2131 ! ENDDO
2132 ! cc
2133
2134 !jyg<
2135
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
480 IF (iflag_wk_pop_dyn >= 1) THEN
2136 DO i = 1, klon
2137 kill_wake(i) = ((wape(i)>=wape2(i)) .AND. (wape2(i)<=wapecut)) .OR. (ktopw(i)<=2) .OR. &
2138 .NOT. ok_qx_qw(i) .OR. (wdens(i) < 2.*wdensmin)
2139 ENDDO
2140 ELSE ! (iflag_wk_pop_dyn >= 1)
2141
2/2
✓ Branch 0 taken 480 times.
✓ Branch 1 taken 477120 times.
477600 DO i = 1, klon
2142 kill_wake(i) = ((wape(i)>=wape2(i)) .AND. (wape2(i)<=wapecut)) .OR. (ktopw(i)<=2) .OR. &
2143
7/8
✓ Branch 0 taken 444913 times.
✓ Branch 1 taken 32207 times.
✓ Branch 2 taken 53182 times.
✓ Branch 3 taken 391731 times.
✓ Branch 4 taken 82565 times.
✓ Branch 5 taken 2824 times.
✓ Branch 6 taken 82565 times.
✗ Branch 7 not taken.
560165 .NOT. ok_qx_qw(i)
2144 ENDDO
2145 ENDIF ! (iflag_wk_pop_dyn >= 1)
2146 !>jyg
2147
2148
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = 1, klev
2149
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
2150 !!jyg IF (((wape(i)>=wape2(i)) .AND. (wape2(i)<=wapecut)) .OR. (ktopw(i)<=2) .OR. &
2151 !!jyg .NOT. ok_qx_qw(i)) THEN
2152
2/2
✓ Branch 0 taken 15387645 times.
✓ Branch 1 taken 3220035 times.
18626400 IF (kill_wake(i)) THEN
2153 ! cc
2154 15387645 dtls(i, k) = 0.
2155 15387645 dqls(i, k) = 0.
2156 15387645 deltatw(i, k) = 0.
2157 15387645 deltaqw(i, k) = 0.
2158 15387645 d_deltatw2(i,k) = -deltatw0(i,k)
2159 15387645 d_deltaqw2(i,k) = -deltaqw0(i,k)
2160 END IF ! (kill_wake(i))
2161 END DO
2162 END DO
2163
2164
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
2165 !!jyg IF (((wape(i)>=wape2(i)) .AND. (wape2(i)<=wapecut)) .OR. (ktopw(i)<=2) .OR. &
2166 !!jyg .NOT. ok_qx_qw(i)) THEN
2167
2/2
✓ Branch 0 taken 394555 times.
✓ Branch 1 taken 82565 times.
477600 IF (kill_wake(i)) THEN
2168 394555 ktopw(i) = 0
2169 394555 wape(i) = 0.
2170 394555 cstar(i) = 0.
2171 !!jyg Outside subroutine "Wake" hw, wdens and sigmaw are zero when there are no wakes
2172 !! hw(i) = hwmin !jyg
2173 !! sigmaw(i) = sigmad !jyg
2174 394555 hw(i) = 0. !jyg
2175 394555 fip(i) = 0.
2176 !! sigmaw(i) = 0. !jyg
2177 sigmaw_targ = 0.
2178 394555 d_sigmaw2(i) = d_sigmaw2(i) + sigmaw_targ - sigmaw(i)
2179 394555 sigmaw(i) = sigmaw_targ
2180
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 394555 times.
394555 IF (iflag_wk_pop_dyn >= 1) THEN
2181 !! awdens(i) = 0.
2182 !! wdens(i) = 0.
2183 wdens_targ = 0.
2184 d_wdens2(i) = wdens_targ - wdens(i)
2185 wdens(i) = wdens_targ
2186 wdens_targ = 0.
2187 d_awdens2(i) = wdens_targ - awdens(i)
2188 awdens(i) = wdens_targ
2189 ENDIF ! (iflag_wk_pop_dyn >= 1)
2190 ELSE ! (kill_wake(i))
2191 82565 wape(i) = wape2(i)
2192 82565 cstar(i) = cstar2(i)
2193 END IF ! (kill_wake(i))
2194 ! c print*,'wape wape2 ktopw OK_qx_qw =',
2195 ! c $ wape(i),wape2(i),ktopw(i),OK_qx_qw(i)
2196 END DO
2197
2198
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
480 IF (prt_level>=10) THEN
2199 PRINT *, 'wake-6, wape wape2 ktopw OK_qx_qw =', &
2200 wape(igout),wape2(igout),ktopw(igout),OK_qx_qw(igout)
2201 ENDIF
2202
2203
2204 ! -----------------------------------------------------------------
2205 ! Get back to tendencies per second
2206
2207
2/2
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
19200 DO k = 1, klev
2208
2/2
✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
18626880 DO i = 1, klon
2209
2210 ! cc nrlmd IF ( wk_adv(i) .AND. k .LE. kupper(i)) THEN
2211 !jyg<
2212 !! IF (ok_qx_qw(i) .AND. k<=kupper(i)) THEN
2213
1/2
✓ Branch 0 taken 18607680 times.
✗ Branch 1 not taken.
18626400 IF (ok_qx_qw(i)) THEN
2214 !>jyg
2215 ! cc
2216 18607680 dtls(i, k) = dtls(i, k)/dtime
2217 18607680 dqls(i, k) = dqls(i, k)/dtime
2218 18607680 d_deltatw2(i, k) = d_deltatw2(i, k)/dtime
2219 18607680 d_deltaqw2(i, k) = d_deltaqw2(i, k)/dtime
2220 18607680 d_deltat_gw(i, k) = d_deltat_gw(i, k)/dtime
2221 ! c print*,'k,dqls,omg,entr,detr',k,dqls(i,k),omg(i,k),entr(i,k)
2222 ! c $ ,death_rate(i)*sigmaw(i)
2223 END IF
2224 END DO
2225 END DO
2226 !jyg<
2227
2/2
✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
477600 DO i = 1, klon
2228 477120 d_sigmaw2(i) = d_sigmaw2(i)/dtime
2229 477120 d_awdens2(i) = d_awdens2(i)/dtime
2230 477600 d_wdens2(i) = d_wdens2(i)/dtime
2231 ENDDO
2232 !>jyg
2233
2234
2235
2236 480 RETURN
2237 END SUBROUTINE wake
2238
2239 4800 SUBROUTINE wake_vec_modulation(nlon, nl, wk_adv, epsilon, qe, d_qe, deltaqw, &
2240 4800 d_deltaqw, sigmaw, d_sigmaw, alpha)
2241 ! ------------------------------------------------------
2242 ! Dtermination du coefficient alpha tel que les tendances
2243 ! corriges alpha*d_G, pour toutes les grandeurs G, correspondent
2244 ! a une humidite positive dans la zone (x) et dans la zone (w).
2245 ! ------------------------------------------------------
2246 IMPLICIT NONE
2247
2248 ! Input
2249 REAL qe(nlon, nl), d_qe(nlon, nl)
2250 REAL deltaqw(nlon, nl), d_deltaqw(nlon, nl)
2251 REAL sigmaw(nlon), d_sigmaw(nlon)
2252 LOGICAL wk_adv(nlon)
2253 INTEGER nl, nlon
2254 ! Output
2255 REAL alpha(nlon)
2256 ! Internal variables
2257 9600 REAL zeta(nlon, nl)
2258 4800 REAL alpha1(nlon)
2259 REAL x, a, b, c, discrim
2260 REAL epsilon
2261 ! DATA epsilon/1.e-15/
2262 INTEGER i,k
2263
2264
2/2
✓ Branch 0 taken 187200 times.
✓ Branch 1 taken 4800 times.
192000 DO k = 1, nl
2265
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186264000 DO i = 1, nlon
2266
1/2
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
186264000 IF (wk_adv(i)) THEN
2267
2/2
✓ Branch 0 taken 178237248 times.
✓ Branch 1 taken 7839552 times.
186076800 IF ((deltaqw(i,k)+d_deltaqw(i,k))>=0.) THEN
2268 178237248 zeta(i, k) = 0.
2269 ELSE
2270 7839552 zeta(i, k) = 1.
2271 END IF
2272 END IF
2273 END DO
2274
2/2
✓ Branch 0 taken 186076800 times.
✓ Branch 1 taken 187200 times.
186268800 DO i = 1, nlon
2275
1/2
✓ Branch 0 taken 186076800 times.
✗ Branch 1 not taken.
186264000 IF (wk_adv(i)) THEN
2276 x = qe(i, k) + (zeta(i,k)-sigmaw(i))*deltaqw(i, k) + d_qe(i, k) + &
2277 (zeta(i,k)-sigmaw(i))*d_deltaqw(i, k) - d_sigmaw(i) * &
2278 186076800 (deltaqw(i,k)+d_deltaqw(i,k))
2279 186076800 a = -d_sigmaw(i)*d_deltaqw(i, k)
2280 b = d_qe(i, k) + (zeta(i,k)-sigmaw(i))*d_deltaqw(i, k) - &
2281 186076800 deltaqw(i, k)*d_sigmaw(i)
2282 186076800 c = qe(i, k) + (zeta(i,k)-sigmaw(i))*deltaqw(i, k) + epsilon
2283 186076800 discrim = b*b - 4.*a*c
2284 ! print*, 'x, a, b, c, discrim', x, a, b, c, discrim
2285
2/2
✓ Branch 0 taken 181842516 times.
✓ Branch 1 taken 4234284 times.
186076800 IF (a+b>=0.) THEN !! Condition suffisante pour la positivit� de ovap
2286 181842516 alpha1(i) = 1.
2287 ELSE
2288
1/2
✓ Branch 0 taken 4234284 times.
✗ Branch 1 not taken.
4234284 IF (x>=0.) THEN
2289 4234284 alpha1(i) = 1.
2290 ELSE
2291 IF (a>0.) THEN
2292 alpha1(i) = 0.9*min( (2.*c)/(-b+sqrt(discrim)), &
2293 (-b+sqrt(discrim))/(2.*a) )
2294 ELSE IF (a==0.) THEN
2295 alpha1(i) = 0.9*(-c/b)
2296 ELSE
2297 ! print*,'a,b,c discrim',a,b,c discrim
2298 alpha1(i) = 0.9*max( (2.*c)/(-b+sqrt(discrim)), &
2299 (-b+sqrt(discrim))/(2.*a))
2300 END IF
2301 END IF
2302 END IF
2303 186076800 alpha(i) = min(alpha(i), alpha1(i))
2304 END IF
2305 END DO
2306 END DO
2307
2308 4800 RETURN
2309 END SUBROUTINE wake_vec_modulation
2310
2311
2312
2313
2314