trunk/dyn3d/advtrac.f

SUBROUTINE advtrac(pbaru,pbarv,p,masse,q,iapptrac,teta,pk)

  ! From dyn3d/advtrac.F,v 1.4 2005/04/13 08:58:34
  !     Auteur :  F. Hourdin

  !     Modif. P. Le Van     (20/12/97)
  !            F. Codron     (10/99)
  !            D. Le Croller (07/2001)
  !            M.A Filiberti (04/2002)

  USE dimens_m
  USE paramet_m
  USE comconst
  USE comvert
  USE conf_gcm_m
  USE logic
  USE comgeom
  USE temps
  USE ener
  USE advtrac_m
  USE comdissip
  IMPLICIT NONE


  !-------------------------------------------------------------------
  !     Arguments
  !-------------------------------------------------------------------
  !     Ajout PPM
  !--------------------------------------------------------
  REAL massebx(ip1jmp1,llm), masseby(ip1jm,llm)
  !--------------------------------------------------------
  INTEGER iapptrac
  REAL pbaru(ip1jmp1,llm), pbarv(ip1jm,llm)
  REAL q(ip1jmp1,llm,nqmx), masse(ip1jmp1,llm)
  REAL p(ip1jmp1,llmp1), teta(ip1jmp1,llm)
  REAL pk(ip1jmp1,llm)

  !-------------------------------------------------------------
  !     Variables locales
  !-------------------------------------------------------------

  REAL pbaruc(ip1jmp1,llm), pbarvc(ip1jm,llm)
  REAL massem(ip1jmp1,llm), zdp(ip1jmp1)
  REAL pbarug(ip1jmp1,llm), pbarvg(ip1jm,llm), wg(ip1jmp1,llm)
  REAL cpuadv(nqmx)
  COMMON /cpuadv/cpuadv

  INTEGER iadvtr
  INTEGER ij, l, iq
  REAL zdpmin, zdpmax
  EXTERNAL minmax
  SAVE iadvtr, massem, pbaruc, pbarvc
  DATA iadvtr/0/
  !----------------------------------------------------------
  !     Rajouts pour PPM
  !----------------------------------------------------------
  INTEGER indice, n
  ! Pas de temps adaptatif pour que CFL<1                
  REAL dtbon
  REAL cflmaxz, aaa, & ! CFL maximum                                
       bbb
  REAL psppm(iim,jjp1) ! pression  au sol                           
  REAL unatppm(iim,jjp1,llm), vnatppm(iim,jjp1,llm)
  REAL qppm(iim*jjp1,llm,nqmx)
  REAL fluxwppm(iim,jjp1,llm)
  REAL apppm(llmp1), bpppm(llmp1)
  LOGICAL dum, fill
  DATA fill/ .TRUE./
  DATA dum/ .TRUE./


  IF (iadvtr==0) THEN
     CALL initial0(ijp1llm,pbaruc)
     CALL initial0(ijmllm,pbarvc)
  END IF

  !   accumulation des flux de masse horizontaux
  DO l = 1, llm
     DO ij = 1, ip1jmp1
        pbaruc(ij,l) = pbaruc(ij,l) + pbaru(ij,l)
     END DO
     DO ij = 1, ip1jm
        pbarvc(ij,l) = pbarvc(ij,l) + pbarv(ij,l)
     END DO
  END DO

  !   selection de la masse instantannee des mailles avant le transport.
  IF (iadvtr==0) THEN

     CALL scopy(ip1jmp1*llm,masse,1,massem,1)
     !cc         CALL filtreg ( massem ,jjp1, llm,-2, 2, .TRUE., 1 )

  END IF

  iadvtr = iadvtr + 1
  iapptrac = iadvtr


  !   Test pour savoir si on advecte a ce pas de temps
  IF (iadvtr==iapp_tracvl) THEN

     !c   ..  Modif P.Le Van  ( 20/12/97 )  ....
     !c

     !   traitement des flux de masse avant advection.
     !     1. calcul de w
     !     2. groupement des mailles pres du pole.

     CALL groupe(massem,pbaruc,pbarvc,pbarug,pbarvg,wg)


     !  test sur l'eventuelle creation de valeurs negatives de la masse
     DO l = 1, llm - 1
        DO ij = iip2 + 1, ip1jm
           zdp(ij) = pbarug(ij-1,l) - pbarug(ij,l) - pbarvg(ij-iip1,l) + &
                pbarvg(ij,l) + wg(ij,l+1) - wg(ij,l)
        END DO
        CALL scopy(jjm-1,zdp(iip1+iip1),iip1,zdp(iip2),iip1)
        DO ij = iip2, ip1jm
           zdp(ij) = zdp(ij)*dtvr/massem(ij,l)
        END DO


        CALL minmax(ip1jm-iip1,zdp(iip2),zdpmin,zdpmax)

        IF (max(abs(zdpmin),abs(zdpmax))>0.5) THEN
           PRINT *, 'WARNING DP/P l=', l, '  MIN:', zdpmin, '   MAX:', zdpmax
        END IF

     END DO

     !-------------------------------------------------------------------
     !   Advection proprement dite (Modification Le Croller (07/2001)
     !-------------------------------------------------------------------

     !----------------------------------------------------
     !        Calcul des moyennes basées sur la masse
     !----------------------------------------------------
     CALL massbar(massem,massebx,masseby)

     !-----------------------------------------------------------
     !     Appel des sous programmes d'advection
     !-----------------------------------------------------------
     DO iq = 1, nqmx
        !        call clock(t_initial)
        IF (iadv(iq)==0) CYCLE
        !   ----------------------------------------------------------------
        !   Schema de Van Leer I MUSCL
        !   ----------------------------------------------------------------
        IF (iadv(iq)==10) THEN
           CALL vlsplt(q(1,1,iq),2.,massem,wg,pbarug,pbarvg,dtvr)


           !   ----------------------------------------------------------------
           !   Schema "pseudo amont" + test sur humidite specifique
           !    pour la vapeur d'eau. F. Codron
           !   ----------------------------------------------------------------
        ELSE IF (iadv(iq)==14) THEN

           CALL vlspltqs(q(1,1,1),2.,massem,wg,pbarug,pbarvg,dtvr,p,pk,teta)
           !   ----------------------------------------------------------------
           !   Schema de Frederic Hourdin
           !   ----------------------------------------------------------------
        ELSE IF (iadv(iq)==12) THEN
           !            Pas de temps adaptatif
           CALL adaptdt(iadv(iq),dtbon,n,pbarug,massem)
           IF (n>1) THEN
              WRITE (*,*) 'WARNING horizontal dt=', dtbon, 'dtvr=', dtvr, &
                   'n=', n
           END IF
           DO indice = 1, n
              CALL advn(q(1,1,iq),massem,wg,pbarug,pbarvg,dtbon,1)
           END DO
        ELSE IF (iadv(iq)==13) THEN
           !            Pas de temps adaptatif
           CALL adaptdt(iadv(iq),dtbon,n,pbarug,massem)
           IF (n>1) THEN
              WRITE (*,*) 'WARNING horizontal dt=', dtbon, 'dtvr=', dtvr, &
                   'n=', n
           END IF
           DO indice = 1, n
              CALL advn(q(1,1,iq),massem,wg,pbarug,pbarvg,dtbon,2)
           END DO
           !   ----------------------------------------------------------------
           !   Schema de pente SLOPES
           !   ----------------------------------------------------------------
        ELSE IF (iadv(iq)==20) THEN
           CALL pentes_ini(q(1,1,iq),wg,massem,pbarug,pbarvg,0)

           !   ----------------------------------------------------------------
           !   Schema de Prather
           !   ----------------------------------------------------------------
        ELSE IF (iadv(iq)==30) THEN
           !            Pas de temps adaptatif
           CALL adaptdt(iadv(iq),dtbon,n,pbarug,massem)
           IF (n>1) THEN
              WRITE (*,*) 'WARNING horizontal dt=', dtbon, 'dtvr=', dtvr, &
                   'n=', n
           END IF
           CALL prather(q(1,1,iq),wg,massem,pbarug,pbarvg,n,dtbon)
           !   ----------------------------------------------------------------
           !   Schemas PPM Lin et Rood
           !   ----------------------------------------------------------------
        ELSE IF (iadv(iq)==11 .OR. (iadv(iq)>=16 .AND. iadv(iq)<=18)) THEN

           !        Test sur le flux horizontal
           !        Pas de temps adaptatif
           CALL adaptdt(iadv(iq),dtbon,n,pbarug,massem)
           IF (n>1) THEN
              WRITE (*,*) 'WARNING horizontal dt=', dtbon, 'dtvr=', dtvr, &
                   'n=', n
           END IF
           !        Test sur le flux vertical
           cflmaxz = 0.
           DO l = 2, llm
              DO ij = iip2, ip1jm
                 aaa = wg(ij,l)*dtvr/massem(ij,l)
                 cflmaxz = max(cflmaxz,aaa)
                 bbb = -wg(ij,l)*dtvr/massem(ij,l-1)
                 cflmaxz = max(cflmaxz,bbb)
              END DO
           END DO
           IF (cflmaxz>=1) THEN
              WRITE (*,*) 'WARNING vertical', 'CFLmaxz=', cflmaxz
           END IF

           !-----------------------------------------------------------
           !        Ss-prg interface LMDZ.4->PPM3d
           !-----------------------------------------------------------

           CALL interpre(q(1,1,iq),qppm(1,1,iq),wg,fluxwppm,massem,apppm, &
                bpppm,massebx,masseby,pbarug,pbarvg,unatppm,vnatppm,psppm)

           DO indice = 1, n
              !----------------------------------------------------------------
              !                         VL (version PPM) horiz. et PPM vert.
              !---------------------------------------------------------------
              IF (iadv(iq)==11) THEN
                 !                  Ss-prg PPM3d de Lin
                 CALL ppm3d(1,qppm(1,1,iq),psppm,psppm,unatppm,vnatppm, &
                      fluxwppm,dtbon,2,2,2,1,iim,jjp1,2,llm,apppm,bpppm,0.01, &
                      6400000,fill,dum,220.)

                 !-----------------------------------------------------------
                 !                           Monotonic PPM
                 !-------------------------------------------------------
              ELSE IF (iadv(iq)==16) THEN
                 !                  Ss-prg PPM3d de Lin
                 CALL ppm3d(1,qppm(1,1,iq),psppm,psppm,unatppm,vnatppm, &
                      fluxwppm,dtbon,3,3,3,1,iim,jjp1,2,llm,apppm,bpppm,0.01, &
                      6400000,fill,dum,220.)
                 !                           Semi Monotonic PPM
              ELSE IF (iadv(iq)==17) THEN
                 !                  Ss-prg PPM3d de Lin
                 CALL ppm3d(1,qppm(1,1,iq),psppm,psppm,unatppm,vnatppm, &
                      fluxwppm,dtbon,4,4,4,1,iim,jjp1,2,llm,apppm,bpppm,0.01, &
                      6400000,fill,dum,220.)
                 !                         Positive Definite PPM
              ELSE IF (iadv(iq)==18) THEN
                 !                  Ss-prg PPM3d de Lin
                 CALL ppm3d(1,qppm(1,1,iq),psppm,psppm,unatppm,vnatppm, &
                      fluxwppm,dtbon,5,5,5,1,iim,jjp1,2,llm,apppm,bpppm,0.01, &
                      6400000,fill,dum,220.)
              END IF
           END DO
           !-----------------------------------------------------------------
           !               Ss-prg interface PPM3d-LMDZ.4
           !-----------------------------------------------------------------
           CALL interpost(q(1,1,iq),qppm(1,1,iq))
        END IF
        !----------------------------------------------------------------------

        !-----------------------------------------------------------------
        ! On impose une seule valeur du traceur au pôle Sud j=jjm+1=jjp1
        ! et Nord j=1
        !-----------------------------------------------------------------

        !                  call traceurpole(q(1,1,iq),massem)

        ! calcul du temps cpu pour un schema donne

        !                  call clock(t_final)
        !ym                  tps_cpu=t_final-t_initial
        !ym                  cpuadv(iq)=cpuadv(iq)+tps_cpu

     END DO


     !------------------------------------------------------------------
     !   on reinitialise a zero les flux de masse cumules
     !---------------------------------------------------
     iadvtr = 0


  END IF
  ! if iadvtr.EQ.iapp_tracvl                                 
  RETURN
END SUBROUTINE advtrac
1	guez	3	SUBROUTINE advtrac(pbaru,pbarv,p,masse,q,iapptrac,teta,pk)
2
3			! From dyn3d/advtrac.F,v 1.4 2005/04/13 08:58:34
4			! Auteur : F. Hourdin
5
6			! Modif. P. Le Van (20/12/97)
7			! F. Codron (10/99)
8			! D. Le Croller (07/2001)
9			! M.A Filiberti (04/2002)
10
11			USE dimens_m
12			USE paramet_m
13			USE comconst
14			USE comvert
15			USE conf_gcm_m
16			USE logic
17			USE comgeom
18			USE temps
19			USE ener
20			USE advtrac_m
21			USE comdissip
22			IMPLICIT NONE
23
24
25			!-------------------------------------------------------------------
26			! Arguments
27			!-------------------------------------------------------------------
28			! Ajout PPM
29			!--------------------------------------------------------
30			REAL massebx(ip1jmp1,llm), masseby(ip1jm,llm)
31			!--------------------------------------------------------
32			INTEGER iapptrac
33			REAL pbaru(ip1jmp1,llm), pbarv(ip1jm,llm)
34			REAL q(ip1jmp1,llm,nqmx), masse(ip1jmp1,llm)
35			REAL p(ip1jmp1,llmp1), teta(ip1jmp1,llm)
36			REAL pk(ip1jmp1,llm)
37
38			!-------------------------------------------------------------
39			! Variables locales
40			!-------------------------------------------------------------
41
42			REAL pbaruc(ip1jmp1,llm), pbarvc(ip1jm,llm)
43			REAL massem(ip1jmp1,llm), zdp(ip1jmp1)
44			REAL pbarug(ip1jmp1,llm), pbarvg(ip1jm,llm), wg(ip1jmp1,llm)
45			REAL cpuadv(nqmx)
46			COMMON /cpuadv/cpuadv
47
48			INTEGER iadvtr
49			INTEGER ij, l, iq
50			REAL zdpmin, zdpmax
51			EXTERNAL minmax
52			SAVE iadvtr, massem, pbaruc, pbarvc
53			DATA iadvtr/0/
54			!----------------------------------------------------------
55			! Rajouts pour PPM
56			!----------------------------------------------------------
57			INTEGER indice, n
58			! Pas de temps adaptatif pour que CFL<1
59			REAL dtbon
60			REAL cflmaxz, aaa, & ! CFL maximum
61			bbb
62			REAL psppm(iim,jjp1) ! pression au sol
63			REAL unatppm(iim,jjp1,llm), vnatppm(iim,jjp1,llm)
64			REAL qppm(iim*jjp1,llm,nqmx)
65			REAL fluxwppm(iim,jjp1,llm)
66			REAL apppm(llmp1), bpppm(llmp1)
67			LOGICAL dum, fill
68			DATA fill/ .TRUE./
69			DATA dum/ .TRUE./
70
71
72			IF (iadvtr==0) THEN
73			CALL initial0(ijp1llm,pbaruc)
74			CALL initial0(ijmllm,pbarvc)
75			END IF
76
77			! accumulation des flux de masse horizontaux
78			DO l = 1, llm
79			DO ij = 1, ip1jmp1
80			pbaruc(ij,l) = pbaruc(ij,l) + pbaru(ij,l)
81			END DO
82			DO ij = 1, ip1jm
83			pbarvc(ij,l) = pbarvc(ij,l) + pbarv(ij,l)
84			END DO
85			END DO
86
87			! selection de la masse instantannee des mailles avant le transport.
88			IF (iadvtr==0) THEN
89
90			CALL scopy(ip1jmp1*llm,masse,1,massem,1)
91			!cc CALL filtreg ( massem ,jjp1, llm,-2, 2, .TRUE., 1 )
92
93			END IF
94
95			iadvtr = iadvtr + 1
96			iapptrac = iadvtr
97
98
99			! Test pour savoir si on advecte a ce pas de temps
100			IF (iadvtr==iapp_tracvl) THEN
101
102			!c .. Modif P.Le Van ( 20/12/97 ) ....
103			!c
104
105			! traitement des flux de masse avant advection.
106			! 1. calcul de w
107			! 2. groupement des mailles pres du pole.
108
109			CALL groupe(massem,pbaruc,pbarvc,pbarug,pbarvg,wg)
110
111
112			! test sur l'eventuelle creation de valeurs negatives de la masse
113			DO l = 1, llm - 1
114			DO ij = iip2 + 1, ip1jm
115			zdp(ij) = pbarug(ij-1,l) - pbarug(ij,l) - pbarvg(ij-iip1,l) + &
116			pbarvg(ij,l) + wg(ij,l+1) - wg(ij,l)
117			END DO
118			CALL scopy(jjm-1,zdp(iip1+iip1),iip1,zdp(iip2),iip1)
119			DO ij = iip2, ip1jm
120			zdp(ij) = zdp(ij)*dtvr/massem(ij,l)
121			END DO
122
123
124			CALL minmax(ip1jm-iip1,zdp(iip2),zdpmin,zdpmax)
125
126			IF (max(abs(zdpmin),abs(zdpmax))>0.5) THEN
127			PRINT *, 'WARNING DP/P l=', l, ' MIN:', zdpmin, ' MAX:', zdpmax
128			END IF
129
130			END DO
131
132			!-------------------------------------------------------------------
133			! Advection proprement dite (Modification Le Croller (07/2001)
134			!-------------------------------------------------------------------
135
136			!----------------------------------------------------
137			! Calcul des moyennes basées sur la masse
138			!----------------------------------------------------
139			CALL massbar(massem,massebx,masseby)
140
141			!-----------------------------------------------------------
142			! Appel des sous programmes d'advection
143			!-----------------------------------------------------------
144			DO iq = 1, nqmx
145			! call clock(t_initial)
146			IF (iadv(iq)==0) CYCLE
147			! ----------------------------------------------------------------
148			! Schema de Van Leer I MUSCL
149			! ----------------------------------------------------------------
150			IF (iadv(iq)==10) THEN
151			CALL vlsplt(q(1,1,iq),2.,massem,wg,pbarug,pbarvg,dtvr)
152
153
154			! ----------------------------------------------------------------
155			! Schema "pseudo amont" + test sur humidite specifique
156			! pour la vapeur d'eau. F. Codron
157			! ----------------------------------------------------------------
158			ELSE IF (iadv(iq)==14) THEN
159
160			CALL vlspltqs(q(1,1,1),2.,massem,wg,pbarug,pbarvg,dtvr,p,pk,teta)
161			! ----------------------------------------------------------------
162			! Schema de Frederic Hourdin
163			! ----------------------------------------------------------------
164			ELSE IF (iadv(iq)==12) THEN
165			! Pas de temps adaptatif
166			CALL adaptdt(iadv(iq),dtbon,n,pbarug,massem)
167			IF (n>1) THEN
168			WRITE (,) 'WARNING horizontal dt=', dtbon, 'dtvr=', dtvr, &
169			'n=', n
170			END IF
171			DO indice = 1, n
172			CALL advn(q(1,1,iq),massem,wg,pbarug,pbarvg,dtbon,1)
173			END DO
174			ELSE IF (iadv(iq)==13) THEN
175			! Pas de temps adaptatif
176			CALL adaptdt(iadv(iq),dtbon,n,pbarug,massem)
177			IF (n>1) THEN
178			WRITE (,) 'WARNING horizontal dt=', dtbon, 'dtvr=', dtvr, &
179			'n=', n
180			END IF
181			DO indice = 1, n
182			CALL advn(q(1,1,iq),massem,wg,pbarug,pbarvg,dtbon,2)
183			END DO
184			! ----------------------------------------------------------------
185			! Schema de pente SLOPES
186			! ----------------------------------------------------------------
187			ELSE IF (iadv(iq)==20) THEN
188			CALL pentes_ini(q(1,1,iq),wg,massem,pbarug,pbarvg,0)
189
190			! ----------------------------------------------------------------
191			! Schema de Prather
192			! ----------------------------------------------------------------
193			ELSE IF (iadv(iq)==30) THEN
194			! Pas de temps adaptatif
195			CALL adaptdt(iadv(iq),dtbon,n,pbarug,massem)
196			IF (n>1) THEN
197			WRITE (,) 'WARNING horizontal dt=', dtbon, 'dtvr=', dtvr, &
198			'n=', n
199			END IF
200			CALL prather(q(1,1,iq),wg,massem,pbarug,pbarvg,n,dtbon)
201			! ----------------------------------------------------------------
202			! Schemas PPM Lin et Rood
203			! ----------------------------------------------------------------
204			ELSE IF (iadv(iq)==11 .OR. (iadv(iq)>=16 .AND. iadv(iq)<=18)) THEN
205
206			! Test sur le flux horizontal
207			! Pas de temps adaptatif
208			CALL adaptdt(iadv(iq),dtbon,n,pbarug,massem)
209			IF (n>1) THEN
210			WRITE (,) 'WARNING horizontal dt=', dtbon, 'dtvr=', dtvr, &
211			'n=', n
212			END IF
213			! Test sur le flux vertical
214			cflmaxz = 0.
215			DO l = 2, llm
216			DO ij = iip2, ip1jm
217			aaa = wg(ij,l)*dtvr/massem(ij,l)
218			cflmaxz = max(cflmaxz,aaa)
219			bbb = -wg(ij,l)*dtvr/massem(ij,l-1)
220			cflmaxz = max(cflmaxz,bbb)
221			END DO
222			END DO
223			IF (cflmaxz>=1) THEN
224			WRITE (,) 'WARNING vertical', 'CFLmaxz=', cflmaxz
225			END IF
226
227			!-----------------------------------------------------------
228			! Ss-prg interface LMDZ.4->PPM3d
229			!-----------------------------------------------------------
230
231			CALL interpre(q(1,1,iq),qppm(1,1,iq),wg,fluxwppm,massem,apppm, &
232			bpppm,massebx,masseby,pbarug,pbarvg,unatppm,vnatppm,psppm)
233
234			DO indice = 1, n
235			!----------------------------------------------------------------
236			! VL (version PPM) horiz. et PPM vert.
237			!---------------------------------------------------------------
238			IF (iadv(iq)==11) THEN
239			! Ss-prg PPM3d de Lin
240			CALL ppm3d(1,qppm(1,1,iq),psppm,psppm,unatppm,vnatppm, &
241			fluxwppm,dtbon,2,2,2,1,iim,jjp1,2,llm,apppm,bpppm,0.01, &
242			6400000,fill,dum,220.)
243
244			!-----------------------------------------------------------
245			! Monotonic PPM
246			!-------------------------------------------------------
247			ELSE IF (iadv(iq)==16) THEN
248			! Ss-prg PPM3d de Lin
249			CALL ppm3d(1,qppm(1,1,iq),psppm,psppm,unatppm,vnatppm, &
250			fluxwppm,dtbon,3,3,3,1,iim,jjp1,2,llm,apppm,bpppm,0.01, &
251			6400000,fill,dum,220.)
252			! Semi Monotonic PPM
253			ELSE IF (iadv(iq)==17) THEN
254			! Ss-prg PPM3d de Lin
255			CALL ppm3d(1,qppm(1,1,iq),psppm,psppm,unatppm,vnatppm, &
256			fluxwppm,dtbon,4,4,4,1,iim,jjp1,2,llm,apppm,bpppm,0.01, &
257			6400000,fill,dum,220.)
258			! Positive Definite PPM
259			ELSE IF (iadv(iq)==18) THEN
260			! Ss-prg PPM3d de Lin
261			CALL ppm3d(1,qppm(1,1,iq),psppm,psppm,unatppm,vnatppm, &
262			fluxwppm,dtbon,5,5,5,1,iim,jjp1,2,llm,apppm,bpppm,0.01, &
263			6400000,fill,dum,220.)
264			END IF
265			END DO
266			!-----------------------------------------------------------------
267			! Ss-prg interface PPM3d-LMDZ.4
268			!-----------------------------------------------------------------
269			CALL interpost(q(1,1,iq),qppm(1,1,iq))
270			END IF
271			!----------------------------------------------------------------------
272
273			!-----------------------------------------------------------------
274			! On impose une seule valeur du traceur au pôle Sud j=jjm+1=jjp1
275			! et Nord j=1
276			!-----------------------------------------------------------------
277
278			! call traceurpole(q(1,1,iq),massem)
279
280			! calcul du temps cpu pour un schema donne
281
282			! call clock(t_final)
283			!ym tps_cpu=t_final-t_initial
284			!ym cpuadv(iq)=cpuadv(iq)+tps_cpu
285
286			END DO
287
288
289			!------------------------------------------------------------------
290			! on reinitialise a zero les flux de masse cumules
291			!---------------------------------------------------
292			iadvtr = 0
293
294
295			END IF
296			! if iadvtr.EQ.iapp_tracvl
297			RETURN
298			END SUBROUTINE advtrac