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

Last change on this file since 5174 was 5160, checked in by abarral, 3 months ago

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