Directory: | ./ |
---|---|
File: | phys/acama_gwd_rando_m.f90 |
Date: | 2022-01-11 19:19:34 |
Exec | Total | Coverage | |
---|---|---|---|
Lines: | 137 | 144 | 95.1% |
Branches: | 204 | 222 | 91.9% |
Line | Branch | Exec | Source |
---|---|---|---|
1 | ! | ||
2 | ! $Id: acama_gwd_rando_m.F90 3977 2021-08-25 17:24:20Z fhourdin $ | ||
3 | ! | ||
4 | module ACAMA_GWD_rando_m | ||
5 | |||
6 | implicit none | ||
7 | |||
8 | contains | ||
9 | |||
10 |
5/10✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 480 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 480 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 480 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 480 times.
|
480 | SUBROUTINE ACAMA_GWD_rando(DTIME, pp, plat, tt, uu, vv, rot, & |
11 |
7/14✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 480 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 480 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 480 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 480 times.
✗ Branch 10 not taken.
✓ Branch 11 taken 480 times.
✗ Branch 12 not taken.
✓ Branch 13 taken 480 times.
|
480 | zustr, zvstr, d_u, d_v,east_gwstress,west_gwstress) |
12 | |||
13 | ! Parametrization of the momentum flux deposition due to a discrete | ||
14 | ! number of gravity waves. | ||
15 | ! Author: F. Lott, A. de la Camara | ||
16 | ! July, 24th, 2014 | ||
17 | ! Gaussian distribution of the source, source is vorticity squared | ||
18 | ! Reference: de la Camara and Lott (GRL, 2015, vol 42, 2071-2078 ) | ||
19 | ! Lott et al (JAS, 2010, vol 67, page 157-170) | ||
20 | ! Lott et al (JAS, 2012, vol 69, page 2134-2151) | ||
21 | |||
22 | ! ONLINE: | ||
23 | use dimphy, only: klon, klev | ||
24 | use assert_m, only: assert | ||
25 | USE ioipsl_getin_p_mod, ONLY : getin_p | ||
26 | USE vertical_layers_mod, ONLY : presnivs | ||
27 | |||
28 | include "YOMCST.h" | ||
29 | include "clesphys.h" | ||
30 | ! OFFLINE: | ||
31 | ! include "dimensions.h" | ||
32 | ! include "dimphy.h" | ||
33 | !END DIFFERENCE | ||
34 | include "YOEGWD.h" | ||
35 | |||
36 | ! 0. DECLARATIONS: | ||
37 | |||
38 | ! 0.1 INPUTS | ||
39 | REAL, intent(in)::DTIME ! Time step of the Physics | ||
40 | REAL, intent(in):: PP(:, :) ! (KLON, KLEV) Pressure at full levels | ||
41 | REAL, intent(in):: ROT(:,:) ! Relative vorticity | ||
42 | REAL, intent(in):: TT(:, :) ! (KLON, KLEV) Temp at full levels | ||
43 | REAL, intent(in):: UU(:, :) ! (KLON, KLEV) Zonal wind at full levels | ||
44 | REAL, intent(in):: VV(:, :) ! (KLON, KLEV) Merid wind at full levels | ||
45 | REAL, intent(in):: PLAT(:) ! (KLON) LATITUDE | ||
46 | |||
47 | ! 0.2 OUTPUTS | ||
48 | REAL, intent(out):: zustr(:), zvstr(:) ! (KLON) Surface Stresses | ||
49 | |||
50 | REAL, intent(inout):: d_u(:, :), d_v(:, :) | ||
51 | REAL, intent(inout):: east_gwstress(:, :) ! Profile of eastward stress | ||
52 | REAL, intent(inout):: west_gwstress(:, :) ! Profile of westward stress | ||
53 | ! (KLON, KLEV) tendencies on winds | ||
54 | |||
55 | ! O.3 INTERNAL ARRAYS | ||
56 | 960 | REAL BVLOW(klon) ! LOW LEVEL BV FREQUENCY | |
57 | 960 | REAL ROTBA(KLON),CORIO(KLON) ! BAROTROPIC REL. VORTICITY AND PLANETARY | |
58 | 960 | REAL UZ(KLON, KLEV + 1) | |
59 | |||
60 | INTEGER II, JJ, LL | ||
61 | |||
62 | ! 0.3.0 TIME SCALE OF THE LIFE CYCLE OF THE WAVES PARAMETERIZED | ||
63 | |||
64 | REAL DELTAT | ||
65 | |||
66 | ! 0.3.1 GRAVITY-WAVES SPECIFICATIONS | ||
67 | |||
68 | INTEGER, PARAMETER:: NK = 2, NP = 2, NO = 2, NW = NK * NP * NO | ||
69 | INTEGER JK, JP, JO, JW | ||
70 | INTEGER, PARAMETER:: NA = 5 !number of realizations to get the phase speed | ||
71 | REAL KMIN, KMAX ! Min and Max horizontal wavenumbers | ||
72 | REAL CMIN, CMAX ! Min and Max absolute ph. vel. | ||
73 | REAL CPHA ! absolute PHASE VELOCITY frequency | ||
74 | 960 | REAL ZK(NW, KLON) ! Horizontal wavenumber amplitude | |
75 | 960 | REAL ZP(NW, KLON) ! Horizontal wavenumber angle | |
76 | 960 | REAL ZO(NW, KLON) ! Absolute frequency ! | |
77 | |||
78 | ! Waves Intr. freq. at the 1/2 lev surrounding the full level | ||
79 | 960 | REAL ZOM(NW, KLON), ZOP(NW, KLON) | |
80 | |||
81 | ! Wave EP-fluxes at the 2 semi levels surrounding the full level | ||
82 | 960 | REAL WWM(NW, KLON), WWP(NW, KLON) | |
83 | |||
84 | 960 | REAL RUW0(NW, KLON) ! Fluxes at launching level | |
85 | |||
86 | 960 | REAL RUWP(NW, KLON), RVWP(NW, KLON) | |
87 | ! Fluxes X and Y for each waves at 1/2 Levels | ||
88 | |||
89 | INTEGER LAUNCH, LTROP ! Launching altitude and tropo altitude | ||
90 | |||
91 | REAL XLAUNCH ! Controle the launching altitude | ||
92 | REAL XTROP ! SORT of Tropopause altitude | ||
93 | 960 | REAL RUW(KLON, KLEV + 1) ! Flux x at semi levels | |
94 | 960 | REAL RVW(KLON, KLEV + 1) ! Flux y at semi levels | |
95 | |||
96 | REAL PRMAX ! Maximum value of PREC, and for which our linear formula | ||
97 | ! for GWs parameterisation apply | ||
98 | |||
99 | ! 0.3.2 PARAMETERS OF WAVES DISSIPATIONS | ||
100 | |||
101 | REAL RDISS, ZOISEC ! COEFF DE DISSIPATION, SECURITY FOR INTRINSIC FREQ | ||
102 | REAL CORSEC ! SECURITY FOR INTRINSIC CORIOLIS | ||
103 | REAL RUWFRT,SATFRT | ||
104 | |||
105 | ! 0.3.3 BACKGROUND FLOW AT 1/2 LEVELS AND VERTICAL COORDINATE | ||
106 | |||
107 | REAL H0 ! Characteristic Height of the atmosphere | ||
108 | REAL DZ ! Characteristic depth of the source! | ||
109 | REAL PR, TR ! Reference Pressure and Temperature | ||
110 | |||
111 | 960 | REAL ZH(KLON, KLEV + 1) ! Log-pressure altitude | |
112 | |||
113 | 960 | REAL UH(KLON, KLEV + 1), VH(KLON, KLEV + 1) ! Winds at 1/2 levels | |
114 | 960 | REAL PH(KLON, KLEV + 1) ! Pressure at 1/2 levels | |
115 | REAL PSEC ! Security to avoid division by 0 pressure | ||
116 | 960 | REAL PHM1(KLON, KLEV + 1) ! 1/Press at 1/2 levels | |
117 | 960 | REAL BV(KLON, KLEV + 1) ! Brunt Vaisala freq. (BVF) at 1/2 levels | |
118 | REAL BVSEC ! Security to avoid negative BVF | ||
119 | |||
120 | 960 | REAL, DIMENSION(klev+1) ::HREF | |
121 | LOGICAL, SAVE :: gwd_reproductibilite_mpiomp=.true. | ||
122 | LOGICAL, SAVE :: firstcall = .TRUE. | ||
123 | !$OMP THREADPRIVATE(firstcall,gwd_reproductibilite_mpiomp) | ||
124 | |||
125 | CHARACTER (LEN=20) :: modname='acama_gwd_rando_m' | ||
126 | CHARACTER (LEN=80) :: abort_message | ||
127 | |||
128 | |||
129 | |||
130 |
2/2✓ Branch 0 taken 1 times.
✓ Branch 1 taken 479 times.
|
480 | IF (firstcall) THEN |
131 | ! Cle introduite pour resoudre un probleme de non reproductibilite | ||
132 | ! Le but est de pouvoir tester de revenir a la version precedenete | ||
133 | ! A eliminer rapidement | ||
134 | 1 | CALL getin_p('gwd_reproductibilite_mpiomp',gwd_reproductibilite_mpiomp) | |
135 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1 times.
|
1 | IF (NW+4*(NA-1)+NA>=KLEV) THEN |
136 | ✗ | abort_message = 'NW+3*NA>=KLEV Probleme pour generation des ondes' | |
137 | ✗ | CALL abort_physic (modname,abort_message,1) | |
138 | ENDIF | ||
139 | 1 | firstcall=.false. | |
140 | ! CALL iophys_ini(dtime) | ||
141 | ENDIF | ||
142 | |||
143 | !----------------------------------------------------------------- | ||
144 | |||
145 | ! 1. INITIALISATIONS | ||
146 | |||
147 | ! 1.1 Basic parameter | ||
148 | |||
149 | ! Are provided from elsewhere (latent heat of vaporization, dry | ||
150 | ! gaz constant for air, gravity constant, heat capacity of dry air | ||
151 | ! at constant pressure, earth rotation rate, pi). | ||
152 | |||
153 | ! 1.2 Tuning parameters of V14 | ||
154 | |||
155 | ! Values for linear in rot (recommended): | ||
156 | ! RUWFRT=0.005 ! As RUWMAX but for frontal waves | ||
157 | ! SATFRT=1.00 ! As SAT but for frontal waves | ||
158 | ! Values when rot^2 is used | ||
159 | ! RUWFRT=0.02 ! As RUWMAX but for frontal waves | ||
160 | ! SATFRT=1.00 ! As SAT but for frontal waves | ||
161 | ! CMAX = 30. ! Characteristic phase speed | ||
162 | ! Values when rot^2*EXP(-pi*sqrt(J)) is used | ||
163 | ! RUWFRT=2.5 ! As RUWMAX but for frontal waves ~ N0*F0/4*DZ | ||
164 | ! SATFRT=0.60 ! As SAT but for frontal waves | ||
165 | 480 | RUWFRT=gwd_front_ruwmax | |
166 | 480 | SATFRT=gwd_front_sat | |
167 | CMAX = 50. ! Characteristic phase speed | ||
168 | ! Phase speed test | ||
169 | ! RUWFRT=0.01 | ||
170 | ! CMAX = 50. ! Characteristic phase speed (TEST) | ||
171 | ! Values when rot^2 and exp(-m^2*dz^2) are used | ||
172 | ! RUWFRT=0.03 ! As RUWMAX but for frontal waves | ||
173 | ! SATFRT=1.00 ! As SAT but for frontal waves | ||
174 | ! CRUCIAL PARAMETERS FOR THE WIND FILTERING | ||
175 | XLAUNCH=0.95 ! Parameter that control launching altitude | ||
176 | RDISS = 0.5 ! Diffusion parameter | ||
177 | |||
178 | ! maximum of rain for which our theory applies (in kg/m^2/s) | ||
179 | |||
180 | DZ = 1000. ! Characteristic depth of the source | ||
181 | XTROP=0.2 ! Parameter that control tropopause altitude | ||
182 | DELTAT=24.*3600. ! Time scale of the waves (first introduced in 9b) | ||
183 | ! DELTAT=DTIME ! No AR-1 Accumulation, OR OFFLINE | ||
184 | |||
185 | KMIN = 2.E-5 | ||
186 | ! minimum horizontal wavenumber (inverse of the subgrid scale resolution) | ||
187 | |||
188 | KMAX = 1.E-3 ! Max horizontal wavenumber | ||
189 | CMIN = 1. ! Min phase velocity | ||
190 | |||
191 | TR = 240. ! Reference Temperature | ||
192 | PR = 101300. ! Reference pressure | ||
193 | 480 | H0 = RD * TR / RG ! Characteristic vertical scale height | |
194 | |||
195 | BVSEC = 5.E-3 ! Security to avoid negative BVF | ||
196 | PSEC = 1.E-6 ! Security to avoid division by 0 pressure | ||
197 | ZOISEC = 1.E-6 ! Security FOR 0 INTRINSIC FREQ | ||
198 | 480 | CORSEC = ROMEGA*2.*SIN(2.*RPI/180.)! Security for CORIO | |
199 | |||
200 | ! ONLINE | ||
201 | call assert(klon == (/size(pp, 1), size(tt, 1), size(uu, 1), & | ||
202 | size(vv, 1), size(rot,1), size(zustr), size(zvstr), size(d_u, 1), & | ||
203 | size(d_v, 1), & | ||
204 | size(east_gwstress,1), size(west_gwstress,1) /), & | ||
205 |
2/2✓ Branch 0 taken 5280 times.
✓ Branch 1 taken 480 times.
|
5760 | "ACAMA_GWD_RANDO klon") |
206 | call assert(klev == (/size(pp, 2), size(tt, 2), size(uu, 2), & | ||
207 | size(vv, 2), size(d_u, 2), size(d_v, 2), & | ||
208 | size(east_gwstress,2), size(west_gwstress,2) /), & | ||
209 |
2/2✓ Branch 0 taken 3840 times.
✓ Branch 1 taken 480 times.
|
4320 | "ACAMA_GWD_RANDO klev") |
210 | ! END ONLINE | ||
211 | |||
212 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
|
480 | IF(DELTAT < DTIME)THEN |
213 | ! PRINT *, 'flott_gwd_rando: deltat < dtime!' | ||
214 | ! STOP 1 | ||
215 | ✗ | abort_message=' deltat < dtime! ' | |
216 | ✗ | CALL abort_physic(modname,abort_message,1) | |
217 | ENDIF | ||
218 | |||
219 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 480 times.
|
480 | IF (KLEV < NW) THEN |
220 | ! PRINT *, 'flott_gwd_rando: you will have problem with random numbers' | ||
221 | ! STOP 1 | ||
222 | ✗ | abort_message=' you will have problem with random numbers' | |
223 | ✗ | CALL abort_physic(modname,abort_message,1) | |
224 | ENDIF | ||
225 | |||
226 | ! 2. EVALUATION OF THE BACKGROUND FLOW AT SEMI-LEVELS | ||
227 | |||
228 | ! Pressure and Inv of pressure | ||
229 |
2/2✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 480 times.
|
18720 | DO LL = 2, KLEV |
230 |
2/2✓ Branch 0 taken 18130560 times.
✓ Branch 1 taken 18240 times.
|
18148800 | PH(:, LL) = EXP((LOG(PP(:, LL)) + LOG(PP(:, LL - 1))) / 2.) |
231 |
2/2✓ Branch 0 taken 18130560 times.
✓ Branch 1 taken 18240 times.
|
18149280 | PHM1(:, LL) = 1. / PH(:, LL) |
232 | end DO | ||
233 | |||
234 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | PH(:, KLEV + 1) = 0. |
235 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | PHM1(:, KLEV + 1) = 1. / PSEC |
236 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | PH(:, 1) = 2. * PP(:, 1) - PH(:, 2) |
237 | |||
238 | ! Launching altitude | ||
239 | |||
240 |
1/2✓ Branch 0 taken 480 times.
✗ Branch 1 not taken.
|
480 | IF (gwd_reproductibilite_mpiomp) THEN |
241 | ! Reprend la formule qui calcule PH en fonction de PP=play | ||
242 |
2/2✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 480 times.
|
18720 | DO LL = 2, KLEV |
243 | 18720 | HREF(LL) = EXP((LOG(presnivs(LL)) + LOG(presnivs(LL - 1))) / 2.) | |
244 | end DO | ||
245 | 480 | HREF(KLEV + 1) = 0. | |
246 | 480 | HREF(1) = 2. * presnivs(1) - HREF(2) | |
247 | ELSE | ||
248 | ✗ | HREF(1:KLEV)=PH(KLON/2,1:KLEV) | |
249 | ENDIF | ||
250 | |||
251 | LAUNCH=0 | ||
252 | LTROP =0 | ||
253 |
2/2✓ Branch 0 taken 480 times.
✓ Branch 1 taken 18720 times.
|
19200 | DO LL = 1, KLEV |
254 |
2/2✓ Branch 0 taken 1920 times.
✓ Branch 1 taken 16800 times.
|
19200 | IF (HREF(LL) / HREF(1) > XLAUNCH) LAUNCH = LL |
255 | ENDDO | ||
256 |
2/2✓ Branch 0 taken 480 times.
✓ Branch 1 taken 18720 times.
|
19200 | DO LL = 1, KLEV |
257 |
2/2✓ Branch 0 taken 8160 times.
✓ Branch 1 taken 10560 times.
|
19200 | IF (HREF(LL) / HREF(1) > XTROP) LTROP = LL |
258 | ENDDO | ||
259 | !LAUNCH=22 ; LTROP=33 | ||
260 | ! print*,'LAUNCH=',LAUNCH,'LTROP=',LTROP | ||
261 | |||
262 | |||
263 | ! PRINT *,'LAUNCH IN ACAMARA:',LAUNCH | ||
264 | |||
265 | ! Log pressure vert. coordinate | ||
266 |
2/2✓ Branch 0 taken 19200 times.
✓ Branch 1 taken 480 times.
|
19680 | DO LL = 1, KLEV + 1 |
267 |
2/2✓ Branch 0 taken 19084800 times.
✓ Branch 1 taken 19200 times.
|
19104480 | ZH(:, LL) = H0 * LOG(PR / (PH(:, LL) + PSEC)) |
268 | end DO | ||
269 | |||
270 | ! BV frequency | ||
271 |
2/2✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 480 times.
|
18720 | DO LL = 2, KLEV |
272 | ! BVSEC: BV Frequency (UH USED IS AS A TEMPORARY ARRAY DOWN TO WINDS) | ||
273 | UH(:, LL) = 0.5 * (TT(:, LL) + TT(:, LL - 1)) & | ||
274 | * RD**2 / RCPD / H0**2 + (TT(:, LL) & | ||
275 |
2/2✓ Branch 0 taken 18130560 times.
✓ Branch 1 taken 18240 times.
|
18149280 | - TT(:, LL - 1)) / (ZH(:, LL) - ZH(:, LL - 1)) * RD / H0 |
276 | end DO | ||
277 | BVLOW = 0.5 * (TT(:, LTROP )+ TT(:, LAUNCH)) & | ||
278 | * RD**2 / RCPD / H0**2 + (TT(:, LTROP ) & | ||
279 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | - TT(:, LAUNCH))/(ZH(:, LTROP )- ZH(:, LAUNCH)) * RD / H0 |
280 | |||
281 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | UH(:, 1) = UH(:, 2) |
282 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | UH(:, KLEV + 1) = UH(:, KLEV) |
283 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | BV(:, 1) = UH(:, 2) |
284 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | BV(:, KLEV + 1) = UH(:, KLEV) |
285 | ! SMOOTHING THE BV HELPS | ||
286 |
2/2✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 480 times.
|
18720 | DO LL = 2, KLEV |
287 |
2/2✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 18130560 times.
|
18149280 | BV(:, LL)=(UH(:, LL+1)+2.*UH(:, LL)+UH(:, LL-1))/4. |
288 | end DO | ||
289 | |||
290 |
4/4✓ Branch 0 taken 19200 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 19084800 times.
✓ Branch 3 taken 19200 times.
|
19104480 | BV=MAX(SQRT(MAX(BV, 0.)), BVSEC) |
291 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | BVLOW=MAX(SQRT(MAX(BVLOW, 0.)), BVSEC) |
292 | |||
293 | ! WINDS | ||
294 |
2/2✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 480 times.
|
18720 | DO LL = 2, KLEV |
295 |
2/2✓ Branch 0 taken 18130560 times.
✓ Branch 1 taken 18240 times.
|
18148800 | UH(:, LL) = 0.5 * (UU(:, LL) + UU(:, LL - 1)) ! Zonal wind |
296 |
2/2✓ Branch 0 taken 18130560 times.
✓ Branch 1 taken 18240 times.
|
18148800 | VH(:, LL) = 0.5 * (VV(:, LL) + VV(:, LL - 1)) ! Meridional wind |
297 | UZ(:, LL) = ABS((SQRT(UU(:, LL)**2+VV(:, LL)**2) & | ||
298 | - SQRT(UU(:,LL-1)**2+VV(:, LL-1)**2)) & | ||
299 |
2/2✓ Branch 0 taken 18130560 times.
✓ Branch 1 taken 18240 times.
|
18149280 | /(ZH(:, LL)-ZH(:, LL-1)) ) |
300 | end DO | ||
301 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | UH(:, 1) = 0. |
302 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | VH(:, 1) = 0. |
303 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | UH(:, KLEV + 1) = UU(:, KLEV) |
304 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | VH(:, KLEV + 1) = VV(:, KLEV) |
305 | |||
306 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | UZ(:, 1) = UZ(:, 2) |
307 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | UZ(:, KLEV + 1) = UZ(:, KLEV) |
308 |
4/4✓ Branch 0 taken 19200 times.
✓ Branch 1 taken 480 times.
✓ Branch 2 taken 19084800 times.
✓ Branch 3 taken 19200 times.
|
19104480 | UZ(:, :) = MAX(UZ(:,:), PSEC) |
309 | |||
310 | ! BAROTROPIC VORTICITY AND INTEGRATED CORIOLIS PARAMETER | ||
311 | |||
312 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | CORIO(:) = MAX(ROMEGA*2.*ABS(SIN(PLAT(:)*RPI/180.)),CORSEC) |
313 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | ROTBA(:)=0. |
314 |
2/2✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 480 times.
|
18720 | DO LL = 1,KLEV-1 |
315 | !ROTBA(:) = ROTBA(:) + (ROT(:,LL)+ROT(:,LL+1))/2./RG*(PP(:,LL)-PP(:,LL+1)) | ||
316 | ! Introducing the complete formula (exp of Richardson number): | ||
317 | ROTBA(:) = ROTBA(:) + & | ||
318 | !((ROT(:,LL)+ROT(:,LL+1))/2.)**2 & | ||
319 | (CORIO(:)*TANH(ABS(ROT(:,LL)+ROT(:,LL+1))/2./CORIO(:)))**2 & | ||
320 | /RG*(PP(:,LL)-PP(:,LL+1)) & | ||
321 | * EXP(-RPI*BV(:,LL+1)/UZ(:,LL+1)) & | ||
322 | ! * DZ*BV(:,LL+1)/4./ABS(CORIO(:)) | ||
323 |
2/2✓ Branch 0 taken 18240 times.
✓ Branch 1 taken 18130560 times.
|
18149280 | * DZ*BV(:,LL+1)/4./1.E-4 ! Changes after 1991 |
324 | !ARRET | ||
325 | ENDDO | ||
326 | ! PRINT *,'MAX ROTBA:',MAXVAL(ROTBA) | ||
327 | ! ROTBA(:)=(1.*ROTBA(:) & ! Testing zone | ||
328 | ! +0.15*CORIO(:)**2 & | ||
329 | ! /(COS(PLAT(:)*RPI/180.)+0.02) & | ||
330 | ! )*DZ*0.01/0.0001/4. ! & ! Testing zone | ||
331 | ! MODIF GWD4 AFTER 1985 | ||
332 | ! *(1.25+SIN(PLAT(:)*RPI/180.))/(1.05+SIN(PLAT(:)*RPI/180.))/1.25 | ||
333 | ! *1./(COS(PLAT(:)*RPI/180.)+0.02) | ||
334 | ! CORIO(:) = MAX(ROMEGA*2.*ABS(SIN(PLAT(:)*RPI/180.)),ZOISEC)/RG*PP(:,1) | ||
335 | |||
336 | ! 3 WAVES CHARACTERISTICS CHOSEN RANDOMLY AT THE LAUNCH ALTITUDE | ||
337 | |||
338 | ! The mod functions of weird arguments are used to produce the | ||
339 | ! waves characteristics in an almost stochastic way | ||
340 | |||
341 | JW = 0 | ||
342 |
2/2✓ Branch 0 taken 3840 times.
✓ Branch 1 taken 480 times.
|
4320 | DO JW = 1, NW |
343 | ! Angle | ||
344 |
2/2✓ Branch 0 taken 3816960 times.
✓ Branch 1 taken 3840 times.
|
3821280 | DO II = 1, KLON |
345 | ! Angle (0 or PI so far) | ||
346 | ! ZP(JW, II) = (SIGN(1., 0.5 - MOD(TT(II, JW) * 10., 1.)) + 1.) & | ||
347 | ! * RPI / 2. | ||
348 | ! Angle between 0 and pi | ||
349 | 3816960 | ZP(JW, II) = MOD(TT(II, JW) * 10., 1.) * RPI | |
350 | ! TEST WITH POSITIVE WAVES ONLY (Part I/II) | ||
351 | ! ZP(JW, II) = 0. | ||
352 | ! Horizontal wavenumber amplitude | ||
353 | 3816960 | ZK(JW, II) = KMIN + (KMAX - KMIN) * MOD(TT(II, JW) * 100., 1.) | |
354 | ! Horizontal phase speed | ||
355 | CPHA = 0. | ||
356 |
2/2✓ Branch 0 taken 19084800 times.
✓ Branch 1 taken 3816960 times.
|
22901760 | DO JJ = 1, NA |
357 | CPHA = CPHA + & | ||
358 | 22901760 | CMAX*2.*(MOD(TT(II, JW+4*(JJ-1)+JJ)**2, 1.)-0.5)*SQRT(3.)/SQRT(NA*1.) | |
359 | END DO | ||
360 |
2/2✓ Branch 0 taken 1908580 times.
✓ Branch 1 taken 1908380 times.
|
3816960 | IF (CPHA.LT.0.) THEN |
361 | 1908580 | CPHA = -1.*CPHA | |
362 | 1908580 | ZP(JW,II) = ZP(JW,II) + RPI | |
363 | ! TEST WITH POSITIVE WAVES ONLY (Part II/II) | ||
364 | ! ZP(JW, II) = 0. | ||
365 | ENDIF | ||
366 | 3816960 | CPHA = CPHA + CMIN !we dont allow |c|<1m/s | |
367 | ! Absolute frequency is imposed | ||
368 | 3816960 | ZO(JW, II) = CPHA * ZK(JW, II) | |
369 | ! Intrinsic frequency is imposed | ||
370 | ZO(JW, II) = ZO(JW, II) & | ||
371 | + ZK(JW, II) * COS(ZP(JW, II)) * UH(II, LAUNCH) & | ||
372 | 3816960 | + ZK(JW, II) * SIN(ZP(JW, II)) * VH(II, LAUNCH) | |
373 | ! Momentum flux at launch lev | ||
374 | ! LAUNCHED RANDOM WAVES WITH LOG-NORMAL AMPLITUDE | ||
375 | ! RIGHT IN THE SH (GWD4 after 1990) | ||
376 | 3816960 | RUW0(JW, II) = 0. | |
377 |
2/2✓ Branch 0 taken 19084800 times.
✓ Branch 1 taken 3816960 times.
|
22901760 | DO JJ = 1, NA |
378 | RUW0(JW, II) = RUW0(JW,II) + & | ||
379 | 22901760 | 2.*(MOD(TT(II, JW+4*(JJ-1)+JJ)**2, 1.)-0.5)*SQRT(3.)/SQRT(NA*1.) | |
380 | END DO | ||
381 | RUW0(JW, II) = RUWFRT & | ||
382 | * EXP(RUW0(JW,II))/1250. & ! 2 mpa at south pole | ||
383 | 3820800 | *((1.05+SIN(PLAT(II)*RPI/180.))/(1.01+SIN(PLAT(II)*RPI/180.))-2.05/2.01) | |
384 | ! RUW0(JW, II) = RUWFRT | ||
385 | ENDDO | ||
386 | end DO | ||
387 | |||
388 | ! 4. COMPUTE THE FLUXES | ||
389 | |||
390 | ! 4.0 | ||
391 | |||
392 | ! 4.1 Vertical velocity at launching altitude to ensure | ||
393 | ! the correct value to the imposed fluxes. | ||
394 | |||
395 |
2/2✓ Branch 0 taken 3840 times.
✓ Branch 1 taken 480 times.
|
4320 | DO JW = 1, NW |
396 | |||
397 | ! Evaluate intrinsic frequency at launching altitude: | ||
398 | ZOP(JW, :) = ZO(JW, :) & | ||
399 | - ZK(JW, :) * COS(ZP(JW, :)) * UH(:, LAUNCH) & | ||
400 |
2/2✓ Branch 0 taken 3816960 times.
✓ Branch 1 taken 3840 times.
|
3820800 | - ZK(JW, :) * SIN(ZP(JW, :)) * VH(:, LAUNCH) |
401 | |||
402 | ! VERSION WITH FRONTAL SOURCES | ||
403 | |||
404 | ! Momentum flux at launch level imposed by vorticity sources | ||
405 | |||
406 | ! tanh limitation for values above CORIO (inertial instability). | ||
407 | ! WWP(JW, :) = RUW0(JW, :) & | ||
408 | WWP(JW, :) = RUWFRT & | ||
409 | ! * (CORIO(:)*TANH(ROTBA(:)/CORIO(:)))**2 & | ||
410 | ! * ABS((CORIO(:)*TANH(ROTBA(:)/CORIO(:)))*CORIO(:)) & | ||
411 | ! CONSTANT FLUX | ||
412 | ! * (CORIO(:)*CORIO(:)) & | ||
413 | ! MODERATION BY THE DEPTH OF THE SOURCE (DZ HERE) | ||
414 | ! *EXP(-BVLOW(:)**2/MAX(ABS(ZOP(JW, :)),ZOISEC)**2 & | ||
415 | ! *ZK(JW, :)**2*DZ**2) & | ||
416 | ! COMPLETE FORMULA: | ||
417 | !* CORIO(:)**2*TANH(ROTBA(:)/CORIO(:)**2) & | ||
418 | * ROTBA(:) & | ||
419 | ! RESTORE DIMENSION OF A FLUX | ||
420 | ! *RD*TR/PR | ||
421 | ! *1. + RUW0(JW, :) | ||
422 |
2/2✓ Branch 0 taken 3816960 times.
✓ Branch 1 taken 3840 times.
|
3820800 | *1. |
423 | |||
424 | ! Factor related to the characteristics of the waves: NONE | ||
425 | |||
426 | ! Moderation by the depth of the source (dz here): NONE | ||
427 | |||
428 | ! Put the stress in the right direction: | ||
429 | |||
430 |
2/2✓ Branch 0 taken 3816960 times.
✓ Branch 1 taken 3840 times.
|
3820800 | RUWP(JW, :) = SIGN(1., ZOP(JW, :))*COS(ZP(JW, :)) * WWP(JW, :) |
431 |
2/2✓ Branch 0 taken 3816960 times.
✓ Branch 1 taken 3840 times.
|
3821280 | RVWP(JW, :) = SIGN(1., ZOP(JW, :))*SIN(ZP(JW, :)) * WWP(JW, :) |
432 | |||
433 | end DO | ||
434 | |||
435 | ! 4.2 Uniform values below the launching altitude | ||
436 | |||
437 |
2/2✓ Branch 0 taken 1920 times.
✓ Branch 1 taken 480 times.
|
2400 | DO LL = 1, LAUNCH |
438 |
2/2✓ Branch 0 taken 1908480 times.
✓ Branch 1 taken 1920 times.
|
1910400 | RUW(:, LL) = 0 |
439 |
2/2✓ Branch 0 taken 1908480 times.
✓ Branch 1 taken 1920 times.
|
1910400 | RVW(:, LL) = 0 |
440 |
2/2✓ Branch 0 taken 15360 times.
✓ Branch 1 taken 1920 times.
|
17760 | DO JW = 1, NW |
441 |
2/2✓ Branch 0 taken 15267840 times.
✓ Branch 1 taken 15360 times.
|
15283200 | RUW(:, LL) = RUW(:, LL) + RUWP(JW, :) |
442 |
2/2✓ Branch 0 taken 15267840 times.
✓ Branch 1 taken 15360 times.
|
15285120 | RVW(:, LL) = RVW(:, LL) + RVWP(JW, :) |
443 | end DO | ||
444 | end DO | ||
445 | |||
446 | ! 4.3 Loop over altitudes, with passage from one level to the next | ||
447 | ! done by i) conserving the EP flux, ii) dissipating a little, | ||
448 | ! iii) testing critical levels, and vi) testing the breaking. | ||
449 | |||
450 |
2/2✓ Branch 0 taken 16800 times.
✓ Branch 1 taken 480 times.
|
17280 | DO LL = LAUNCH, KLEV - 1 |
451 | ! Warning: all the physics is here (passage from one level | ||
452 | ! to the next) | ||
453 |
2/2✓ Branch 0 taken 134400 times.
✓ Branch 1 taken 16800 times.
|
151200 | DO JW = 1, NW |
454 |
2/2✓ Branch 0 taken 133593600 times.
✓ Branch 1 taken 134400 times.
|
133728000 | ZOM(JW, :) = ZOP(JW, :) |
455 |
2/2✓ Branch 0 taken 133593600 times.
✓ Branch 1 taken 134400 times.
|
133728000 | WWM(JW, :) = WWP(JW, :) |
456 | ! Intrinsic Frequency | ||
457 | ZOP(JW, :) = ZO(JW, :) - ZK(JW, :) * COS(ZP(JW, :)) * UH(:, LL + 1) & | ||
458 |
2/2✓ Branch 0 taken 133593600 times.
✓ Branch 1 taken 134400 times.
|
133728000 | - ZK(JW, :) * SIN(ZP(JW, :)) * VH(:, LL + 1) |
459 | |||
460 | ! No breaking (Eq.6) | ||
461 | ! Dissipation (Eq. 8) | ||
462 | WWP(JW, :) = WWM(JW, :) * EXP(- 4. * RDISS * PR / (PH(:, LL + 1) & | ||
463 | + PH(:, LL)) * ((BV(:, LL + 1) + BV(:, LL)) / 2.)**3 & | ||
464 | / MAX(ABS(ZOP(JW, :) + ZOM(JW, :)) / 2., ZOISEC)**4 & | ||
465 |
2/2✓ Branch 0 taken 133593600 times.
✓ Branch 1 taken 134400 times.
|
133728000 | * ZK(JW, :)**3 * (ZH(:, LL + 1) - ZH(:, LL))) |
466 | |||
467 | ! Critical levels (forced to zero if intrinsic frequency changes sign) | ||
468 | ! Saturation (Eq. 12) | ||
469 | WWP(JW, :) = min(WWP(JW, :), MAX(0., & | ||
470 | SIGN(1., ZOP(JW, :) * ZOM(JW, :))) * ABS(ZOP(JW, :))**3 & | ||
471 | ! / BV(:, LL + 1) * EXP(- ZH(:, LL + 1) / H0) * SATFRT**2 * KMIN**2 & | ||
472 | / BV(:, LL + 1) * EXP(- ZH(:, LL + 1) / H0) * KMIN**2 & | ||
473 | ! *(SATFRT*(2.5+1.5*TANH((ZH(:,LL+1)/H0-8.)/2.)))**2 & | ||
474 | *SATFRT**2 & | ||
475 |
2/2✓ Branch 0 taken 133593600 times.
✓ Branch 1 taken 134400 times.
|
133744800 | / ZK(JW, :)**4) |
476 | end DO | ||
477 | |||
478 | ! Evaluate EP-flux from Eq. 7 and give the right orientation to | ||
479 | ! the stress | ||
480 | |||
481 |
2/2✓ Branch 0 taken 134400 times.
✓ Branch 1 taken 16800 times.
|
151200 | DO JW = 1, NW |
482 |
2/2✓ Branch 0 taken 133593600 times.
✓ Branch 1 taken 134400 times.
|
133728000 | RUWP(JW, :) = SIGN(1., ZOP(JW, :))*COS(ZP(JW, :)) * WWP(JW, :) |
483 |
2/2✓ Branch 0 taken 133593600 times.
✓ Branch 1 taken 134400 times.
|
133744800 | RVWP(JW, :) = SIGN(1., ZOP(JW, :))*SIN(ZP(JW, :)) * WWP(JW, :) |
484 | end DO | ||
485 | |||
486 |
2/2✓ Branch 0 taken 16699200 times.
✓ Branch 1 taken 16800 times.
|
16716000 | RUW(:, LL + 1) = 0. |
487 |
2/2✓ Branch 0 taken 16699200 times.
✓ Branch 1 taken 16800 times.
|
16716000 | RVW(:, LL + 1) = 0. |
488 | |||
489 |
2/2✓ Branch 0 taken 16800 times.
✓ Branch 1 taken 134400 times.
|
151680 | DO JW = 1, NW |
490 |
2/2✓ Branch 0 taken 133593600 times.
✓ Branch 1 taken 134400 times.
|
133728000 | RUW(:, LL + 1) = RUW(:, LL + 1) + RUWP(JW, :) |
491 |
2/2✓ Branch 0 taken 133593600 times.
✓ Branch 1 taken 134400 times.
|
133728000 | RVW(:, LL + 1) = RVW(:, LL + 1) + RVWP(JW, :) |
492 |
2/2✓ Branch 0 taken 133593600 times.
✓ Branch 1 taken 134400 times.
|
133728000 | EAST_GWSTRESS(:, LL)=EAST_GWSTRESS(:, LL)+MAX(0.,RUWP(JW,:))/FLOAT(NW) |
493 |
2/2✓ Branch 0 taken 133593600 times.
✓ Branch 1 taken 134400 times.
|
133744800 | WEST_GWSTRESS(:, LL)=WEST_GWSTRESS(:, LL)+MIN(0.,RUWP(JW,:))/FLOAT(NW) |
494 | end DO | ||
495 | end DO | ||
496 | |||
497 | ! 5 CALCUL DES TENDANCES: | ||
498 | |||
499 | ! 5.1 Rectification des flux au sommet et dans les basses couches | ||
500 | |||
501 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | RUW(:, KLEV + 1) = 0. |
502 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | RVW(:, KLEV + 1) = 0. |
503 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | RUW(:, 1) = RUW(:, LAUNCH) |
504 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | RVW(:, 1) = RVW(:, LAUNCH) |
505 |
2/2✓ Branch 0 taken 1920 times.
✓ Branch 1 taken 480 times.
|
2400 | DO LL = 1, LAUNCH |
506 |
2/2✓ Branch 0 taken 1908480 times.
✓ Branch 1 taken 1920 times.
|
1910400 | RUW(:, LL) = RUW(:, LAUNCH+1) |
507 |
2/2✓ Branch 0 taken 1908480 times.
✓ Branch 1 taken 1920 times.
|
1910400 | RVW(:, LL) = RVW(:, LAUNCH+1) |
508 |
2/2✓ Branch 0 taken 1908480 times.
✓ Branch 1 taken 1920 times.
|
1910400 | EAST_GWSTRESS(:, LL)=EAST_GWSTRESS(:, LAUNCH) |
509 |
2/2✓ Branch 0 taken 1908480 times.
✓ Branch 1 taken 1920 times.
|
1910880 | WEST_GWSTRESS(:, LL)=WEST_GWSTRESS(:, LAUNCH) |
510 | end DO | ||
511 | |||
512 | ! AR-1 RECURSIVE FORMULA (13) IN VERSION 4 | ||
513 |
2/2✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
|
19200 | DO LL = 1, KLEV |
514 | D_U(:, LL) = (1.-DTIME/DELTAT) * D_U(:, LL) + DTIME/DELTAT/REAL(NW) * & | ||
515 | RG * (RUW(:, LL + 1) - RUW(:, LL)) & | ||
516 |
2/2✓ Branch 0 taken 18607680 times.
✓ Branch 1 taken 18720 times.
|
18626400 | / (PH(:, LL + 1) - PH(:, LL)) * DTIME |
517 | ! NO AR1 FOR MERIDIONAL TENDENCIES | ||
518 | ! D_V(:, LL) = (1.-DTIME/DELTAT) * D_V(:, LL) + DTIME/DELTAT/REAL(NW) * & | ||
519 | D_V(:, LL) = 1./REAL(NW) * & | ||
520 | RG * (RVW(:, LL + 1) - RVW(:, LL)) & | ||
521 |
2/2✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 18607680 times.
|
18626880 | / (PH(:, LL + 1) - PH(:, LL)) * DTIME |
522 | ENDDO | ||
523 | |||
524 | ! Cosmetic: evaluation of the cumulated stress | ||
525 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | ZUSTR = 0. |
526 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | ZVSTR = 0. |
527 |
2/2✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 480 times.
|
19200 | DO LL = 1, KLEV |
528 |
2/2✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 18607680 times.
|
18626880 | ZUSTR = ZUSTR + D_U(:, LL) / RG * (PH(:, LL + 1) - PH(:, LL))/DTIME |
529 | ! ZVSTR = ZVSTR + D_V(:, LL) / RG * (PH(:, LL + 1) - PH(:, LL))/DTIME | ||
530 | ENDDO | ||
531 | ! COSMETICS TO VISUALIZE ROTBA | ||
532 |
2/2✓ Branch 0 taken 477120 times.
✓ Branch 1 taken 480 times.
|
477600 | ZVSTR = ROTBA |
533 | |||
534 | 480 | END SUBROUTINE ACAMA_GWD_RANDO | |
535 | |||
536 | end module ACAMA_GWD_rando_m | ||
537 |