source: LMDZ6/branches/Amaury_dev/libf/dynphy_lonlat/lmdz_calfis_loc.F90 @ 5137

Last change on this file since 5137 was 5136, checked in by abarral, 8 weeks ago

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