Context Navigation

advect_new_p.F @ 3281

Last change on this file since 3281 was 2622, checked in by Ehouarn Millour, 8 years ago
Some code tidying: turn ener.h into ener_mod.F90 EM
Property copyright set to Name of program: LMDZ Creation date: 1984 Version: LMDZ5 License: CeCILL version 2 Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539 See the license file in the root directory Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 7.1 KB

Line
1	!
2	! $Header$
3	!
4	SUBROUTINE advect_new_p(ucov,vcov,teta,w,massebx,masseby,
5	& du,dv,dteta)
6	USE parallel_lmdz
7	USE write_field_p
8	USE comconst_mod, ONLY: daysec
9	USE logic_mod, ONLY: conser
10
11	IMPLICIT NONE
12	c=======================================================================
13	c
14	c Auteurs: P. Le Van , Fr. Hourdin .
15	c -------
16	c
17	c Objet:
18	c ------
19	c
20	c *************************************************************
21	c .... calcul des termes d'advection vertic.pour u,v,teta,q ...
22	c *************************************************************
23	c ces termes sont ajoutes a du,dv,dteta et dq .
24	c Modif F.Forget 03/94 : on retire q de advect
25	c
26	c=======================================================================
27	c-----------------------------------------------------------------------
28	c Declarations:
29	c -------------
30
31	#include "dimensions.h"
32	#include "paramet.h"
33	#include "comgeom.h"
34
35	c Arguments:
36	c ----------
37
38	REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm),teta(ip1jmp1,llm)
39	REAL massebx(ip1jmp1,llm),masseby(ip1jm,llm),w(ip1jmp1,llm)
40	REAL dv(ip1jm,llm),du(ip1jmp1,llm),dteta(ip1jmp1,llm)
41	REAL,SAVE :: dv1(ip1jm,llm),du1(ip1jmp1,llm),dteta1(ip1jmp1,llm)
42	REAL,SAVE :: dv2(ip1jm,llm),du2(ip1jmp1,llm),dteta2(ip1jmp1,llm)
43	c Local:
44	c ------
45
46	REAL,SAVE :: uav(ip1jmp1,llm),vav(ip1jm,llm)
47	REAL wsur2(ip1jmp1)
48	REAL unsaire2(ip1jmp1), ge(ip1jmp1)
49	REAL deuxjour, ww, gt, uu, vv
50
51	INTEGER ij,l,ijb,ije
52
53	EXTERNAL SSUM
54	REAL SSUM
55
56	c-----------------------------------------------------------------------
57	c 2. Calculs preliminaires:
58	c -------------------------
59
60	IF (conser) THEN
61	deuxjour = 2. * daysec
62
63	DO 1 ij = 1, ip1jmp1
64	unsaire2(ij) = unsaire(ij) * unsaire(ij)
65	1 CONTINUE
66	END IF
67
68
69	c------------------ -yy ----------------------------------------------
70	c . Calcul de u
71
72	c$OMP MASTER
73	ijb=ij_begin
74	ije=ij_end
75	if (pole_nord) ijb=ijb+iip1
76	if (pole_sud) ije=ije-iip1
77
78	DO ij=ijb,ije
79	du2(ij,1)=0.
80	ENDDO
81
82	ijb=ij_begin
83	ije=ij_end
84	if (pole_sud) ije=ij_end-iip1
85
86	DO ij=ijb,ije
87	dv2(ij,1)=0.
88	ENDDO
89
90	ijb=ij_begin
91	ije=ij_end
92
93	DO ij=ijb,ije
94	dteta2(ij,1)=0.
95	ENDDO
96	c$OMP END MASTER
97
98
99	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
100	DO l=1,llm
101
102	ijb=ij_begin
103	ije=ij_end
104	if (pole_nord) ijb=ijb+iip1
105	if (pole_sud) ije=ije-iip1
106
107	c DO ij = iip2, ip1jmp1
108	c uav(ij,l) = 0.25 * ( ucov(ij,l) + ucov(ij-iip1,l) )
109	c ENDDO
110
111	c DO ij = iip2, ip1jm
112	c uav(ij,l) = uav(ij,l) + uav(ij+iip1,l)
113	c ENDDO
114
115	DO ij = ijb, ije
116
117	uav(ij,l)=0.25*(ucov(ij,l)+ucov(ij-iip1,l))
118	. +0.25*(ucov(ij+iip1,l)+ucov(ij,l))
119	ENDDO
120
121	if (pole_nord) then
122	DO ij = 1, iip1
123	uav(ij ,l) = 0.
124	ENDDO
125	endif
126
127	if (pole_sud) then
128	DO ij = 1, iip1
129	uav(ip1jm+ij,l) = 0.
130	ENDDO
131	endif
132
133	ENDDO
134	c$OMP END DO
135	c call write_field3d_p('uav',reshape(uav,(/iip1,jjp1,llm/)))
136
137	c------------------ -xx ----------------------------------------------
138	c . Calcul de v
139
140	ijb=ij_begin
141	ije=ij_end
142	if (pole_sud) ije=ij_end-iip1
143
144	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
145	DO l=1,llm
146
147	DO ij = ijb+1, ije
148	vav(ij,l) = 0.25 * ( vcov(ij,l) + vcov(ij-1,l) )
149	ENDDO
150
151	DO ij = ijb,ije,iip1
152	vav(ij,l) = vav(ij+iim,l)
153	ENDDO
154
155
156	DO ij = ijb, ije-1
157	vav(ij,l) = vav(ij,l) + vav(ij+1,l)
158	ENDDO
159
160	DO ij = ijb, ije, iip1
161	vav(ij+iim,l) = vav(ij,l)
162	ENDDO
163
164	ENDDO
165	c$OMP END DO
166	c call write_field3d_p('vav',reshape(vav,(/iip1,jjm,llm/)))
167
168	c-----------------------------------------------------------------------
169	c$OMP BARRIER
170
171	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
172	DO 20 l = 1, llmm1
173
174
175	c ...... calcul de - w/2. au niveau l+1 .......
176	ijb=ij_begin
177	ije=ij_end+iip1
178	if (pole_sud) ije=ij_end
179
180	DO 5 ij = ijb, ije
181	wsur2( ij ) = - 0.5 * w( ij,l+1 )
182	5 CONTINUE
183
184
185	c ..................... calcul pour du ..................
186
187	ijb=ij_begin
188	ije=ij_end
189	if (pole_nord) ijb=ijb+iip1
190	if (pole_sud) ije=ije-iip1
191
192	DO 6 ij = ijb ,ije-1
193	ww = wsur2 ( ij ) + wsur2( ij+1 )
194	uu = 0.5 * ( ucov(ij,l) + ucov(ij,l+1) )
195	du1(ij,l) = ww * ( uu - uav(ij, l ) )/massebx(ij, l )
196	du2(ij,l+1)= ww * ( uu - uav(ij,l+1) )/massebx(ij,l+1)
197	6 CONTINUE
198
199	c ................. calcul pour dv .....................
200	ijb=ij_begin
201	ije=ij_end
202	if (pole_sud) ije=ij_end-iip1
203
204	DO 8 ij = ijb, ije
205	ww = wsur2( ij+iip1 ) + wsur2( ij )
206	vv = 0.5 * ( vcov(ij,l) + vcov(ij,l+1) )
207	dv1(ij,l) = ww * (vv - vav(ij, l ) )/masseby(ij, l )
208	dv2(ij,l+1)= ww * (vv - vav(ij,l+1) )/masseby(ij,l+1)
209	8 CONTINUE
210
211	c
212
213	c ............................................................
214	c ............... calcul pour dh ...................
215	c ............................................................
216
217	c ---z
218	c calcul de - d( teta * w ) qu'on ajoute a dh
219	c ...............
220	ijb=ij_begin
221	ije=ij_end
222
223	DO 15 ij = ijb, ije
224	ww = wsur2(ij) * (teta(ij,l) + teta(ij,l+1) )
225	dteta1(ij, l ) = ww
226	dteta2(ij,l+1) = ww
227	15 CONTINUE
228
229	c ym ---> conser a voir plus tard
230
231	c IF( conser) THEN
232	c
233	c DO 17 ij = 1,ip1jmp1
234	c ge(ij) = wsur2(ij) * wsur2(ij) * unsaire2(ij)
235	c 17 CONTINUE
236	c gt = SSUM( ip1jmp1,ge,1 )
237	c gtot(l) = deuxjour * SQRT( gt/ip1jmp1 )
238	c END IF
239
240	20 CONTINUE
241	c$OMP END DO
242
243	ijb=ij_begin
244	ije=ij_end
245	if (pole_nord) ijb=ijb+iip1
246	if (pole_sud) ije=ije-iip1
247
248	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
249	DO l=1,llm
250	DO ij=ijb,ije-1
251	du(ij,l)=du(ij,l)+du2(ij,l)-du1(ij,l)
252	ENDDO
253
254	DO ij = ijb+iip1-1, ije, iip1
255	du( ij, l ) = du( ij -iim, l )
256	ENDDO
257	ENDDO
258	c$OMP END DO NOWAIT
259
260	ijb=ij_begin
261	ije=ij_end
262	if (pole_sud) ije=ij_end-iip1
263
264	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
265	DO l=1,llm
266	DO ij=ijb,ije
267	dv(ij,l)=dv(ij,l)+dv2(ij,l)-dv1(ij,l)
268	ENDDO
269	ENDDO
270	c$OMP END DO NOWAIT
271	ijb=ij_begin
272	ije=ij_end
273
274	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
275	DO l=1,llm
276	DO ij=ijb,ije
277	dteta(ij,l)=dteta(ij,l)+dteta2(ij,l)-dteta1(ij,l)
278	ENDDO
279	ENDDO
280	c$OMP END DO NOWAIT
281
282	RETURN
283	END

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

Download in other formats:

Original Format