Context Navigation

integrd_loc.F @ 5353

Last change on this file since 5353 was 4143, checked in by dcugnet, 3 years ago
Some variables are renamed or replaced by direct equivalents: iso_indnum -> tracers(:)%iso_iName niso_possibles -> niso iqiso -> iqIsoPha ; index_trac -> itZonIso ok_iso_verif -> isoCheck ntraceurs_zone -> nzone ; ntraciso -> ntiso qperemin -> min_qparent ; masseqmin -> min_qmass ; ratiomin -> min_ratio Some renamed variables are only aliased with the older name (using USE <module>, ONLY: <oldName> => <newName>) in routines where they are repeated many times. Few hard-coded indexes are now computed (examples: ilic, iso, ivap, irneb, iq_vap, iq_liq, iso_H2O, iso_HDO, iso_HTO, iso_O17, iso_O18). The IF(isoCheck) test is now embedded in the check_isotopes_seq and check_isotopes_loc routines (lighter calling).
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
File size: 11.5 KB

Line
1	!
2	! $Id: integrd_p.F 1299 2010-01-20 14:27:21Z fairhead $
3	!
4	SUBROUTINE integrd_loc
5	$ ( nq,vcovm1,ucovm1,tetam1,psm1,massem1,
6	$ dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps0,masse,phis) !,finvmaold)
7	USE parallel_lmdz
8	USE control_mod
9	USE mod_filtreg_p
10	USE write_field_loc
11	USE write_field
12	USE integrd_mod
13	USE comconst_mod, ONLY: pi
14	USE logic_mod, ONLY: leapf
15	USE comvert_mod, ONLY: ap, bp
16	USE temps_mod, ONLY: dt
17
18	IMPLICIT NONE
19
20
21	c=======================================================================
22	c
23	c Auteur: P. Le Van
24	c -------
25	c
26	c objet:
27	c ------
28	c
29	c Incrementation des tendances dynamiques
30	c
31	c=======================================================================
32	c-----------------------------------------------------------------------
33	c Declarations:
34	c -------------
35
36	include "dimensions.h"
37	include "paramet.h"
38	include "comgeom.h"
39	include "iniprint.h"
40
41	c Arguments:
42	c ----------
43
44	INTEGER,intent(in) :: nq ! number of tracers to handle in this routine
45
46	REAL,INTENT(INOUT) :: vcov(ijb_v:ije_v,llm) ! covariant meridional wind
47	REAL,INTENT(INOUT) :: ucov(ijb_u:ije_u,llm) ! covariant zonal wind
48	REAL,INTENT(INOUT) :: teta(ijb_u:ije_u,llm) ! potential temperature
49	REAL,INTENT(INOUT) :: q(ijb_u:ije_u,llm,nq) ! advected tracers
50	REAL,INTENT(INOUT) :: ps0(ijb_u:ije_u) ! surface pressure
51	REAL,INTENT(INOUT) :: masse(ijb_u:ije_u,llm) ! atmospheric mass
52	REAL,INTENT(INOUT) :: phis(ijb_u:ije_u) ! ground geopotential !!! unused
53	! values at previous time step
54	REAL,INTENT(INOUT) :: vcovm1(ijb_v:ije_v,llm)
55	REAL,INTENT(INOUT) :: ucovm1(ijb_u:ije_u,llm)
56	REAL,INTENT(INOUT) :: tetam1(ijb_u:ije_u,llm)
57	REAL,INTENT(INOUT) :: psm1(ijb_u:ije_u)
58	REAL,INTENT(INOUT) :: massem1(ijb_u:ije_u,llm)
59	! the tendencies to add
60	REAL,INTENT(INOUT) :: dv(ijb_v:ije_v,llm)
61	REAL,INTENT(INOUT) :: du(ijb_u:ije_u,llm)
62	REAL,INTENT(INOUT) :: dteta(ijb_u:ije_u,llm)
63	REAL,INTENT(INOUT) :: dp(ijb_u:ije_u)
64	REAL,INTENT(INOUT) :: dq(ijb_u:ije_u,llm,nq) !!! unused
65	! REAL,INTENT(INOUT) ::finvmaold(ijb_u:ije_u,llm) !!! unused
66
67	c Local:
68	c ------
69
70	REAL vscr( ijb_v:ije_v ),uscr( ijb_u:ije_u )
71	REAL hscr( ijb_u:ije_u ),pscr(ijb_u:ije_u)
72	REAL massescr( ijb_u:ije_u,llm )
73	! REAL finvmasse(ijb_u:ije_u,llm)
74	REAL tpn,tps,tppn(iim),tpps(iim)
75	REAL qpn,qps,qppn(iim),qpps(iim)
76
77	INTEGER l,ij,iq,i,j
78
79	REAL SSUM
80	EXTERNAL SSUM
81	INTEGER ijb,ije,jjb,jje
82	LOGICAL :: checksum
83	LOGICAL,SAVE :: checksum_all=.TRUE.
84	INTEGER :: stop_it
85	INTEGER :: ierr
86
87	!write(,) 'integrd 88: entree, nq=',nq
88	c-----------------------------------------------------------------------
89
90	c$OMP BARRIER
91	if (pole_nord) THEN
92	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
93	DO l = 1,llm
94	DO ij = 1,iip1
95	ucov( ij , l) = 0.
96	uscr( ij ) = 0.
97	ENDDO
98	ENDDO
99	c$OMP END DO NOWAIT
100	ENDIF
101
102	if (pole_sud) THEN
103	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
104	DO l = 1,llm
105	DO ij = 1,iip1
106	ucov( ij +ip1jm, l) = 0.
107	uscr( ij +ip1jm ) = 0.
108	ENDDO
109	ENDDO
110	c$OMP END DO NOWAIT
111	ENDIF
112
113	c ............ integration de ps ..............
114
115	c CALL SCOPY(ip1jmp1*llm, masse, 1, massescr, 1)
116
117	ijb=ij_begin
118	ije=ij_end
119	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
120	DO l = 1,llm
121	massescr(ijb:ije,l)=masse(ijb:ije,l)
122	ENDDO
123	c$OMP END DO NOWAIT
124
125	c$OMP DO SCHEDULE(STATIC)
126	DO 2 ij = ijb,ije
127	pscr (ij) = ps0(ij)
128	ps (ij) = psm1(ij) + dt * dp(ij)
129
130	2 CONTINUE
131
132	c$OMP END DO
133	c$OMP BARRIER
134	c --> ici synchro OPENMP pour ps
135
136	checksum=.TRUE.
137	stop_it=0
138
139	c$OMP MASTER
140	!c$OMP DO SCHEDULE(STATIC)
141	DO ij = ijb,ije
142	IF( ps(ij).LT.0. ) THEN
143	IF (checksum) stop_it=ij
144	checksum=.FALSE.
145	ENDIF
146	ENDDO
147	!c$OMP END DO NOWAIT
148
149	! CALL MPI_ALLREDUCE(checksum,checksum_all,1,
150	! & MPI_LOGICAL,MPI_LOR,COMM_LMDZ,ierr)
151	IF( .NOT. checksum ) THEN
152	write(lunout,*) "integrd: ps = ", ps(stop_it)
153	write(lunout,*) " at node ij =", stop_it
154	! since ij=j+(i-1)*jjp1 , we have
155	j=modulo(stop_it,jjp1)
156	i=1+(stop_it-j)/jjp1
157	write(lunout,) " lon = ",rlonv(i)180./pi, " deg",
158	& " lat = ",rlatu(j)*180./pi, " deg"
159	call abort_gcm("integrd_loc", "negative surface pressure", 1)
160	ENDIF
161
162	c$OMP END MASTER
163	c$OMP BARRIER
164	!write(,) 'integrd 170'
165	IF (.NOT. Checksum_all) THEN
166	call WriteField_v('int_vcov',vcov)
167	call WriteField_u('int_ucov',ucov)
168	call WriteField_u('int_teta',teta)
169	call WriteField_u('int_ps0',ps0)
170	call WriteField_u('int_masse',masse)
171	call WriteField_u('int_phis',phis)
172	call WriteField_v('int_vcovm1',vcovm1)
173	call WriteField_u('int_ucovm1',ucovm1)
174	call WriteField_u('int_tetam1',tetam1)
175	call WriteField_u('int_psm1',psm1)
176	call WriteField_u('int_massem1',massem1)
177
178	call WriteField_v('int_dv',dv)
179	call WriteField_u('int_du',du)
180	call WriteField_u('int_dteta',dteta)
181	call WriteField_u('int_dp',dp)
182	! call WriteField_u('int_finvmaold',finvmaold)
183	do j=1,nq
184	call WriteField_u('int_q'//trim(int2str(j)),
185	. q(:,:,j))
186	call WriteField_u('int_dq'//trim(int2str(j)),
187	. dq(:,:,j))
188	enddo
189	call abort_gcm("integrd_loc", "", 1)
190	ENDIF
191
192
193	c
194	!write(,) 'integrd 200'
195	C$OMP MASTER
196	if (pole_nord) THEN
197
198	DO ij = 1, iim
199	tppn(ij) = aire( ij ) * ps( ij )
200	ENDDO
201	tpn = SSUM(iim,tppn,1)/apoln
202	DO ij = 1, iip1
203	ps( ij ) = tpn
204	ENDDO
205
206	ENDIF
207
208	if (pole_sud) THEN
209
210	DO ij = 1, iim
211	tpps(ij) = aire(ij+ip1jm) * ps(ij+ip1jm)
212	ENDDO
213	tps = SSUM(iim,tpps,1)/apols
214	DO ij = 1, iip1
215	ps(ij+ip1jm) = tps
216	ENDDO
217
218	ENDIF
219	c$OMP END MASTER
220	c$OMP BARRIER
221	!write(,) 'integrd 217'
222	c
223	c ... Calcul de la nouvelle masse d'air au dernier temps integre t+1 ...
224	c
225
226	CALL pression_loc ( ip1jmp1, ap, bp, ps, p )
227
228	c$OMP BARRIER
229	CALL massdair_loc ( p , masse )
230
231	! Ehouarn : we don't use/need finvmaold and finvmasse,
232	! so might as well not compute them
233	!c CALL SCOPY( ijp1llm , masse, 1, finvmasse, 1 )
234	! ijb=ij_begin
235	! ije=ij_end
236	!
237	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
238	! DO l = 1,llm
239	! finvmasse(ijb:ije,l)=masse(ijb:ije,l)
240	! ENDDO
241	!c$OMP END DO NOWAIT
242
243	! jjb=jj_begin
244	! jje=jj_end
245	! CALL filtreg_p( finvmasse,jjb_u,jje_u,jjb,jje, jjp1, llm,
246	! & -2, 2, .TRUE., 1 )
247	c
248
249	c ............ integration de ucov, vcov, h ..............
250
251	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
252	DO 10 l = 1,llm
253
254	ijb=ij_begin
255	ije=ij_end
256	if (pole_nord) ijb=ij_begin+iip1
257	if (pole_sud) ije=ij_end-iip1
258
259	DO 4 ij = ijb,ije
260	uscr( ij ) = ucov( ij,l )
261	ucov( ij,l ) = ucovm1( ij,l ) + dt * du( ij,l )
262	4 CONTINUE
263
264	ijb=ij_begin
265	ije=ij_end
266	if (pole_sud) ije=ij_end-iip1
267
268	DO 5 ij = ijb,ije
269	vscr( ij ) = vcov( ij,l )
270	vcov( ij,l ) = vcovm1( ij,l ) + dt * dv( ij,l )
271	5 CONTINUE
272
273	ijb=ij_begin
274	ije=ij_end
275
276	DO 6 ij = ijb,ije
277	hscr( ij ) = teta(ij,l)
278	teta ( ij,l ) = tetam1(ij,l) * massem1(ij,l) / masse(ij,l)
279	$ + dt * dteta(ij,l) / masse(ij,l)
280	6 CONTINUE
281
282	c .... Calcul de la valeur moyenne, unique aux poles pour teta ......
283	c
284	c
285	!write(,) 'integrd 291'
286	IF (pole_nord) THEN
287
288	DO ij = 1, iim
289	tppn(ij) = aire( ij ) * teta( ij ,l)
290	ENDDO
291	tpn = SSUM(iim,tppn,1)/apoln
292
293	DO ij = 1, iip1
294	teta( ij ,l) = tpn
295	ENDDO
296
297	ENDIF
298
299	IF (pole_sud) THEN
300
301	DO ij = 1, iim
302	tpps(ij) = aire(ij+ip1jm) * teta(ij+ip1jm,l)
303	ENDDO
304	tps = SSUM(iim,tpps,1)/apols
305
306	DO ij = 1, iip1
307	teta(ij+ip1jm,l) = tps
308	ENDDO
309
310	ENDIF
311	c
312
313	IF(leapf) THEN
314	c CALL SCOPY ( ip1jmp1, uscr(1), 1, ucovm1(1, l), 1 )
315	c CALL SCOPY ( ip1jm, vscr(1), 1, vcovm1(1, l), 1 )
316	c CALL SCOPY ( ip1jmp1, hscr(1), 1, tetam1(1, l), 1 )
317	ijb=ij_begin
318	ije=ij_end
319	ucovm1(ijb:ije,l)=uscr(ijb:ije)
320	tetam1(ijb:ije,l)=hscr(ijb:ije)
321	if (pole_sud) ije=ij_end-iip1
322	vcovm1(ijb:ije,l)=vscr(ijb:ije)
323
324	END IF
325
326	10 CONTINUE
327	c$OMP END DO NOWAIT
328
329	c
330	c ....... integration de q ......
331	c
332	ijb=ij_begin
333	ije=ij_end
334
335	if (planet_type.eq."earth") then
336	! Earth-specific treatment of first 2 tracers (water)
337	c$OMP BARRIER
338	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
339	DO l = 1, llm
340	DO ij = ijb, ije
341	deltap(ij,l) = p(ij,l) - p(ij,l+1)
342	ENDDO
343	ENDDO
344
345	c$OMP END DO NOWAIT
346	c$OMP BARRIER
347
348	call check_isotopes(q,ijb,ije,'integrd 342')
349
350	!write(,) 'integrd 341'
351	CALL qminimum_loc( q, nq, deltap )
352	!write(,) 'integrd 343'
353
354	call check_isotopes(q,ijb,ije,'integrd 346')
355	c
356	c ..... Calcul de la valeur moyenne, unique aux poles pour q .....
357	c
358	c$OMP BARRIER
359	IF (pole_nord) THEN
360
361	DO iq = 1, nq
362
363	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
364	DO l = 1, llm
365
366	DO ij = 1, iim
367	qppn(ij) = aire( ij ) * q( ij ,l,iq)
368	ENDDO
369	qpn = SSUM(iim,qppn,1)/apoln
370
371	DO ij = 1, iip1
372	q( ij ,l,iq) = qpn
373	ENDDO
374
375	ENDDO
376	c$OMP END DO NOWAIT
377
378	ENDDO
379
380	ENDIF
381
382	IF (pole_sud) THEN
383
384	DO iq = 1, nq
385
386	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
387	DO l = 1, llm
388
389	DO ij = 1, iim
390	qpps(ij) = aire(ij+ip1jm) * q(ij+ip1jm,l,iq)
391	ENDDO
392	qps = SSUM(iim,qpps,1)/apols
393
394	DO ij = 1, iip1
395	q(ij+ip1jm,l,iq) = qps
396	ENDDO
397
398	ENDDO
399	c$OMP END DO NOWAIT
400
401	ENDDO
402
403	ENDIF
404
405	call check_isotopes(q,ijb,ije,'integrd 409')
406
407	! Ehouarn: forget about finvmaold
408	!c CALL SCOPY( ijp1llm , finvmasse, 1, finvmaold, 1 )
409
410	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
411	! DO l = 1, llm
412	! finvmaold(ijb:ije,l)=finvmasse(ijb:ije,l)
413	! ENDDO
414	!c$OMP END DO NOWAIT
415
416	endif ! of if (planet_type.eq."earth")
417
418	c
419	c
420	c ..... FIN de l'integration de q .......
421
422	15 continue
423	!write(,) 'integrd 410'
424
425	c$OMP DO SCHEDULE(STATIC)
426	DO ij=ijb,ije
427	ps0(ij)=ps(ij)
428	ENDDO
429	c$OMP END DO NOWAIT
430
431	c .................................................................
432
433
434	IF( leapf ) THEN
435	c CALL SCOPY ( ip1jmp1 , pscr , 1, psm1 , 1 )
436	c CALL SCOPY ( ip1jmp1*llm, massescr, 1, massem1, 1 )
437	c$OMP DO SCHEDULE(STATIC)
438	DO ij=ijb,ije
439	psm1(ij)=pscr(ij)
440	ENDDO
441	c$OMP END DO NOWAIT
442
443	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
444	DO l = 1, llm
445	massem1(ijb:ije,l)=massescr(ijb:ije,l)
446	ENDDO
447	c$OMP END DO NOWAIT
448	END IF
449	c$OMP BARRIER
450	RETURN
451	END

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

Download in other formats:

Original Format