source: LMDZ6/trunk/libf/dynphy_lonlat/calfis_loc.F @ 4127

Last change on this file since 4127 was 4056, checked in by dcugnet, 3 years ago

Most of the changes are intended to help to eventually remove the constraints about the tracers assumptions, in particular water tracers.

  • Remove index tables itr_indice and niadv, replaced by tracers(:)%isAdvected and tracers(:)%isH2OFamily. Most of the loops are now from 1 to nqtot:
    • DO iq=nqo+1,nqtot loops are replaced with: DO iq=1,nqtot

IF(tracers(iq)%isH2Ofamily) CYCLE

  • DO it=1,nbtr; iq=niadv(it+nqo)

and DO it=1,nqtottr; iq=itr_indice(it) loops are replaced with:

it = 0
DO iq = 1, nqtot

IF(.NOT.tracers(iq)%isAdvected .OR. tracers(iq)%isH2Ofamily) CYCLE
it = it+1

  • Move some StratAer? related code from infotrac to infotrac_phy
  • Remove "nqperes" variable:

DO iq=1,nqpere loops are replaced with:
DO iq=1,nqtot

IF(tracers(iq)%parent/='air') CYCLE

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