source: LMDZ5/trunk/libf/dynlmdz_phylmd/calfis_p.F @ 2239

Last change on this file since 2239 was 2239, checked in by Ehouarn Millour, 10 years ago

Reorganizing physics/dynamics interface:

  • what is related to dynamics-physics interface is now in a seperate directory: dynlmdz_phy* for physics in phy*
  • 1d model and related dependencies (including a couple from "dynamics", set up as symbolic links) is now in subdirectory "dyn1d" of phy*.
  • "bibio" directory is now "misc" and should only contain autonomous utilities.
  • "cosp" is now a subdirectory of phylmd.

EM

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