source: LMDZ4/branches/LMDZ4-dev/libf/dyn3dpar/integrd_p.F @ 1237

Last change on this file since 1237 was 1222, checked in by Ehouarn Millour, 15 years ago

Changes and cleanups to enable compiling without physics
and without ioipsl.

IOIPSL related cleanups:

  • bibio/writehist.F encapsulate the routine (which needs IOIPSL to function)

with #ifdef IOIPSL flag.

  • dyn3d/abort_gcm.F, dyn3dpar/abort_gcm.F and dyn3dpar/getparam.F90: use ioipsl_getincom module when not compiling with IOIPSL library, in order to always be able to use getin() routine.
  • removed unused "use IOIPSL" in dyn3dpar/guide_p_mod.F90
  • calendar related issue: Initialize day_ref and annee_ref in iniacademic.F (i.e. when they are not read from start.nc file)

Earth-specific programs/routines/modules:
create_etat0.F, fluxstokenc.F, limit_netcdf.F, startvar.F
(versions in dyn3d and dyn3dpar)
These routines and modules, which by design and porpose are made to function with
Earth physics are encapsulated with #CPP_EARTH cpp flag.

Earth-specific instructions:

  • calls to qminimum (specific treatment of first 2 tracers, i.e. water) in dyn3d/caladvtrac.F, dyn3d/integrd.F, dyn3dpar/caladvtrac_p.F, dyn3dpar/integrd_p.F only if (planet_type == 'earth')

Interaction with parallel physics:

  • routine dyn3dpar/parallel.F90 uses "surface_data" module (which is in the physics ...) to know value of "type_ocean" . Encapsulated that with #ifdef CPP_EARTH and set to a default type_ocean="dummy" otherwise.
  • So far, only Earth physics are parallelized, so all the interaction between parallel dynamics and parallel physics are encapsulated with #ifdef CCP_EARTH (this way we can run parallel without any physics). The (dyn3dpar) routines which contains such interaction are: bands.F90, gr_dyn_fi_p.F, gr_fi_dyn_p.F, mod_interface_dyn_phys.F90 This should later (when improving dyn/phys interface) be encapsulated with a more general and appropriate #ifdef CPP_PHYS cpp flag.

I checked that these changes do not alter results (on the simple
32x24x11 bench) on Ciclad (seq & mpi), Brodie (seq, mpi & omp) and
Vargas (seq, mpi & omp).

EM
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 8.4 KB
Line 
1!
2! $Id: integrd_p.F 1222 2009-08-07 11:48:33Z fairhead $
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      IMPLICIT NONE
9
10
11c=======================================================================
12c
13c   Auteur:  P. Le Van
14c   -------
15c
16c   objet:
17c   ------
18c
19c   Incrementation des tendances dynamiques
20c
21c=======================================================================
22c-----------------------------------------------------------------------
23c   Declarations:
24c   -------------
25
26#include "dimensions.h"
27#include "paramet.h"
28#include "comconst.h"
29#include "comgeom.h"
30#include "comvert.h"
31#include "logic.h"
32#include "temps.h"
33#include "serre.h"
34#include "control.h"
35
36c   Arguments:
37c   ----------
38
39      INTEGER nq
40
41      REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm),teta(ip1jmp1,llm)
42      REAL q(ip1jmp1,llm,nq)
43      REAL ps0(ip1jmp1),masse(ip1jmp1,llm),phis(ip1jmp1)
44
45      REAL vcovm1(ip1jm,llm),ucovm1(ip1jmp1,llm)
46      REAL tetam1(ip1jmp1,llm),psm1(ip1jmp1),massem1(ip1jmp1,llm)
47
48      REAL dv(ip1jm,llm),du(ip1jmp1,llm)
49      REAL dteta(ip1jmp1,llm),dp(ip1jmp1)
50      REAL dq(ip1jmp1,llm,nq), finvmaold(ip1jmp1,llm)
51
52c   Local:
53c   ------
54
55      REAL vscr( ip1jm ),uscr( ip1jmp1 ),hscr( ip1jmp1 ),pscr(ip1jmp1)
56      REAL massescr( ip1jmp1,llm ), finvmasse(ip1jmp1,llm)
57      REAL,SAVE :: p(ip1jmp1,llmp1)
58      REAL tpn,tps,tppn(iim),tpps(iim)
59      REAL qpn,qps,qppn(iim),qpps(iim)
60      REAL,SAVE :: deltap( ip1jmp1,llm )
61
62      INTEGER  l,ij,iq
63
64      REAL SSUM
65      EXTERNAL SSUM
66      INTEGER ijb,ije,jjb,jje
67      REAL,SAVE :: ps(ip1jmp1)
68      LOGICAL :: checksum
69      INTEGER :: stop_it
70c-----------------------------------------------------------------------
71c$OMP BARRIER     
72      if (pole_nord) THEN
73c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
74        DO  l = 1,llm
75          DO  ij = 1,iip1
76           ucov(    ij    , l) = 0.
77           uscr(     ij      ) = 0.
78           ENDDO
79        ENDDO
80c$OMP END DO NOWAIT       
81      ENDIF
82
83      if (pole_sud) THEN
84c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
85        DO  l = 1,llm
86          DO  ij = 1,iip1
87           ucov( ij +ip1jm, l) = 0.
88           uscr( ij +ip1jm   ) = 0.
89          ENDDO
90        ENDDO
91c$OMP END DO NOWAIT     
92      ENDIF
93
94c    ............    integration  de       ps         ..............
95
96c      CALL SCOPY(ip1jmp1*llm, masse, 1, massescr, 1)
97
98      ijb=ij_begin
99      ije=ij_end
100c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
101      DO  l = 1,llm
102        massescr(ijb:ije,l)=masse(ijb:ije,l)
103      ENDDO
104c$OMP END DO NOWAIT
105
106c$OMP DO SCHEDULE(STATIC)
107      DO 2 ij = ijb,ije
108       pscr (ij)    = ps0(ij)
109       ps (ij)      = psm1(ij) + dt * dp(ij)
110   2  CONTINUE
111c$OMP END DO 
112c$OMP BARRIER
113c --> ici synchro OPENMP pour ps
114       
115      checksum=.TRUE.
116      stop_it=0
117
118c$OMP DO SCHEDULE(STATIC)
119      DO ij = ijb,ije
120         IF( ps(ij).LT.0. ) THEN
121           IF (checksum) stop_it=ij
122           checksum=.FALSE.
123         ENDIF
124       ENDDO
125c$OMP END DO NOWAIT
126       
127        IF( .NOT. checksum ) THEN
128         PRINT*,' Au point ij = ',stop_it, ' , pression sol neg. '
129     &         , ps(stop_it)
130         STOP' dans integrd'
131        ENDIF
132
133c
134C$OMP MASTER
135      if (pole_nord) THEN
136     
137        DO  ij    = 1, iim
138         tppn(ij) = aire(   ij   ) * ps(  ij    )
139        ENDDO
140         tpn      = SSUM(iim,tppn,1)/apoln
141        DO ij   = 1, iip1
142         ps(   ij   )  = tpn
143        ENDDO
144     
145      ENDIF
146     
147      if (pole_sud) THEN
148     
149        DO  ij    = 1, iim
150         tpps(ij) = aire(ij+ip1jm) * ps(ij+ip1jm)
151        ENDDO
152         tps      = SSUM(iim,tpps,1)/apols
153        DO ij   = 1, iip1
154         ps(ij+ip1jm)  = tps
155        ENDDO
156     
157      ENDIF
158c$OMP END MASTER
159c$OMP BARRIER
160c
161c  ... Calcul  de la nouvelle masse d'air au dernier temps integre t+1 ...
162c
163
164      CALL pression_p ( ip1jmp1, ap, bp, ps, p )
165c$OMP BARRIER
166      CALL massdair_p (     p  , masse         )
167
168c      CALL   SCOPY( ijp1llm  , masse, 1, finvmasse,  1      )
169      ijb=ij_begin
170      ije=ij_end
171     
172c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
173      DO  l = 1,llm
174        finvmasse(ijb:ije,l)=masse(ijb:ije,l)
175      ENDDO
176c$OMP END DO NOWAIT
177
178      jjb=jj_begin
179      jje=jj_end
180      CALL filtreg_p( finvmasse,jjb,jje, jjp1, llm, -2, 2, .TRUE., 1  )
181c
182
183c    ............   integration  de  ucov, vcov,  h     ..............
184
185c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
186      DO 10 l = 1,llm
187     
188      ijb=ij_begin
189      ije=ij_end
190      if (pole_nord) ijb=ij_begin+iip1
191      if (pole_sud)  ije=ij_end-iip1
192     
193      DO 4 ij = ijb,ije
194      uscr( ij )   =  ucov( ij,l )
195      ucov( ij,l ) = ucovm1( ij,l ) + dt * du( ij,l )
196   4  CONTINUE
197
198      ijb=ij_begin
199      ije=ij_end
200      if (pole_sud)  ije=ij_end-iip1
201     
202      DO 5 ij = ijb,ije
203      vscr( ij )   =  vcov( ij,l )
204      vcov( ij,l ) = vcovm1( ij,l ) + dt * dv( ij,l )
205   5  CONTINUE
206     
207      ijb=ij_begin
208      ije=ij_end
209     
210      DO 6 ij = ijb,ije
211      hscr( ij )    =  teta(ij,l)
212      teta ( ij,l ) = tetam1(ij,l) *  massem1(ij,l) / masse(ij,l)
213     $                + dt * dteta(ij,l) / masse(ij,l)
214   6  CONTINUE
215
216c   ....  Calcul de la valeur moyenne, unique  aux poles pour  teta    ......
217c
218c
219      IF (pole_nord) THEN
220       
221        DO  ij   = 1, iim
222          tppn(ij) = aire(   ij   ) * teta(  ij    ,l)
223        ENDDO
224          tpn      = SSUM(iim,tppn,1)/apoln
225
226        DO ij   = 1, iip1
227          teta(   ij   ,l)  = tpn
228        ENDDO
229     
230      ENDIF
231     
232      IF (pole_sud) THEN
233       
234        DO  ij   = 1, iim
235          tpps(ij) = aire(ij+ip1jm) * teta(ij+ip1jm,l)
236        ENDDO
237          tps      = SSUM(iim,tpps,1)/apols
238
239        DO ij   = 1, iip1
240          teta(ij+ip1jm,l)  = tps
241        ENDDO
242     
243      ENDIF
244c
245
246      IF(leapf)  THEN
247c         CALL SCOPY ( ip1jmp1, uscr(1), 1, ucovm1(1, l), 1 )
248c         CALL SCOPY (   ip1jm, vscr(1), 1, vcovm1(1, l), 1 )
249c         CALL SCOPY ( ip1jmp1, hscr(1), 1, tetam1(1, l), 1 )
250        ijb=ij_begin
251        ije=ij_end
252        ucovm1(ijb:ije,l)=uscr(ijb:ije)
253        tetam1(ijb:ije,l)=hscr(ijb:ije)
254        if (pole_sud) ije=ij_end-iip1
255        vcovm1(ijb:ije,l)=vscr(ijb:ije)
256     
257      END IF
258
259  10  CONTINUE
260c$OMP END DO NOWAIT
261
262c
263c   .......  integration de   q   ......
264c
265      ijb=ij_begin
266      ije=ij_end
267
268         if (planet_type.eq."earth") then
269! Earth-specific treatment of first 2 tracers (water)
270c$OMP BARRIER
271c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
272          DO l = 1, llm
273           DO ij = ijb, ije
274            deltap(ij,l) =  p(ij,l) - p(ij,l+1)
275           ENDDO
276          ENDDO
277c$OMP END DO NOWAIT
278c$OMP BARRIER
279
280          CALL qminimum_p( q, nq, deltap )
281         endif ! of if (planet_type.eq."earth")
282c
283c    .....  Calcul de la valeur moyenne, unique  aux poles pour  q .....
284c
285c$OMP BARRIER
286      IF (pole_nord) THEN
287     
288        DO iq = 1, nq
289       
290c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
291          DO l = 1, llm
292 
293             DO ij = 1, iim
294               qppn(ij) = aire(   ij   ) * q(   ij   ,l,iq)
295             ENDDO
296               qpn  =  SSUM(iim,qppn,1)/apoln
297     
298             DO ij = 1, iip1
299               q(   ij   ,l,iq)  = qpn
300             ENDDO   
301 
302          ENDDO
303c$OMP END DO NOWAIT
304
305        ENDDO
306     
307      ENDIF
308
309      IF (pole_sud) THEN
310     
311        DO iq = 1, nq
312
313c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
314          DO l = 1, llm
315 
316             DO ij = 1, iim
317               qpps(ij) = aire(ij+ip1jm) * q(ij+ip1jm,l,iq)
318             ENDDO
319               qps  =  SSUM(iim,qpps,1)/apols
320 
321             DO ij = 1, iip1
322               q(ij+ip1jm,l,iq)  = qps
323             ENDDO   
324 
325          ENDDO
326c$OMP END DO NOWAIT
327
328        ENDDO
329     
330      ENDIF
331     
332c         CALL  SCOPY( ijp1llm , finvmasse, 1, finvmaold, 1 )
333
334c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
335      DO l = 1, llm     
336        finvmaold(ijb:ije,l)=finvmasse(ijb:ije,l)       
337      ENDDO
338c$OMP END DO NOWAIT
339c
340c
341c     .....   FIN  de l'integration  de   q    .......
342
34315    continue
344
345c$OMP DO SCHEDULE(STATIC)
346      DO ij=ijb,ije 
347        ps0(ij)=ps(ij)
348      ENDDO
349c$OMP END DO NOWAIT
350
351c    .................................................................
352
353
354      IF( leapf )  THEN
355c       CALL SCOPY (    ip1jmp1 ,  pscr   , 1,   psm1  , 1 )
356c       CALL SCOPY ( ip1jmp1*llm, massescr, 1,  massem1, 1 )
357c$OMP DO SCHEDULE(STATIC)
358      DO ij=ijb,ije 
359        psm1(ij)=pscr(ij)
360      ENDDO
361c$OMP END DO NOWAIT
362
363c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
364          DO l = 1, llm
365            massem1(ijb:ije,l)=massescr(ijb:ije,l)
366          ENDDO
367c$OMP END DO NOWAIT       
368      END IF
369c$OMP BARRIER
370      RETURN
371      END
Note: See TracBrowser for help on using the repository browser.