source: LMDZ6/branches/Amaury_dev/libf/dyn3dmem/integrd_loc.f90 @ 5133

Last change on this file since 5133 was 5128, checked in by abarral, 5 months ago

Correct bug in vlspltqs_loc.f90 from r2270 where we call SSUM with incorrect arguments.
Merge the three different versions of abort_gcm into one
Fix seq, para 3D compilation broken from r5107 onwards
(lint) usual + Remove uneeded fixed-form continuations

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