source: LMDZ4/trunk/libf/dyn3dpar/calfis_p.F @ 1076

Last change on this file since 1076 was 1000, checked in by Laurent Fairhead, 16 years ago
  • Modifs sur le parallelisme: masquage dans la physique
  • Inclusion strato
  • mise en coherence etat0
  • le mode offline fonctionne maintenant en parallele,
  • les fichiers de la dynamiques sont correctement sortis et peuvent etre reconstruit avec rebuild
  • la version parallele de la dynamique peut s'executer sans MPI (sur 1 proc)
  • L'OPENMP fonctionne maintenant sans la parallelisation MPI.

YM
LF

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 28.1 KB
Line 
1!
2! $Header$
3!
4C
5C
6      SUBROUTINE calfis_p(nq,
7     $                  lafin,
8     $                  rdayvrai,
9     $                  heure,
10     $                  pucov,
11     $                  pvcov,
12     $                  pteta,
13     $                  pq,
14     $                  pmasse,
15     $                  pps,
16     $                  pp,
17     $                  ppk,
18     $                  pphis,
19     $                  pphi,
20     $                  pducov,
21     $                  pdvcov,
22     $                  pdteta,
23     $                  pdq,
24     $                  flxw,
25     $                  clesphy0,
26     $                  pdufi,
27     $                  pdvfi,
28     $                  pdhfi,
29     $                  pdqfi,
30     $                  pdpsfi)
31c
32c    Auteur :  P. Le Van, F. Hourdin
33c   .........
34      USE dimphy
35      USE mod_phys_lmdz_para, mpi_root_xx=>mpi_root
36      USE parallel, ONLY : omp_chunk, using_mpi
37      USE mod_interface_dyn_phys
38      USE Write_Field
39      Use Write_field_p
40      USE Times
41      USE IOPHY
42      IMPLICIT NONE
43c=======================================================================
44c
45c   1. rearrangement des tableaux et transformation
46c      variables dynamiques  >  variables physiques
47c   2. calcul des termes physiques
48c   3. retransformation des tendances physiques en tendances dynamiques
49c
50c   remarques:
51c   ----------
52c
53c    - les vents sont donnes dans la physique par leurs composantes
54c      naturelles.
55c    - la variable thermodynamique de la physique est une variable
56c      intensive :   T
57c      pour la dynamique on prend    T * ( preff / p(l) ) **kappa
58c    - les deux seules variables dependant de la geometrie necessaires
59c      pour la physique sont la latitude pour le rayonnement et
60c      l'aire de la maille quand on veut integrer une grandeur
61c      horizontalement.
62c    - les points de la physique sont les points scalaires de la
63c      la dynamique; numerotation:
64c          1 pour le pole nord
65c          (jjm-1)*iim pour l'interieur du domaine
66c          ngridmx pour le pole sud
67c      ---> ngridmx=2+(jjm-1)*iim
68c
69c     Input :
70c     -------
71c       ecritphy        frequence d'ecriture (en jours)de histphy
72c       pucov           covariant zonal velocity
73c       pvcov           covariant meridional velocity
74c       pteta           potential temperature
75c       pps             surface pressure
76c       pmasse          masse d'air dans chaque maille
77c       pts             surface temperature  (K)
78c       callrad         clef d'appel au rayonnement
79c
80c    Output :
81c    --------
82c        pdufi          tendency for the natural zonal velocity (ms-1)
83c        pdvfi          tendency for the natural meridional velocity
84c        pdhfi          tendency for the potential temperature
85c        pdtsfi         tendency for the surface temperature
86c
87c        pdtrad         radiative tendencies  \  both input
88c        pfluxrad       radiative fluxes      /  and output
89c
90c=======================================================================
91c
92c-----------------------------------------------------------------------
93c
94c    0.  Declarations :
95c    ------------------
96
97#include "dimensions.h"
98#include "paramet.h"
99#include "temps.h"
100#include "advtrac.h"
101
102      INTEGER ngridmx,nq
103      PARAMETER( ngridmx = 2+(jjm-1)*iim - 1/jjm   )
104
105#include "comconst.h"
106#include "comvert.h"
107#include "comgeom2.h"
108#include "control.h"
109#ifdef CPP_MPI
110      include 'mpif.h'
111#endif
112c    Arguments :
113c    -----------
114      LOGICAL  lafin
115      REAL heure
116
117      REAL pvcov(iip1,jjm,llm)
118      REAL pucov(iip1,jjp1,llm)
119      REAL pteta(iip1,jjp1,llm)
120      REAL pmasse(iip1,jjp1,llm)
121      REAL pq(iip1,jjp1,llm,nqmx)
122      REAL pphis(iip1,jjp1)
123      REAL pphi(iip1,jjp1,llm)
124c
125      REAL pdvcov(iip1,jjm,llm)
126      REAL pducov(iip1,jjp1,llm)
127      REAL pdteta(iip1,jjp1,llm)
128      REAL pdq(iip1,jjp1,llm,nqmx)
129c
130      REAL pps(iip1,jjp1)
131      REAL pp(iip1,jjp1,llmp1)
132      REAL ppk(iip1,jjp1,llm)
133c
134      REAL pdvfi(iip1,jjm,llm)
135      REAL pdufi(iip1,jjp1,llm)
136      REAL pdhfi(iip1,jjp1,llm)
137      REAL pdqfi(iip1,jjp1,llm,nqmx)
138      REAL pdpsfi(iip1,jjp1)
139
140      INTEGER        longcles
141      PARAMETER    ( longcles = 20 )
142      REAL clesphy0( longcles )
143
144
145c    Local variables :
146c    -----------------
147
148      INTEGER i,j,l,ig0,ig,iq,iiq
149      REAL,ALLOCATABLE,SAVE :: zpsrf(:)
150      REAL,ALLOCATABLE,SAVE :: zplev(:,:),zplay(:,:)
151      REAL,ALLOCATABLE,SAVE :: zphi(:,:),zphis(:)
152c
153      REAL,ALLOCATABLE,SAVE :: zufi(:,:), zvfi(:,:)
154      REAL,ALLOCATABLE,SAVE :: ztfi(:,:),zqfi(:,:,:)
155c
156      REAL,ALLOCATABLE,SAVE :: pcvgu(:,:), pcvgv(:,:)
157      REAL,ALLOCATABLE,SAVE :: pcvgt(:,:), pcvgq(:,:,:)
158c
159c      REAL,ALLOCATABLE,SAVE :: pvervel(:,:)
160c
161      REAL,ALLOCATABLE,SAVE :: zdufi(:,:),zdvfi(:,:)
162      REAL,ALLOCATABLE,SAVE :: zdtfi(:,:),zdqfi(:,:,:)
163      REAL,ALLOCATABLE,SAVE :: zdpsrf(:)
164      REAL,SAVE,ALLOCATABLE ::  flxwfi(:,:)     ! Flux de masse verticale sur la grille physiq
165
166c
167      REAL,ALLOCATABLE,SAVE :: zplev_omp(:,:)
168      REAL,ALLOCATABLE,SAVE :: zplay_omp(:,:)
169      REAL,ALLOCATABLE,SAVE :: zphi_omp(:,:)
170      REAL,ALLOCATABLE,SAVE :: zphis_omp(:)
171      REAL,ALLOCATABLE,SAVE :: presnivs_omp(:)
172      REAL,ALLOCATABLE,SAVE :: zufi_omp(:,:)
173      REAL,ALLOCATABLE,SAVE :: zvfi_omp(:,:)
174      REAL,ALLOCATABLE,SAVE :: ztfi_omp(:,:)
175      REAL,ALLOCATABLE,SAVE :: zqfi_omp(:,:,:)
176c      REAL,ALLOCATABLE,SAVE :: pvervel_omp(:,:)
177      REAL,ALLOCATABLE,SAVE :: zdufi_omp(:,:)
178      REAL,ALLOCATABLE,SAVE :: zdvfi_omp(:,:)
179      REAL,ALLOCATABLE,SAVE :: zdtfi_omp(:,:)
180      REAL,ALLOCATABLE,SAVE :: zdqfi_omp(:,:,:)
181      REAL,ALLOCATABLE,SAVE :: zdpsrf_omp(:)
182      REAL,SAVE,ALLOCATABLE ::  flxwfi_omp(:,:)     ! Flux de masse verticale sur la grille physiq
183
184c$OMP THREADPRIVATE(zplev_omp,zplay_omp,zphi_omp,zphis_omp,
185c$OMP+                 presnivs_omp,zufi_omp,zvfi_omp,ztfi_omp,
186c$OMP+                 zqfi_omp,zdufi_omp,zdvfi_omp,
187c$OMP+                 zdtfi_omp,zdqfi_omp,zdpsrf_omp,flxwfi_omp)       
188
189      LOGICAL,SAVE :: first_omp=.true.
190c$OMP THREADPRIVATE(first_omp)
191     
192      REAL zsin(iim),zcos(iim),z1(iim)
193      REAL zsinbis(iim),zcosbis(iim),z1bis(iim)
194      REAL unskap, pksurcp
195c
196cIM diagnostique PVteta, Amip2
197      INTEGER ntetaSTD
198      PARAMETER(ntetaSTD=3)
199      REAL rtetaSTD(ntetaSTD)
200      DATA rtetaSTD/350., 380., 405./
201      REAL PVteta(klon,ntetaSTD)
202     
203      REAL flxw(iip1,jjp1,llm)  ! Flux de masse verticale sur la grille dynamique
204     
205      REAL SSUM
206
207      LOGICAL firstcal, debut
208      DATA firstcal/.true./
209      SAVE firstcal,debut
210c$OMP THREADPRIVATE(firstcal,debut)
211      REAL rdayvrai
212     
213      REAL,SAVE,dimension(1:iim,1:llm):: du_send,du_recv,dv_send,dv_recv
214      INTEGER :: ierr
215#ifdef CPP_MPI
216      INTEGER,dimension(MPI_STATUS_SIZE,4) :: Status
217#else
218      INTEGER,dimension(1,4) :: Status
219#endif
220      INTEGER, dimension(4) :: Req
221      REAL,ALLOCATABLE,SAVE:: zdufi2(:,:),zdvfi2(:,:)
222      integer :: k,kstart,kend
223      INTEGER :: offset 
224c
225c-----------------------------------------------------------------------
226c
227c    1. Initialisations :
228c    --------------------
229c
230
231      klon=klon_mpi
232     
233      PVteta(:,:)=0.
234           
235      IF (ngridmx.NE.2+(jjm-1)*iim) THEN
236         PRINT*,'STOP dans calfis'
237         PRINT*,'La dimension ngridmx doit etre egale a 2 + (jjm-1)*iim'
238         PRINT*,'  ngridmx  jjm   iim   '
239         PRINT*,ngridmx,jjm,iim
240         STOP
241      ENDIF
242
243c-----------------------------------------------------------------------
244c   latitude, longitude et aires des mailles pour la physique:
245c   ----------------------------------------------------------
246
247c
248      IF ( firstcal )  THEN
249          debut = .TRUE.
250c$OMP MASTER
251      ALLOCATE(zpsrf(klon))
252      ALLOCATE(zplev(klon,llm+1),zplay(klon,llm))
253      ALLOCATE(zphi(klon,llm),zphis(klon))
254      ALLOCATE(zufi(klon,llm), zvfi(klon,llm))
255      ALLOCATE(ztfi(klon,llm),zqfi(klon,llm,nqmx))
256      ALLOCATE(pcvgu(klon,llm), pcvgv(klon,llm))
257      ALLOCATE(pcvgt(klon,llm), pcvgq(klon,llm,2))
258c      ALLOCATE(pvervel(klon,llm))
259      ALLOCATE(zdufi(klon,llm),zdvfi(klon,llm))
260      ALLOCATE(zdtfi(klon,llm),zdqfi(klon,llm,nqmx))
261      ALLOCATE(zdpsrf(klon))
262      ALLOCATE(zdufi2(klon+iim,llm),zdvfi2(klon+iim,llm))
263      ALLOCATE(flxwfi(klon,llm))
264c$OMP END MASTER
265c$OMP BARRIER     
266      ELSE
267          debut = .FALSE.
268      ENDIF
269
270c
271c
272c-----------------------------------------------------------------------
273c   40. transformation des variables dynamiques en variables physiques:
274c   ---------------------------------------------------------------
275
276c   41. pressions au sol (en Pascals)
277c   ----------------------------------
278
279c$OMP MASTER
280      call start_timer(timer_physic)
281c$OMP END MASTER
282
283c$OMP MASTER             
284      do ig0=1,klon
285        i=index_i(ig0)
286        j=index_j(ig0)
287        zpsrf(ig0)=pps(i,j)
288      enddo
289c$OMP END MASTER
290
291
292c   42. pression intercouches :
293c
294c   -----------------------------------------------------------------
295c     .... zplev  definis aux (llm +1) interfaces des couches  ....
296c     .... zplay  definis aux (  llm )    milieux des couches  ....
297c   -----------------------------------------------------------------
298
299c    ...    Exner = cp * ( p(l) / preff ) ** kappa     ....
300c
301       unskap   = 1./ kappa
302c
303c      print *,omp_rank,'klon--->',klon
304c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
305      DO l = 1, llmp1
306        do ig0=1,klon
307          i=index_i(ig0)
308          j=index_j(ig0)
309          zplev( ig0,l ) = pp(i,j,l)
310        enddo
311      ENDDO
312c$OMP END DO NOWAIT
313c
314c
315
316c   43. temperature naturelle (en K) et pressions milieux couches .
317c   ---------------------------------------------------------------
318c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
319      DO l=1,llm
320
321        do ig0=1,klon
322          i=index_i(ig0)
323          j=index_j(ig0)
324          pksurcp        = ppk(i,j,l) / cpp
325          zplay(ig0,l)   = preff * pksurcp ** unskap
326          ztfi(ig0,l)    = pteta(i,j,l)  * pksurcp
327c          pcvgt(ig0,l)   = pdteta(i,j,l) * pksurcp / pmasse(i,j,l)
328        enddo
329
330      ENDDO
331c$OMP END DO NOWAIT
332
333c   43.bis traceurs
334c   ---------------
335c
336
337      DO iq=1,nq
338         iiq=niadv(iq)
339c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
340         DO l=1,llm
341           do ig0=1,klon
342             i=index_i(ig0)
343             j=index_j(ig0)
344             zqfi(ig0,l,iq)  = pq(i,j,l,iiq)
345           enddo
346         ENDDO
347c$OMP END DO NOWAIT     
348      ENDDO
349
350c   convergence dynamique pour les traceurs "EAU"
351
352      DO iq=1,2
353c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
354         DO l=1,llm
355           do ig0=1,klon
356             i=index_i(ig0)
357             j=index_j(ig0)
358c             pcvgq(ig0,l,iq) = pdq(i,j,l,iq) / pmasse(i,j,l)
359           enddo
360         ENDDO
361c$OMP END DO NOWAIT     
362      ENDDO
363
364
365
366c   Geopotentiel calcule par rapport a la surface locale:
367c   -----------------------------------------------------
368
369      CALL gr_dyn_fi_p(llm,iip1,jjp1,klon,pphi,zphi)
370
371c$OMP MASTER
372      CALL gr_dyn_fi_p(1,iip1,jjp1,klon,pphis,zphis)
373c$OMP END MASTER
374c$OMP BARRIER
375
376c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
377      DO l=1,llm
378         DO ig=1,klon
379           zphi(ig,l)=zphi(ig,l)-zphis(ig)
380         ENDDO
381      ENDDO
382c$OMP END DO NOWAIT
383     
384c   ....  Calcul de la vitesse  verticale  ( en Pa*m*s  ou Kg/s )  ....
385c JG : ancien calcule de omega utilise dans physiq.F. Maintenant le flux
386c      de masse est calclue dans advtrac_p.F 
387c
388cc$OMP DO SCHEDULE(STATIC,OMP_CHUNK)       
389c      DO l=1,llm
390c        do ig0=1,klon
391c           i=index_i(ig0)
392c           j=index_j(ig0)
393c           pvervel(ig0,l) = pw(i,j,l)*g* unsaire(i,j)
394c        enddo
395c       if (is_north_pole) pvervel(1,l)=pw(1,1,l)*g /apoln
396c       if (is_south_pole) pvervel(klon,l)=pw(1,jjp1,l)*g/apols
397c      ENDDO
398cc$OMP END DO NOWAIT
399
400c
401c   45. champ u:
402c   ------------
403
404      kstart=1
405      kend=klon
406     
407      if (is_north_pole) kstart=2
408      if (is_south_pole) kend=klon-1
409     
410c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
411      DO l=1,llm
412        do ig0=kstart,kend
413          i=index_i(ig0)
414          j=index_j(ig0)
415          if (i==1) then
416            zufi(ig0,l)= 0.5 *(  pucov(iim,j,l)/cu(iim,j)
417     $                         + pucov(1,j,l)/cu(1,j) )
418c            pcvgu(ig0,l)= 0.5*(  pducov(iim,j,l)/cu(iim,j)
419c     $                         + pducov(1,j,l)/cu(1,j) )
420          else
421            zufi(ig0,l)= 0.5*(  pucov(i-1,j,l)/cu(i-1,j)
422     $                       + pucov(i,j,l)/cu(i,j) )
423c            pcvgu(ig0,l)= 0.5*(  pducov(i-1,j,l)/cu(i-1,j)
424c     $                        + pducov(i,j,l)/cu(i,j) )
425          endif
426        enddo
427      ENDDO
428c$OMP END DO NOWAIT
429c   46.champ v:
430c   -----------
431c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
432      DO l=1,llm
433        DO ig0=kstart,kend
434          i=index_i(ig0)
435          j=index_j(ig0)
436          zvfi(ig0,l)= 0.5 *(  pvcov(i,j-1,l)/cv(i,j-1)
437     $                       + pvcov(i,j,l)/cv(i,j) )
438   
439c          pcvgv(ig0+i,l)= 0.5 * (  pdvcov(i,j-1,l)/cv(i,j-1)
440c     $                           + pdvcov(i,j,l)/cv(i,j) )
441         ENDDO
442      ENDDO
443c$OMP END DO NOWAIT
444
445c   47. champs de vents aux pole nord   
446c   ------------------------------
447c        U = 1 / pi  *  integrale [ v * cos(long) * d long ]
448c        V = 1 / pi  *  integrale [ v * sin(long) * d long ]
449
450      if (is_north_pole) then
451c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
452        DO l=1,llm
453
454           z1(1)   =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,1,l)/cv(1,1)
455c           z1bis(1)=(rlonu(1)-rlonu(iim)+2.*pi)*pdvcov(1,1,l)/cv(1,1)
456           DO i=2,iim
457              z1(i)   =(rlonu(i)-rlonu(i-1))*pvcov(i,1,l)/cv(i,1)
458c              z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,1,l)/cv(i,1)
459           ENDDO
460 
461           DO i=1,iim
462              zcos(i)   = COS(rlonv(i))*z1(i)
463c              zcosbis(i)= COS(rlonv(i))*z1bis(i)
464              zsin(i)   = SIN(rlonv(i))*z1(i)
465c              zsinbis(i)= SIN(rlonv(i))*z1bis(i)
466           ENDDO
467 
468           zufi(1,l)  = SSUM(iim,zcos,1)/pi
469c           pcvgu(1,l) = SSUM(iim,zcosbis,1)/pi
470           zvfi(1,l)  = SSUM(iim,zsin,1)/pi
471c           pcvgv(1,l) = SSUM(iim,zsinbis,1)/pi
472 
473        ENDDO
474c$OMP END DO NOWAIT     
475      endif
476
477
478c   48. champs de vents aux pole sud:
479c   ---------------------------------
480c        U = 1 / pi  *  integrale [ v * cos(long) * d long ]
481c        V = 1 / pi  *  integrale [ v * sin(long) * d long ]
482
483      if (is_south_pole) then
484c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
485        DO l=1,llm
486 
487         z1(1)   =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,jjm,l)/cv(1,jjm)
488c         z1bis(1)=(rlonu(1)-rlonu(iim)+2.*pi)*pdvcov(1,jjm,l)/cv(1,jjm)
489           DO i=2,iim
490           z1(i)   =(rlonu(i)-rlonu(i-1))*pvcov(i,jjm,l)/cv(i,jjm)
491c           z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,jjm,l)/cv(i,jjm)
492           ENDDO
493 
494           DO i=1,iim
495              zcos(i)    = COS(rlonv(i))*z1(i)
496c              zcosbis(i) = COS(rlonv(i))*z1bis(i)
497              zsin(i)    = SIN(rlonv(i))*z1(i)
498c              zsinbis(i) = SIN(rlonv(i))*z1bis(i)
499           ENDDO
500 
501           zufi(klon,l)  = SSUM(iim,zcos,1)/pi
502c           pcvgu(klon,l) = SSUM(iim,zcosbis,1)/pi
503           zvfi(klon,l)  = SSUM(iim,zsin,1)/pi
504c           pcvgv(klon,l) = SSUM(iim,zsinbis,1)/pi
505
506        ENDDO
507c$OMP END DO NOWAIT       
508      endif
509
510
511      IF (is_sequential) THEN
512c
513cIM calcul PV a teta=350, 380, 405K
514        CALL PVtheta(ngridmx,llm,pucov,pvcov,pteta,
515     $           ztfi,zplay,zplev,
516     $           ntetaSTD,rtetaSTD,PVteta)
517c
518      ENDIF
519
520c On change de grille, dynamique vers physiq, pour le flux de masse verticale
521      CALL gr_dyn_fi_p(llm,iip1,jjp1,klon,flxw,flxwfi)
522
523c-----------------------------------------------------------------------
524c   Appel de la physique:
525c   ---------------------
526
527cc$OMP  PARALLEL DEFAULT(NONE)
528cc$OMP+ PRIVATE(i,l,offset,iq)
529cc$OMP+ SHARED(klon_omp_nb,nq,klon_omp_begin,
530cc$OMP+        debut,lafin,rdayvrai,heure,dtphys,zplev,zplay,
531cc$OMP+        zphi,zphis,presnivs,clesphy0,zufi,zvfi,ztfi,
532cc$OMP+        zqfi,pvervel,zdufi,zdvfi,zdtfi,zdqfi,zdpsrf)
533
534c PRIVATE(zplev_omp,zplay_omp,zphi_omp,zphis_omp,
535c c$OMP+                 presnivs_omp,zufi_omp,zvfi_omp,ztfi_omp,
536c c$OMP+                 zqfi_omp,pvervel_omp,zdufi_omp,zdvfi_omp,
537c c$OMP+                 zdtfi_omp,zdqfi_omp,zdpsrf_omp)
538
539c$OMP BARRIER
540      if (first_omp) then
541        klon=klon_omp
542
543        allocate(zplev_omp(klon,llm+1))
544        allocate(zplay_omp(klon,llm))
545        allocate(zphi_omp(klon,llm))
546        allocate(zphis_omp(klon))
547        allocate(presnivs_omp(llm))
548        allocate(zufi_omp(klon,llm))
549        allocate(zvfi_omp(klon,llm))
550        allocate(ztfi_omp(klon,llm))
551        allocate(zqfi_omp(klon,llm,nq))
552c        allocate(pvervel_omp(klon,llm))
553        allocate(zdufi_omp(klon,llm))
554        allocate(zdvfi_omp(klon,llm))
555        allocate(zdtfi_omp(klon,llm))
556        allocate(zdqfi_omp(klon,llm,nq))
557        allocate(zdpsrf_omp(klon))
558        allocate(flxwfi_omp(klon,llm))
559        first_omp=.false.
560      endif
561       
562           
563      klon=klon_omp
564      offset=klon_omp_begin-1
565     
566      do l=1,llm+1
567        do i=1,klon
568          zplev_omp(i,l)=zplev(offset+i,l)
569        enddo
570      enddo
571         
572       do l=1,llm
573        do i=1,klon 
574          zplay_omp(i,l)=zplay(offset+i,l)
575        enddo
576      enddo
577       
578      do l=1,llm
579        do i=1,klon
580          zphi_omp(i,l)=zphi(offset+i,l)
581        enddo
582      enddo
583       
584      do i=1,klon
585        zphis_omp(i)=zphis(offset+i)
586      enddo
587     
588       
589      do l=1,llm
590        presnivs_omp(l)=presnivs(l)
591      enddo
592       
593      do l=1,llm
594        do i=1,klon
595          zufi_omp(i,l)=zufi(offset+i,l)
596        enddo
597      enddo
598       
599      do l=1,llm
600        do i=1,klon
601          zvfi_omp(i,l)=zvfi(offset+i,l)
602        enddo
603      enddo
604       
605      do l=1,llm
606        do i=1,klon
607          ztfi_omp(i,l)=ztfi(offset+i,l)
608        enddo
609      enddo
610       
611      do iq=1,nq
612        do l=1,llm
613          do i=1,klon
614            zqfi_omp(i,l,iq)=zqfi(offset+i,l,iq)
615          enddo
616        enddo
617      enddo
618       
619c      do l=1,llm
620c        do i=1,klon
621c         pvervel_omp(i,l)=pvervel(offset+i,l)
622c       enddo
623c      enddo
624       
625      do l=1,llm
626        do i=1,klon
627          zdufi_omp(i,l)=zdufi(offset+i,l)
628        enddo
629      enddo
630       
631      do l=1,llm
632        do i=1,klon
633          zdvfi_omp(i,l)=zdvfi(offset+i,l)
634        enddo
635      enddo
636       
637      do l=1,llm
638        do i=1,klon
639          zdtfi_omp(i,l)=zdtfi(offset+i,l)
640        enddo
641      enddo
642       
643      do iq=1,nq
644        do l=1,llm
645          do i=1,klon
646            zdqfi_omp(i,l,iq)=zdqfi(offset+i,l,iq)
647          enddo
648        enddo
649      enddo
650       
651      do i=1,klon
652        zdpsrf_omp(i)=zdpsrf(offset+i)
653      enddo
654
655      do l=1,llm
656        do i=1,klon
657          flxwfi_omp(i,l)=flxwfi(offset+i,l)
658        enddo
659      enddo
660     
661c$OMP BARRIER
662cym      call WriteField_phy_p('zdtfi_omp',zdtfi_omp(:,:),llm)
663     
664      CALL physiq (klon,
665     .             llm,
666     .             nq,
667     .             debut,
668     .             lafin,
669     .             rdayvrai,
670     .             heure,
671     .             dtphys,
672     .             zplev_omp,
673     .             zplay_omp,
674     .             zphi_omp,
675     .             zphis_omp,
676     .             presnivs_omp,
677     .             clesphy0,
678     .             zufi_omp,
679     .             zvfi_omp,
680     .             ztfi_omp,
681     .             zqfi_omp,
682c     .             pvervel_omp,
683c#ifdef INCA
684     .             flxwfi_omp,
685c#endif
686     .             zdufi_omp,
687     .             zdvfi_omp,
688     .             zdtfi_omp,
689     .             zdqfi_omp,
690     .             zdpsrf_omp,
691cIM diagnostique PVteta, Amip2         
692     .             pducov,
693     .             PVteta)
694
695cym      call WriteField_phy_p('zdtfi_omp',zdtfi_omp(:,:),llm)
696
697c$OMP BARRIER
698
699      do l=1,llm+1
700        do i=1,klon
701          zplev(offset+i,l)=zplev_omp(i,l)
702        enddo
703      enddo
704         
705       do l=1,llm
706        do i=1,klon 
707          zplay(offset+i,l)=zplay_omp(i,l)
708        enddo
709      enddo
710       
711      do l=1,llm
712        do i=1,klon
713          zphi(offset+i,l)=zphi_omp(i,l)
714        enddo
715      enddo
716       
717
718      do i=1,klon
719        zphis(offset+i)=zphis_omp(i)
720      enddo
721     
722       
723      do l=1,llm
724        presnivs(l)=presnivs_omp(l)
725      enddo
726       
727      do l=1,llm
728        do i=1,klon
729          zufi(offset+i,l)=zufi_omp(i,l)
730        enddo
731      enddo
732       
733      do l=1,llm
734        do i=1,klon
735          zvfi(offset+i,l)=zvfi_omp(i,l)
736        enddo
737      enddo
738       
739      do l=1,llm
740        do i=1,klon
741          ztfi(offset+i,l)=ztfi_omp(i,l)
742        enddo
743      enddo
744       
745      do iq=1,nq
746        do l=1,llm
747          do i=1,klon
748            zqfi(offset+i,l,iq)=zqfi_omp(i,l,iq)
749          enddo
750        enddo
751      enddo
752       
753c      do l=1,llm
754c        do i=1,klon
755c         pvervel(offset+i,l)=pvervel_omp(i,l)
756c       enddo
757c      enddo
758       
759      do l=1,llm
760        do i=1,klon
761          zdufi(offset+i,l)=zdufi_omp(i,l)
762        enddo
763      enddo
764       
765      do l=1,llm
766        do i=1,klon
767          zdvfi(offset+i,l)=zdvfi_omp(i,l)
768        enddo
769      enddo
770       
771      do l=1,llm
772        do i=1,klon
773          zdtfi(offset+i,l)=zdtfi_omp(i,l)
774        enddo
775      enddo
776       
777      do iq=1,nq
778        do l=1,llm
779          do i=1,klon
780            zdqfi(offset+i,l,iq)=zdqfi_omp(i,l,iq)
781          enddo
782        enddo
783      enddo
784       
785      do i=1,klon
786        zdpsrf(offset+i)=zdpsrf_omp(i)
787      enddo
788     
789
790cc$OMP END PARALLEL
791      klon=klon_mpi
792500   CONTINUE
793c$OMP BARRIER
794
795c$OMP MASTER
796cym      call WriteField_phy('zdtfi',zdtfi(:,:),llm)
797      call stop_timer(timer_physic)
798c$OMP END MASTER
799
800      IF (using_mpi) THEN
801           
802      if (MPI_rank>0) then
803
804c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)       
805       DO l=1,llm     
806        du_send(1:iim,l)=zdufi(1:iim,l)
807        dv_send(1:iim,l)=zdvfi(1:iim,l)
808       ENDDO
809c$OMP END DO NOWAIT       
810
811c$OMP BARRIER
812#ifdef CPP_MPI
813c$OMP MASTER
814!$OMP CRITICAL (MPI)
815        call MPI_ISSEND(du_send,iim*llm,MPI_REAL8,MPI_Rank-1,401,
816     &                   COMM_LMDZ,Req(1),ierr)
817        call MPI_ISSEND(dv_send,iim*llm,MPI_REAL8,MPI_Rank-1,402,
818     &                  COMM_LMDZ,Req(2),ierr)
819!$OMP END CRITICAL (MPI)
820c$OMP END MASTER
821#endif
822c$OMP BARRIER
823     
824      endif
825   
826      if (MPI_rank<MPI_Size-1) then
827c$OMP BARRIER
828#ifdef CPP_MPI
829c$OMP MASTER     
830!$OMP CRITICAL (MPI)
831        call MPI_IRECV(du_recv,iim*llm,MPI_REAL8,MPI_Rank+1,401,
832     &                 COMM_LMDZ,Req(3),ierr)
833        call MPI_IRECV(dv_recv,iim*llm,MPI_REAL8,MPI_Rank+1,402,
834     &                 COMM_LMDZ,Req(4),ierr)
835!$OMP END CRITICAL (MPI)
836c$OMP END MASTER
837#endif
838      endif
839
840c$OMP BARRIER
841
842
843#ifdef CPP_MPI
844c$OMP MASTER   
845!$OMP CRITICAL (MPI)
846      if (MPI_rank>0 .and. MPI_rank< MPI_Size-1) then
847        call MPI_WAITALL(4,Req(1),Status,ierr)
848      else if (MPI_rank>0) then
849        call MPI_WAITALL(2,Req(1),Status,ierr)
850      else if (MPI_rank <MPI_Size-1) then
851        call MPI_WAITALL(2,Req(3),Status,ierr)
852      endif
853!$OMP END CRITICAL (MPI)
854c$OMP END MASTER
855#endif
856
857c$OMP BARRIER     
858
859      ENDIF ! using_mpi
860     
861     
862c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
863      DO l=1,llm
864           
865        zdufi2(1:klon,l)=zdufi(1:klon,l)
866        zdufi2(klon+1:klon+iim,l)=du_recv(1:iim,l)
867           
868        zdvfi2(1:klon,l)=zdvfi(1:klon,l)
869        zdvfi2(klon+1:klon+iim,l)=dv_recv(1:iim,l)
870 
871         pdhfi(:,jj_begin,l)=0
872         pdqfi(:,jj_begin,l,:)=0
873         pdufi(:,jj_begin,l)=0
874         pdvfi(:,jj_begin,l)=0
875         
876         if (.not. is_south_pole) then
877           pdhfi(:,jj_end,l)=0
878           pdqfi(:,jj_end,l,:)=0
879           pdufi(:,jj_end,l)=0
880           pdvfi(:,jj_end,l)=0
881         endif
882     
883       ENDDO
884c$OMP END DO NOWAIT
885
886c$OMP MASTER
887       pdpsfi(:,jj_begin)=0   
888       if (.not. is_south_pole) then
889         pdpsfi(:,jj_end)=0
890       endif
891c$OMP END MASTER
892c-----------------------------------------------------------------------
893c   transformation des tendances physiques en tendances dynamiques:
894c   ---------------------------------------------------------------
895
896c  tendance sur la pression :
897c  -----------------------------------
898c$OMP MASTER
899      CALL gr_fi_dyn_p(1,klon,iip1,jjp1,zdpsrf,pdpsfi)
900c$OMP END MASTER
901c
902c   62. enthalpie potentielle
903c   ---------------------
904     
905      kstart=1
906      kend=klon
907
908      if (is_north_pole) kstart=2
909      if (is_south_pole)  kend=klon-1
910
911c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
912      DO l=1,llm
913
914!!cdir NODEP
915        do ig0=kstart,kend
916          i=index_i(ig0)
917          j=index_j(ig0)
918          pdhfi(i,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l)
919          if (i==1) pdhfi(iip1,j,l) =  cpp * zdtfi(ig0,l) / ppk(i,j,l)
920         enddo         
921
922        if (is_north_pole) then
923            DO i=1,iip1
924              pdhfi(i,1,l)    = cpp *  zdtfi(1,l)      / ppk(i, 1  ,l)
925            enddo
926        endif
927       
928        if (is_south_pole) then
929            DO i=1,iip1
930              pdhfi(i,jjp1,l) = cpp *  zdtfi(klon,l)/ ppk(i,jjp1,l)
931            ENDDO
932        endif
933      ENDDO
934c$OMP END DO NOWAIT
935     
936c   62. humidite specifique
937c   ---------------------
938
939      DO iq=1,nqmx
940c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
941         DO l=1,llm
942!!cdir NODEP
943           do ig0=kstart,kend
944             i=index_i(ig0)
945             j=index_j(ig0)
946             pdqfi(i,j,l,iq) = zdqfi(ig0,l,iq)
947             if (i==1) pdqfi(iip1,j,l,iq) = zdqfi(ig0,l,iq)
948           enddo
949           
950           if (is_north_pole) then
951             do i=1,iip1
952               pdqfi(i,1,l,iq)    = zdqfi(1,l,iq)             
953             enddo
954           endif
955           
956           if (is_south_pole) then
957             do i=1,iip1
958               pdqfi(i,jjp1,l,iq) = zdqfi(klon,l,iq)
959             enddo
960           endif
961           
962         ENDDO
963c$OMP END DO NOWAIT
964      ENDDO
965
966c   63. traceurs
967c   ------------
968C     initialisation des tendances
969
970c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
971      DO l=1,llm
972        pdqfi(:,:,l,:)=0.
973      ENDDO
974c$OMP END DO NOWAIT     
975
976C
977
978      DO iq=1,nq
979         iiq=niadv(iq)
980c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
981         DO l=1,llm
982
983!!cdir NODEP           
984             DO ig0=kstart,kend
985              i=index_i(ig0)
986              j=index_j(ig0)
987              pdqfi(i,j,l,iiq) = zdqfi(ig0,l,iq)
988              if (i==1) pdqfi(iip1,j,l,iiq) = zdqfi(ig0,l,iq)
989            ENDDO
990           
991            IF (is_north_pole) then
992              DO i=1,iip1
993                pdqfi(i,1,l,iiq)    = zdqfi(1,l,iq)
994              ENDDO
995            ENDIF
996           
997            IF (is_south_pole) then
998              DO i=1,iip1
999                pdqfi(i,jjp1,l,iiq) = zdqfi(klon,l,iq)
1000              ENDDO
1001            ENDIF
1002           
1003         ENDDO
1004c$OMP END DO NOWAIT     
1005      ENDDO
1006     
1007c   65. champ u:
1008c   ------------
1009c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1010      DO l=1,llm
1011!!cdir NODEP
1012         do ig0=kstart,kend
1013           i=index_i(ig0)
1014           j=index_j(ig0)
1015           
1016           if (i/=iim) then
1017             pdufi(i,j,l)=0.5*(zdufi2(ig0,l)+zdufi2(ig0+1,l))*cu(i,j)
1018           endif
1019           
1020           if (i==1) then
1021              pdufi(iim,j,l)=0.5*(  zdufi2(ig0,l)
1022     $                            + zdufi2(ig0+iim-1,l))*cu(iim,j)
1023             pdufi(iip1,j,l)=0.5*(zdufi2(ig0,l)+zdufi2(ig0+1,l))*cu(i,j)
1024           endif
1025         
1026         enddo
1027         
1028         if (is_north_pole) then
1029           DO i=1,iip1
1030            pdufi(i,1,l)    = 0.
1031           ENDDO
1032         endif
1033         
1034         if (is_south_pole) then
1035           DO i=1,iip1
1036            pdufi(i,jjp1,l) = 0.
1037           ENDDO
1038         endif
1039         
1040      ENDDO
1041c$OMP END DO NOWAIT
1042
1043c   67. champ v:
1044c   ------------
1045
1046      kstart=1
1047      kend=klon
1048
1049      if (is_north_pole) kstart=2
1050      if (is_south_pole)  kend=klon-1-iim
1051     
1052c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
1053      DO l=1,llm
1054!!cdir NODEP
1055        do ig0=kstart,kend
1056           i=index_i(ig0)
1057           j=index_j(ig0)
1058           pdvfi(i,j,l)=0.5*(zdvfi2(ig0,l)+zdvfi2(ig0+iim,l))*cv(i,j)
1059           if (i==1) pdvfi(iip1,j,l) = 0.5*(zdvfi2(ig0,l)+
1060     $                                      zdvfi2(ig0+iim,l))
1061     $                                    *cv(i,j)
1062        enddo
1063         
1064      ENDDO
1065c$OMP END DO NOWAIT
1066
1067
1068c   68. champ v pres des poles:
1069c   ---------------------------
1070c      v = U * cos(long) + V * SIN(long)
1071
1072      if (is_north_pole) then
1073
1074c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
1075        DO l=1,llm
1076
1077          DO i=1,iim
1078            pdvfi(i,1,l)=
1079     $      zdufi(1,l)*COS(rlonv(i))+zdvfi(1,l)*SIN(rlonv(i))
1080       
1081            pdvfi(i,1,l)=
1082     $      0.5*(pdvfi(i,1,l)+zdvfi(i+1,l))*cv(i,1)
1083          ENDDO
1084
1085          pdvfi(iip1,1,l)  = pdvfi(1,1,l)
1086
1087        ENDDO
1088c$OMP END DO NOWAIT
1089
1090      endif   
1091     
1092      if (is_south_pole) then
1093
1094c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
1095         DO l=1,llm
1096 
1097           DO i=1,iim
1098              pdvfi(i,jjm,l)=zdufi(klon,l)*COS(rlonv(i))
1099     $        +zdvfi(klon,l)*SIN(rlonv(i))
1100
1101              pdvfi(i,jjm,l)=
1102     $        0.5*(pdvfi(i,jjm,l)+zdvfi(klon-iip1+i,l))*cv(i,jjm)
1103           ENDDO
1104
1105           pdvfi(iip1,jjm,l)= pdvfi(1,jjm,l)
1106
1107        ENDDO
1108c$OMP END DO NOWAIT
1109     
1110      endif
1111c-----------------------------------------------------------------------
1112
1113700   CONTINUE
1114 
1115      firstcal = .FALSE.
1116
1117      RETURN
1118      END
Note: See TracBrowser for help on using the repository browser.