New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
trcadv_tvd.F90 in trunk/NEMO/TOP_SRC/TRP – NEMO

source: trunk/NEMO/TOP_SRC/TRP/trcadv_tvd.F90 @ 941

Last change on this file since 941 was 941, checked in by cetlod, 16 years ago

phasing the passive tracer transport module to the new version of NEMO, see ticket 143
  • Property svn:executable set to *
File size: 20.9 KB
RevLine 
[941]1MODULE trcadv_tvd
2   !!==============================================================================
3   !!                       ***  MODULE  trcadv_tvd  ***
4   !! Ocean passive tracers:  horizontal & vertical advective trend
5   !!==============================================================================
6#if defined key_top
7   !!----------------------------------------------------------------------
8   !!   'key_top'                                                TOP models
9   !!----------------------------------------------------------------------
10   !!   trc_adv_tvd  : update the passive tracer trend with the horizontal
11   !!                  and vertical advection trends using a TVD scheme
12   !!   nonosc       : compute monotonic tracer fluxes by a nonoscillatory
13   !!                  algorithm
14   !!----------------------------------------------------------------------
15   USE oce_trc             ! ocean dynamics and active tracers variables
16   USE trc                 ! ocean passive tracers variables
17   USE lbclnk              ! ocean lateral boundary conditions (or mpp link)
18   USE trcbbl              ! advective passive tracers in the BBL
19   USE prtctl_trc      ! Print control for debbuging
20
21   IMPLICIT NONE
22   PRIVATE
23
24   !! * Accessibility
25   PUBLIC trc_adv_tvd    ! routine called by trcstp.F90
26
27   !! * Substitutions
28#  include "top_substitute.h90"
29   !!----------------------------------------------------------------------
30   !!   TOP 1.0 , LOCEAN-IPSL (2005)
31   !! $Header$
32   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
33   !!----------------------------------------------------------------------
34
35CONTAINS
36
37   SUBROUTINE trc_adv_tvd( kt )
38      !!----------------------------------------------------------------------
39      !!                  ***  ROUTINE trc_adv_tvd  ***
40      !!
41      !! **  Purpose :   Compute the now trend due to total advection of
42      !!       tracers and add it to the general trend of tracer equations
43      !!
44      !! **  Method  :   TVD scheme, i.e. 2nd order centered scheme with
45      !!       corrected flux (monotonic correction)
46      !!       note: - this advection scheme needs a leap-frog time scheme
47      !!
48      !! ** Action : - update tra with the now advective tracer trends
49      !!             - save the trends in trtrd ('key_trc_diatrd)
50      !!
51      !! History :
52      !!        !  95-12  (L. Mortier)  Original code
53      !!        !  00-01  (H. Loukos)  adapted to ORCA
54      !!        !  00-10  (MA Foujols E.Kestenare)  include file not routine
55      !!        !  00-12  (E. Kestenare M. Levy)  fix bug in trtrd indexes
56      !!        !  01-07  (E. Durand G. Madec)  adaptation to ORCA config
57      !!   9.0  !  02-06  (C. Ethe, G. Madec)  F90: Free form and module
58      !!----------------------------------------------------------------------
59      !! * Modules used
60#if defined key_trcbbl_adv
61      USE oce_trc            , zun => ua,  &  ! use ua as workspace
62         &                     zvn => va      ! use va as workspace
63      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zwn
64#else
65      USE oce_trc            , zun => un,  &  ! When no bbl, zun == un
66                               zvn => vn,  &  !             zvn == vn
67                               zwn => wn      !             zwn == wn
68#endif
69      !! * Arguments
70      INTEGER, INTENT( in ) ::   kt         ! ocean time-step
71
72      !! * Local declarations
73      INTEGER  ::   ji, jj, jk,jn           ! dummy loop indices
74
75      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   &
76         zti, ztu, ztv, ztw                ! temporary workspace
77
78      REAL(wp) ::   &
79         z2dtt, zbtr, zeu, zev, zew, z2, &  ! temporary scalar
80         zfp_ui, zfp_vj, zfp_wk,         &  !    "         "
81         zfm_ui, zfm_vj, zfm_wk             !    "         "
82
83#if defined key_trc_diatrd
84       REAL(wp) :: &
85          zgm, zgz
86#endif
87
88      CHARACTER (len=22) :: charout
89      !!----------------------------------------------------------------------
90
91      zti(:,:,:) = 0.e0
92
93      IF( kt == nittrc000  .AND. lwp ) THEN
94         WRITE(numout,*)
95         WRITE(numout,*) 'trc_adv_tvd : TVD advection scheme'
96         WRITE(numout,*) '~~~~~~~~~~~'
97      ENDIF
98
99      IF( neuler == 0 .AND. kt == nittrc000 ) THEN
100         z2=1.
101      ELSE
102         z2=2.
103      ENDIF
104
105#if defined key_trcbbl_adv
106      ! Advective Bottom boundary layer: add the velocity
107      ! -------------------------------------------------
108      zun(:,:,:) = un (:,:,:) - u_trc_bbl(:,:,:)
109      zvn(:,:,:) = vn (:,:,:) - v_trc_bbl(:,:,:)
110      zwn(:,:,:) = wn (:,:,:) + w_trc_bbl(:,:,:)
111#endif
112
113      DO jn = 1, jptra
114
115         ! 1. Bottom value : flux set to zero
116         ! ---------------
117         ztu(:,:,jpk) = 0.e0
118         ztv(:,:,jpk) = 0.e0
119         ztw(:,:,jpk) = 0.e0
120         zti(:,:,jpk) = 0.e0
121
122
123         ! 2. upstream advection with initial mass fluxes & intermediate update
124         ! --------------------------------------------------------------------
125         ! upstream tracer flux in the i and j direction
126         DO jk = 1, jpkm1
127            DO jj = 1, jpjm1
128               DO ji = 1, fs_jpim1   ! vector opt.
129                  zeu = 0.5 * e2u(ji,jj) * fse3u(ji,jj,jk) * zun(ji,jj,jk)
130                  zev = 0.5 * e1v(ji,jj) * fse3v(ji,jj,jk) * zvn(ji,jj,jk)
131                  ! upstream scheme
132                  zfp_ui = zeu + ABS( zeu )
133                  zfm_ui = zeu - ABS( zeu )
134                  zfp_vj = zev + ABS( zev )
135                  zfm_vj = zev - ABS( zev )
136                  ztu(ji,jj,jk) = zfp_ui * trb(ji,jj,jk,jn) + zfm_ui * trb(ji+1,jj  ,jk,jn)
137                  ztv(ji,jj,jk) = zfp_vj * trb(ji,jj,jk,jn) + zfm_vj * trb(ji  ,jj+1,jk,jn)
138               END DO
139            END DO
140         END DO
141
142         ! upstream tracer flux in the k direction
143         ! Surface value
144         IF( lk_dynspg_rl ) THEN   ! rigid lid : flux set to zero
145            ztw(:,:,1) = 0.e0
146         ELSE                      ! free surface
147            DO jj = 1, jpj
148               DO ji = 1, jpi
149                  zew = e1t(ji,jj) * e2t(ji,jj) * zwn(ji,jj,1)
150                  ztw(ji,jj,1) = zew * trb(ji,jj,1,jn)
151               END DO
152            END DO
153         ENDIF
154
155         ! Interior value
156         DO jk = 2, jpkm1
157            DO jj = 1, jpj
158               DO ji = 1, jpi
159                  zew = 0.5 * e1t(ji,jj) * e2t(ji,jj) * zwn(ji,jj,jk)
160                  zfp_wk = zew + ABS( zew )
161                  zfm_wk = zew - ABS( zew )
162                  ztw(ji,jj,jk) = zfp_wk * trb(ji,jj,jk,jn) + zfm_wk * trb(ji,jj,jk-1,jn)
163               END DO
164            END DO
165         END DO
166
167         ! total advective trend
168         DO jk = 1, jpkm1
169            DO jj = 2, jpjm1
170               DO ji = fs_2, fs_jpim1   ! vector opt.
171                  zbtr = 1./ ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
172                  zti(ji,jj,jk) = - ( ztu(ji,jj,jk) - ztu(ji-1,jj  ,jk  )   &
173                     &              + ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk  )   &
174                     &              + ztw(ji,jj,jk) - ztw(ji  ,jj  ,jk+1) ) * zbtr
175
176#if defined key_trc_diatrd
177                  IF ( luttrd(jn) ) &
178                     trtrd(ji,jj,jk,ikeep(jn),1) = trtrd(ji,jj,jk,ikeep(jn),1) -  &
179                        &                          zbtr * ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) )                     
180                  IF ( luttrd(jn) ) &
181                     trtrd(ji,jj,jk,ikeep(jn),2) = trtrd(ji,jj,jk,ikeep(jn),2) -  &
182                        &                          zbtr * ( ztv(ji,jj,jk) - ztv(ji,jj-1,jk) )                     
183                  IF ( luttrd(jn) ) &
184                     trtrd(ji,jj,jk,ikeep(jn),3) = trtrd(ji,jj,jk,ikeep(jn),3) -  &
185                        &                          zbtr * ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) )
186#endif
187               END DO
188            END DO
189         END DO
190
191
192         ! update and guess with monotonic sheme
193         DO jk = 1, jpkm1
194            z2dtt = z2 * rdttra(jk) * FLOAT(ndttrc)
195            DO jj = 2, jpjm1
196               DO ji = fs_2, fs_jpim1   ! vector opt.
197                  tra(ji,jj,jk,jn) =  tra(ji,jj,jk,jn) + zti(ji,jj,jk)
198                  zti (ji,jj,jk) = ( trb(ji,jj,jk,jn) + z2dtt * zti(ji,jj,jk) ) * tmask(ji,jj,jk)
199               END DO
200            END DO
201         END DO
202
203         ! Lateral boundary conditions on zti, zsi   (unchanged sign)
204         CALL lbc_lnk( zti, 'T', 1. )
205
206         ! 3. antidiffusive flux : high order minus low order
207         ! --------------------------------------------------
208         ! antidiffusive flux on i and j
209         DO jk = 1, jpkm1
210            DO jj = 1, jpjm1
211               DO ji = 1, fs_jpim1   ! vector opt.
212                  zeu = 0.5 * e2u(ji,jj) * fse3u(ji,jj,jk) * zun(ji,jj,jk)
213                  zev = 0.5 * e1v(ji,jj) * fse3v(ji,jj,jk) * zvn(ji,jj,jk)
214                  ztu(ji,jj,jk) = zeu * ( trn(ji,jj,jk,jn) + trn(ji+1,jj,jk,jn) ) - ztu(ji,jj,jk)
215                  ztv(ji,jj,jk) = zev * ( trn(ji,jj,jk,jn) + trn(ji,jj+1,jk,jn) ) - ztv(ji,jj,jk)
216               END DO
217            END DO
218         END DO
219
220         ! antidiffusive flux on k
221         ! Surface value
222         ztw(:,:,1) = 0.
223
224         ! Interior value
225         DO jk = 2, jpkm1
226            DO jj = 1, jpj
227               DO ji = 1, jpi
228                  zew = 0.5 * e1t(ji,jj) * e2t(ji,jj) * zwn(ji,jj,jk)
229                  ztw(ji,jj,jk) = zew * ( trn(ji,jj,jk,jn) + trn(ji,jj,jk-1,jn) ) - ztw(ji,jj,jk)
230               END DO
231            END DO
232         END DO
233
234         ! Lateral bondary conditions
235         CALL lbc_lnk( ztu, 'U', -1. )
236         CALL lbc_lnk( ztv, 'V', -1. )
237         CALL lbc_lnk( ztw, 'W',  1. )
238
239         ! 4. monotonicity algorithm
240         ! -------------------------
241         CALL nonosc( trb(:,:,:,jn), ztu, ztv, ztw, zti, z2 )
242
243
244         ! 5. final trend with corrected fluxes
245         ! ------------------------------------
246         DO jk = 1, jpkm1
247            DO jj = 2, jpjm1
248               DO ji = fs_2, fs_jpim1   ! vector opt.
249                  zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
250#if defined key_trc_diatrd
251                  IF ( luttrd(jn) ) &
252                     trtrd(ji,jj,jk,ikeep(jn),1) = trtrd(ji,jj,jk,ikeep(jn),1) -  &
253                        &                          zbtr * ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) )                     
254                  IF ( luttrd(jn) ) &
255                     trtrd(ji,jj,jk,ikeep(jn),2) = trtrd(ji,jj,jk,ikeep(jn),2) -  &
256                        &                          zbtr * ( ztv(ji,jj,jk) - ztv(ji,jj-1,jk) )                     
257                  IF ( luttrd(jn) ) &
258                     trtrd(ji,jj,jk,ikeep(jn),3) = trtrd(ji,jj,jk,ikeep(jn),3) -  &
259                        &                          zbtr * ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) )
260#endif
261                  tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn)   &
262                     &         - ( ztu(ji,jj,jk) - ztu(ji-1,jj  ,jk  )   &
263                     &           + ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk  )   &
264                     &           + ztw(ji,jj,jk) - ztw(ji  ,jj  ,jk+1) ) * zbtr
265               END DO
266            END DO
267         END DO
268         ! 6.0 convert the transport trend into advection trend
269         ! ----------------------------------------------------
270         
271#if defined key_trc_diatrd
272         DO jk = 1,jpk
273            DO jj = 2,jpjm1
274               DO  ji = 2,jpim1
275                  zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
276                  zgm = zbtr * trn(ji,jj,jk,jn) *     &
277                     &         (  zun(ji  ,jj,jk) * e2u(ji  ,jj) * fse3u(ji  ,jj,jk)    &
278                     &          - zun(ji-1,jj,jk) * e2u(ji-1,jj) * fse3u(ji-1,jj,jk) )
279                 
280                  zgz = zbtr * trn(ji,jj,jk,jn) *     &
281                     &         (  zvn(ji,jj  ,jk) * e1v(ji,jj  ) * fse3v(ji,jj  ,jk)    &
282                     &          - zvn(ji,jj-1,jk) * e1v(ji,jj-1) * fse3v(ji,jj-1,jk) )
283                 
284                  IF (luttrd(jn)) trtrd(ji,jj,jk,ikeep(jn),1) = trtrd(ji,jj,jk,ikeep(jn),1) + zgm
285                  IF (luttrd(jn)) trtrd(ji,jj,jk,ikeep(jn),2) = trtrd(ji,jj,jk,ikeep(jn),2) + zgz
286                  IF (luttrd(jn)) trtrd(ji,jj,jk,ikeep(jn),3) = trtrd(ji,jj,jk,ikeep(jn),3)    &
287                     &            - trn(ji,jj,jk,jn) * hdivn(ji,jj,jk)
288               END DO
289            END DO
290         END DO
291         
292         ! Lateral boundary conditions on trtrd:
293         
294         IF (luttrd(jn)) CALL lbc_lnk( trtrd(:,:,:,ikeep(jn),1), 'T', 1. )
295         IF (luttrd(jn)) CALL lbc_lnk( trtrd(:,:,:,ikeep(jn),2), 'T', 1. )
296         IF (luttrd(jn)) CALL lbc_lnk( trtrd(:,:,:,ikeep(jn),3), 'T', 1. )
297#endif
298
299      END DO
300
301      IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
302         WRITE(charout, FMT="('tvd - adv')")
303         CALL prt_ctl_trc_info(charout)
304         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd')
305      ENDIF
306
307   END SUBROUTINE trc_adv_tvd
308
309
310   SUBROUTINE nonosc( pbef, paa, pbb, pcc, paft, prdt )
311      !!---------------------------------------------------------------------
312      !!                    ***  ROUTINE nonosc  ***
313      !!     
314      !! **  Purpose :   compute monotonic tracer fluxes from the upstream
315      !!       scheme and the before field by a nonoscillatory algorithm
316      !!
317      !! **  Method  :   ... ???
318      !!       warning : pbef and paft must be masked, but the boundaries
319      !!       conditions on the fluxes are not necessary zalezak (1979)
320      !!       drange (1995) multi-dimensional forward-in-time and upstream-
321      !!       in-space based differencing for fluid
322      !!
323      !! History :
324      !!        !  97-04  (L. Mortier) Original code
325      !!        !  00-02  (H. Loukos)  rewritting for opa8
326      !!        !  00-10  (M.A Foujols, E. Kestenare)  lateral b.c.
327      !!        !  01-07  (E. Durand G. Madec)  adapted for T & S
328      !!   8.5  !  02-06  (G. Madec)  F90: Free form and module
329      !!----------------------------------------------------------------------
330      !! * Arguments
331      REAL(wp), INTENT( in ) ::   &
332         prdt                               ! ???
333      REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT( inout ) ::   &
334         pbef,                            & ! before field
335         paft,                            & ! after field
336         paa,                             & ! monotonic flux in the i direction
337         pbb,                             & ! monotonic flux in the j direction
338         pcc                                ! monotonic flux in the k direction
339
340      !! * Local declarations
341      INTEGER ::   ji, jj, jk               ! dummy loop indices
342      INTEGER ::   ikm1
343      REAL(wp), DIMENSION (jpi,jpj,jpk) ::   zbetup, zbetdo
344      REAL(wp) ::   zpos, zneg, zbt, za, zb, zc, zbig, zrtrn, z2dtt
345      !!----------------------------------------------------------------------
346
347      zbig = 1.e+40
348      zrtrn = 1.e-15
349      zbetup(:,:,:) = 0.e0   ;   zbetdo(:,:,:) = 0.e0
350
351      ! Search local extrema
352      ! --------------------
353      ! large negative value (-zbig) inside land
354      ! large negative value (-zbig) inside land
355      pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) )
356      paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) - zbig * ( 1.e0 - tmask(:,:,:) )
357      ! search maximum in neighbourhood
358      DO jk = 1, jpkm1
359         ikm1 = MAX(jk-1,1)
360         DO jj = 2, jpjm1
361            DO ji = fs_2, fs_jpim1   ! vector opt.
362               zbetup(ji,jj,jk) = MAX(  pbef(ji  ,jj  ,jk  ), paft(ji  ,jj  ,jk  ),   &
363                  &                     pbef(ji-1,jj  ,jk  ), pbef(ji+1,jj  ,jk  ),   &
364                  &                     paft(ji-1,jj  ,jk  ), paft(ji+1,jj  ,jk  ),   &
365                  &                     pbef(ji  ,jj-1,jk  ), pbef(ji  ,jj+1,jk  ),   &
366                  &                     paft(ji  ,jj-1,jk  ), paft(ji  ,jj+1,jk  ),   &
367                  &                     pbef(ji  ,jj  ,ikm1), pbef(ji  ,jj  ,jk+1),   &
368                  &                     paft(ji  ,jj  ,ikm1), paft(ji  ,jj  ,jk+1)  )
369            END DO
370         END DO
371      END DO
372      ! large positive value (+zbig) inside land
373      pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) )
374      paft(:,:,:) = paft(:,:,:) * tmask(:,:,:) + zbig * ( 1.e0 - tmask(:,:,:) )
375      ! search minimum in neighbourhood
376      DO jk = 1, jpkm1
377         ikm1 = MAX(jk-1,1)
378         DO jj = 2, jpjm1
379            DO ji = fs_2, fs_jpim1   ! vector opt.
380               zbetdo(ji,jj,jk) = MIN(  pbef(ji  ,jj  ,jk  ), paft(ji  ,jj  ,jk  ),   &
381                  &                     pbef(ji-1,jj  ,jk  ), pbef(ji+1,jj  ,jk  ),   &
382                  &                     paft(ji-1,jj  ,jk  ), paft(ji+1,jj  ,jk  ),   &
383                  &                     pbef(ji  ,jj-1,jk  ), pbef(ji  ,jj+1,jk  ),   &
384                  &                     paft(ji  ,jj-1,jk  ), paft(ji  ,jj+1,jk  ),   &
385                  &                     pbef(ji  ,jj  ,ikm1), pbef(ji  ,jj  ,jk+1),   &
386                  &                     paft(ji  ,jj  ,ikm1), paft(ji  ,jj  ,jk+1)  )
387            END DO
388         END DO
389      END DO
390
391      ! restore masked values to zero
392      pbef(:,:,:) = pbef(:,:,:) * tmask(:,:,:)
393      paft(:,:,:) = paft(:,:,:) * tmask(:,:,:)
394 
395
396      ! 2. Positive and negative part of fluxes and beta terms
397      ! ------------------------------------------------------
398
399      DO jk = 1, jpkm1
400         z2dtt = prdt * rdttra(jk) * FLOAT(ndttrc)
401         DO jj = 2, jpjm1
402            DO ji = fs_2, fs_jpim1   ! vector opt.
403               ! positive & negative part of the flux
404               zpos = MAX( 0., paa(ji-1,jj  ,jk  ) ) - MIN( 0., paa(ji  ,jj  ,jk  ) )   &
405                  & + MAX( 0., pbb(ji  ,jj-1,jk  ) ) - MIN( 0., pbb(ji  ,jj  ,jk  ) )   &
406                  & + MAX( 0., pcc(ji  ,jj  ,jk+1) ) - MIN( 0., pcc(ji  ,jj  ,jk  ) )
407               zneg = MAX( 0., paa(ji  ,jj  ,jk  ) ) - MIN( 0., paa(ji-1,jj  ,jk  ) )   &
408                  & + MAX( 0., pbb(ji  ,jj  ,jk  ) ) - MIN( 0., pbb(ji  ,jj-1,jk  ) )   &
409                  & + MAX( 0., pcc(ji  ,jj  ,jk  ) ) - MIN( 0., pcc(ji  ,jj  ,jk+1) )
410               ! up & down beta terms
411               zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) / z2dtt
412               zbetup(ji,jj,jk) = ( zbetup(ji,jj,jk) - paft(ji,jj,jk) ) / (zpos+zrtrn) * zbt
413               zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zbetdo(ji,jj,jk) ) / (zneg+zrtrn) * zbt
414            END DO
415         END DO
416      END DO
417
418      ! lateral boundary condition on zbetup & zbetdo   (unchanged sign)
419      CALL lbc_lnk( zbetup, 'T', 1. )
420      CALL lbc_lnk( zbetdo, 'T', 1. )
421
422
423      ! 3. monotonic flux in the i direction, i.e. paa
424      ! ----------------------------------------------
425      DO jk = 1, jpkm1
426         DO jj = 2, jpjm1
427            DO ji = fs_2, fs_jpim1   ! vector opt.
428               zc = paa(ji,jj,jk)
429               IF( zc >= 0. ) THEN
430                  za = MIN( 1., zbetdo(ji,jj,jk), zbetup(ji+1,jj,jk) )
431                  paa(ji,jj,jk) = za * zc
432               ELSE
433                  zb = MIN( 1., zbetup(ji,jj,jk), zbetdo(ji+1,jj,jk) )
434                  paa(ji,jj,jk) = zb * zc
435               ENDIF
436            END DO
437         END DO
438      END DO
439
440      ! lateral boundary condition on paa   (changed sign)
441      CALL lbc_lnk( paa, 'U', -1. )
442
443
444      ! 4. monotonic flux in the j direction, i.e. pbb
445      ! ----------------------------------------------
446      DO jk = 1, jpkm1
447         DO jj = 2, jpjm1
448            DO ji = fs_2, fs_jpim1   ! vector opt.
449               zc = pbb(ji,jj,jk)
450               IF( zc >= 0. ) THEN
451                  za = MIN( 1., zbetdo(ji,jj,jk), zbetup(ji,jj+1,jk) )
452                  pbb(ji,jj,jk) = za * zc
453               ELSE
454                  zb = MIN( 1., zbetup(ji,jj,jk), zbetdo(ji,jj+1,jk) )
455                  pbb(ji,jj,jk) = zb * zc
456               ENDIF
457            END DO
458         END DO
459      END DO
460
461      ! lateral boundary condition on pbb   (changed sign)
462      CALL lbc_lnk( pbb, 'V', -1. )
463
464
465      ! monotonic flux in the k direction, i.e. pcc
466      ! -------------------------------------------
467      DO jk = 2, jpkm1
468         DO jj = 2, jpjm1
469            DO ji = fs_2, fs_jpim1   ! vector opt.
470               zc = pcc(ji,jj,jk)
471               IF( zc >= 0. ) THEN
472                  za = MIN( 1., zbetdo(ji,jj,jk), zbetup(ji,jj,jk-1) )
473                  pcc(ji,jj,jk) = za * zc
474               ELSE
475                  zb = MIN( 1., zbetup(ji,jj,jk), zbetdo(ji,jj,jk-1) )
476                  pcc(ji,jj,jk) = zb * zc
477               ENDIF
478            END DO
479         END DO
480      END DO
481
482      ! lateral boundary condition on pcc   (unchanged sign)
483      CALL lbc_lnk( pcc, 'W', 1. )
484
485   END SUBROUTINE nonosc
486
487#else
488   !!----------------------------------------------------------------------
489   !!   Default option                                         Empty module
490   !!----------------------------------------------------------------------
491CONTAINS
492   SUBROUTINE trc_adv_tvd( kt ) 
493      INTEGER, INTENT(in) :: kt
494      WRITE(*,*) 'trc_adv_tvd: You should not have seen this print! error?', kt
495   END SUBROUTINE trc_adv_tvd
496#endif
497
498   !!======================================================================
499END MODULE trcadv_tvd
Note: See TracBrowser for help on using the repository browser.