My Project
 All Classes Files Functions Variables Macros
integrd_p.F
Go to the documentation of this file.
1 !
2 ! $Id: integrd_p.F 1616 2012-02-17 11:59:00Z emillour $
3 !
4  SUBROUTINE integrd_p
5  $ ( nq,vcovm1,ucovm1,tetam1,psm1,massem1,
6  $ dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps0,masse,phis) !,finvmaold)
7  USE parallel
8  USE control_mod, only : planet_type
9  IMPLICIT NONE
10 
11 
12 c=======================================================================
13 c
14 c Auteur: P. Le Van
15 c -------
16 c
17 c objet:
18 c ------
19 c
20 c Incrementation des tendances dynamiques
21 c
22 c=======================================================================
23 c-----------------------------------------------------------------------
24 c Declarations:
25 c -------------
26 
27 #include "dimensions.h"
28 #include "paramet.h"
29 #include "comconst.h"
30 #include "comgeom.h"
31 #include "comvert.h"
32 #include "logic.h"
33 #include "temps.h"
34 #include "serre.h"
35 #include "iniprint.h"
36 
37 c Arguments:
38 c ----------
39 
40  integer,intent(in) :: nq ! number of tracers to handle in this routine
41  real,intent(inout) :: vcov(ip1jm,llm) ! covariant meridional wind
42  real,intent(inout) :: ucov(ip1jmp1,llm) ! covariant zonal wind
43  real,intent(inout) :: teta(ip1jmp1,llm) ! potential temperature
44  real,intent(inout) :: q(ip1jmp1,llm,nq) ! advected tracers
45  real,intent(inout) :: ps0(ip1jmp1) ! surface pressure
46  real,intent(inout) :: masse(ip1jmp1,llm) ! atmospheric mass
47  real,intent(in) :: phis(ip1jmp1) ! ground geopotential !!! unused
48  ! values at previous time step
49  real,intent(inout) :: vcovm1(ip1jm,llm)
50  real,intent(inout) :: ucovm1(ip1jmp1,llm)
51  real,intent(inout) :: tetam1(ip1jmp1,llm)
52  real,intent(inout) :: psm1(ip1jmp1)
53  real,intent(inout) :: massem1(ip1jmp1,llm)
54  ! the tendencies to add
55  real,intent(in) :: dv(ip1jm,llm)
56  real,intent(in) :: du(ip1jmp1,llm)
57  real,intent(in) :: dteta(ip1jmp1,llm)
58  real,intent(in) :: dp(ip1jmp1)
59  real,intent(in) :: dq(ip1jmp1,llm,nq) !!! unused
60 ! real,intent(out) :: finvmaold(ip1jmp1,llm) !!! unused
61 
62 c Local:
63 c ------
64 
65  REAL vscr( ip1jm ),uscr( ip1jmp1 ),hscr( ip1jmp1 ),pscr(ip1jmp1)
66  REAL massescr( ip1jmp1,llm )
67 ! REAL finvmasse(ip1jmp1,llm)
68  REAL,SAVE :: p(ip1jmp1,llmp1)
69  REAL tpn,tps,tppn(iim),tpps(iim)
70  REAL qpn,qps,qppn(iim),qpps(iim)
71  REAL,SAVE :: deltap( ip1jmp1,llm )
72 
73  INTEGER l,ij,iq,i,j
74 
75  REAL ssum
76  EXTERNAL ssum
77  INTEGER ijb,ije,jjb,jje
78  REAL,SAVE :: ps(ip1jmp1)
79  LOGICAL :: checksum
80  INTEGER :: stop_it
81 c-----------------------------------------------------------------------
82 c$OMP BARRIER
83  if (pole_nord) THEN
84 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
85  DO l = 1,llm
86  DO ij = 1,iip1
87  ucov( ij , l) = 0.
88  uscr( ij ) = 0.
89  ENDDO
90  ENDDO
91 c$OMP END DO NOWAIT
92  ENDIF
93 
94  if (pole_sud) THEN
95 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
96  DO l = 1,llm
97  DO ij = 1,iip1
98  ucov( ij +ip1jm, l) = 0.
99  uscr( ij +ip1jm ) = 0.
100  ENDDO
101  ENDDO
102 c$OMP END DO NOWAIT
103  ENDIF
104 
105 c ............ integration de ps ..............
106 
107 c CALL SCOPY(ip1jmp1*llm, masse, 1, massescr, 1)
108 
109  ijb=ij_begin
110  ije=ij_end
111 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
112  DO l = 1,llm
113  massescr(ijb:ije,l)=masse(ijb:ije,l)
114  ENDDO
115 c$OMP END DO NOWAIT
116 
117 c$OMP DO SCHEDULE(STATIC)
118  DO 2 ij = ijb,ije
119  pscr(ij) = ps0(ij)
120  ps(ij) = psm1(ij) + dt * dp(ij)
121  2 CONTINUE
122 c$OMP END DO
123 c$OMP BARRIER
124 c --> ici synchro OPENMP pour ps
125 
126  checksum=.true.
127  stop_it=0
128 
129 c$OMP DO SCHEDULE(STATIC)
130  DO ij = ijb,ije
131  IF( ps(ij).LT.0. ) THEN
132  IF (checksum) stop_it=ij
133  checksum=.false.
134  ENDIF
135  ENDDO
136 c$OMP END DO NOWAIT
137 
138  IF( .NOT. checksum ) THEN
139  write(lunout,*) "integrd: negative surface pressure ",
140  & ps(stop_it)
141  write(lunout,*) " at node ij =", stop_it
142  ! since ij=j+(i-1)*jjp1 , we have
143  j=modulo(stop_it,jjp1)
144  i=1+(stop_it-j)/jjp1
145  write(lunout,*) " lon = ",rlonv(i)*180./pi, " deg",
146  & " lat = ",rlatu(j)*180./pi, " deg"
147  ENDIF
148 
149 c
150 C$OMP MASTER
151  if (pole_nord) THEN
152 
153  DO ij = 1, iim
154  tppn(ij) = aire( ij ) * ps( ij )
155  ENDDO
156  tpn = ssum(iim,tppn,1)/apoln
157  DO ij = 1, iip1
158  ps( ij ) = tpn
159  ENDDO
160 
161  ENDIF
162 
163  if (pole_sud) THEN
164 
165  DO ij = 1, iim
166  tpps(ij) = aire(ij+ip1jm) * ps(ij+ip1jm)
167  ENDDO
168  tps = ssum(iim,tpps,1)/apols
169  DO ij = 1, iip1
170  ps(ij+ip1jm) = tps
171  ENDDO
172 
173  ENDIF
174 c$OMP END MASTER
175 c$OMP BARRIER
176 c
177 c ... Calcul de la nouvelle masse d'air au dernier temps integre t+1 ...
178 c
179 
180  CALL pression_p( ip1jmp1, ap, bp, ps, p )
181 c$OMP BARRIER
182  CALL massdair_p( p , masse )
183 
184 ! Ehouarn : we don't use/need finvmaold and finvmasse,
185 ! so might as well not compute them
186 !c CALL SCOPY( ijp1llm , masse, 1, finvmasse, 1 )
187 ! ijb=ij_begin
188 ! ije=ij_end
189 !
190 !c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
191 ! DO l = 1,llm
192 ! finvmasse(ijb:ije,l)=masse(ijb:ije,l)
193 ! ENDDO
194 !c$OMP END DO NOWAIT
195 !
196 ! jjb=jj_begin
197 ! jje=jj_end
198 ! CALL filtreg_p( finvmasse,jjb,jje, jjp1, llm, -2, 2, .TRUE., 1 )
199 c
200 
201 c ............ integration de ucov, vcov, h ..............
202 
203 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
204  DO 10 l = 1,llm
205 
206  ijb=ij_begin
207  ije=ij_end
208  if (pole_nord) ijb=ij_begin+iip1
209  if (pole_sud) ije=ij_end-iip1
210 
211  DO 4 ij = ijb,ije
212  uscr( ij ) = ucov( ij,l )
213  ucov( ij,l ) = ucovm1( ij,l ) + dt * du( ij,l )
214  4 CONTINUE
215 
216  ijb=ij_begin
217  ije=ij_end
218  if (pole_sud) ije=ij_end-iip1
219 
220  DO 5 ij = ijb,ije
221  vscr( ij ) = vcov( ij,l )
222  vcov( ij,l ) = vcovm1( ij,l ) + dt * dv( ij,l )
223  5 CONTINUE
224 
225  ijb=ij_begin
226  ije=ij_end
227 
228  DO 6 ij = ijb,ije
229  hscr( ij ) = teta(ij,l)
230  teta( ij,l ) = tetam1(ij,l) * massem1(ij,l) / masse(ij,l)
231  $ + dt * dteta(ij,l) / masse(ij,l)
232  6 CONTINUE
233 
234 c .... Calcul de la valeur moyenne, unique aux poles pour teta ......
235 c
236 c
237  IF (pole_nord) THEN
238 
239  DO ij = 1, iim
240  tppn(ij) = aire( ij ) * teta( ij ,l)
241  ENDDO
242  tpn = ssum(iim,tppn,1)/apoln
243 
244  DO ij = 1, iip1
245  teta( ij ,l) = tpn
246  ENDDO
247 
248  ENDIF
249 
250  IF (pole_sud) THEN
251 
252  DO ij = 1, iim
253  tpps(ij) = aire(ij+ip1jm) * teta(ij+ip1jm,l)
254  ENDDO
255  tps = ssum(iim,tpps,1)/apols
256 
257  DO ij = 1, iip1
258  teta(ij+ip1jm,l) = tps
259  ENDDO
260 
261  ENDIF
262 c
263 
264  IF(leapf) THEN
265 c CALL SCOPY ( ip1jmp1, uscr(1), 1, ucovm1(1, l), 1 )
266 c CALL SCOPY ( ip1jm, vscr(1), 1, vcovm1(1, l), 1 )
267 c CALL SCOPY ( ip1jmp1, hscr(1), 1, tetam1(1, l), 1 )
268  ijb=ij_begin
269  ije=ij_end
270  ucovm1(ijb:ije,l)=uscr(ijb:ije)
271  tetam1(ijb:ije,l)=hscr(ijb:ije)
272  if (pole_sud) ije=ij_end-iip1
273  vcovm1(ijb:ije,l)=vscr(ijb:ije)
274 
275  END IF
276 
277  10 CONTINUE
278 c$OMP END DO NOWAIT
279 
280 c
281 c ....... integration de q ......
282 c
283  ijb=ij_begin
284  ije=ij_end
285 
286  if (planet_type.eq."earth") then
287 ! Earth-specific treatment of first 2 tracers (water)
288 c$OMP BARRIER
289 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
290  DO l = 1, llm
291  DO ij = ijb, ije
292  deltap(ij,l) = p(ij,l) - p(ij,l+1)
293  ENDDO
294  ENDDO
295 c$OMP END DO NOWAIT
296 c$OMP BARRIER
297 
298  CALL qminimum_p( q, nq, deltap )
299 c
300 c ..... Calcul de la valeur moyenne, unique aux poles pour q .....
301 c
302 c$OMP BARRIER
303  IF (pole_nord) THEN
304 
305  DO iq = 1, nq
306 
307 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
308  DO l = 1, llm
309 
310  DO ij = 1, iim
311  qppn(ij) = aire( ij ) * q( ij ,l,iq)
312  ENDDO
313  qpn = ssum(iim,qppn,1)/apoln
314 
315  DO ij = 1, iip1
316  q( ij ,l,iq) = qpn
317  ENDDO
318 
319  ENDDO
320 c$OMP END DO NOWAIT
321 
322  ENDDO
323 
324  ENDIF
325 
326  IF (pole_sud) THEN
327 
328  DO iq = 1, nq
329 
330 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
331  DO l = 1, llm
332 
333  DO ij = 1, iim
334  qpps(ij) = aire(ij+ip1jm) * q(ij+ip1jm,l,iq)
335  ENDDO
336  qps = ssum(iim,qpps,1)/apols
337 
338  DO ij = 1, iip1
339  q(ij+ip1jm,l,iq) = qps
340  ENDDO
341 
342  ENDDO
343 c$OMP END DO NOWAIT
344 
345  ENDDO
346 
347  ENDIF
348 
349 ! Ehouarn: forget about finvmaold
350 !c CALL SCOPY( ijp1llm , finvmasse, 1, finvmaold, 1 )
351 !
352 !c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
353 ! DO l = 1, llm
354 ! finvmaold(ijb:ije,l)=finvmasse(ijb:ije,l)
355 ! ENDDO
356 !c$OMP END DO NOWAIT
357 
358  endif ! of if (planet_type.eq."earth")
359 
360 c
361 c
362 c ..... FIN de l'integration de q .......
363 
364 15 continue
365 
366 c$OMP DO SCHEDULE(STATIC)
367  DO ij=ijb,ije
368  ps0(ij)=ps(ij)
369  ENDDO
370 c$OMP END DO NOWAIT
371 
372 c .................................................................
373 
374 
375  IF( leapf ) THEN
376 c CALL SCOPY ( ip1jmp1 , pscr , 1, psm1 , 1 )
377 c CALL SCOPY ( ip1jmp1*llm, massescr, 1, massem1, 1 )
378 c$OMP DO SCHEDULE(STATIC)
379  DO ij=ijb,ije
380  psm1(ij)=pscr(ij)
381  ENDDO
382 c$OMP END DO NOWAIT
383 
384 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
385  DO l = 1, llm
386  massem1(ijb:ije,l)=massescr(ijb:ije,l)
387  ENDDO
388 c$OMP END DO NOWAIT
389  END IF
390 c$OMP BARRIER
391  RETURN
392  END