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