source: LMDZ6/branches/contrails/libf/dyn3dmem/integrd_loc.f90 @ 5445

Last change on this file since 5445 was 5285, checked in by abarral, 8 weeks ago

As discussed internally, remove generic ONLY: ... for new _mod_h modules

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