Context Navigation

← Previous Revision
Next Revision →
Blame
Revision Log

lmdz_gwd_front.f90

Last change on this file was 6064, checked in by evignon, 5 weeks ago
suppression de print restants
File size: 19.3 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format