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.
trczdf_iso_vopt.F90 in trunk/NEMO/TOP_SRC/TRP – NEMO

source: trunk/NEMO/TOP_SRC/TRP/trczdf_iso_vopt.F90 @ 403

Last change on this file since 403 was 349, checked in by opalod, 18 years ago

nemo_v1_update_031 : CT : change header names
  • Property svn:eol-style set to native
  • Property svn:executable set to *
  • Property svn:keywords set to Author Date Id Revision
File size: 24.4 KB
Line 
1MODULE trczdf_iso_vopt
2   !!==============================================================================
3   !!                 ***  MODULE  trczdf_iso_vopt  ***
4   !! Ocean passive tracers:  vertical component of the tracer mixing trend
5   !!==============================================================================
6#if defined key_passivetrc && ( defined key_ldfslp   ||   defined key_esopa )
7   !!----------------------------------------------------------------------
8   !!   'key_ldfslp'                  rotation of the lateral mixing tensor
9   !!----------------------------------------------------------------------
10   !!   trc_zdf_iso_vopt : Update the tracer trend with the vertical part of
11   !!                  the isopycnal or geopotential s-coord. operator and
12   !!                  the vertical diffusion. vector optimization, use
13   !!                  k-j-i loops.
14   !!   trc_zdf_iso  :
15   !!   trc_zdf_zdf  :
16   !!----------------------------------------------------------------------
17   !! * Modules used
18   USE oce_trc         ! ocean dynamics and tracers variables
19   USE trc             ! ocean passive tracers variables
20   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
21   USE trctrp_lec      ! passive tracers transport
22   USE prtctl_trc          ! Print control for debbuging
23
24   IMPLICIT NONE
25   PRIVATE
26
27   !! * Routine accessibility
28   PUBLIC trc_zdf_iso_vopt   !  routine called by step.F90
29
30   !! * Module variables
31   REAL(wp), DIMENSION(jpk) ::  &
32      rdttrc                          ! vertical profile of 2 x time-step
33
34   !! * Substitutions
35#  include "passivetrc_substitute.h90"
36   !!----------------------------------------------------------------------
37   !!   TOP 1.0 , LOCEAN-IPSL (2005)
38   !! $Header$
39   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
40   !!----------------------------------------------------------------------
41
42CONTAINS
43   
44   SUBROUTINE trc_zdf_iso_vopt( kt )
45      !!----------------------------------------------------------------------
46      !!                  ***  ROUTINE trc_zdf_iso_vopt  ***
47      !!
48      !! ** Purpose :
49      !! ** Method  :
50      !! ** Action  :
51      !!
52      !! History :
53      !!   8.5  !  02-08  (G. Madec)  F90: Free form and module
54      !!   9.0  !  04-03  (C. Ethe)   adapted for passive tracers
55      !!---------------------------------------------------------------------
56      !! * Arguments
57      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
58      CHARACTER (len=22) :: charout
59      !!---------------------------------------------------------------------
60
61      IF( kt == nittrc000 ) THEN
62         IF(lwp)WRITE(numout,*)
63         IF(lwp)WRITE(numout,*) 'trc_zdf_iso_vopt : vertical mixing computation'
64         IF(lwp)WRITE(numout,*) '~~~~~~~~~~~~~~~~   is  iso-neutral diffusion : implicit vertical time stepping'
65#if defined key_trcldf_eiv && defined key_diaeiv 
66         w_trc_eiv(:,:,:) = 0.e0
67#endif
68      ENDIF
69
70
71      ! I. vertical extra-diagonal part of the rotated tensor
72      ! -----------------------------------------------------
73
74      CALL trc_zdf_iso
75
76      IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
77         WRITE(charout, FMT="('zdf - 1')")
78         CALL prt_ctl_trc_info(charout)
79         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd')
80      ENDIF
81
82      ! II. vertical diffusion (including the vertical diagonal part of the rotated tensor)
83      ! ----------------------
84
85      CALL trc_zdf_zdf( kt )
86
87      IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
88         WRITE(charout, FMT="('zdf - 2')")
89         CALL prt_ctl_trc_info(charout)
90         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd')
91      ENDIF
92
93   END SUBROUTINE trc_zdf_iso_vopt
94
95
96   SUBROUTINE trc_zdf_zdf( kt )
97      !!----------------------------------------------------------------------
98      !!                 ***  ROUTINE trc_zdf_zdf  ***
99      !!                   
100      !! ** Purpose :   Compute the trend due to the vertical tracer diffusion
101      !!     including the vertical component of lateral mixing (only for 2nd
102      !!     order operator, for fourth order it is already computed and add
103      !!     to the general trend in traldf.F) and add it to the general trend
104      !!     of the tracer equations.
105      !!
106      !! ** Method  :   The vertical component of the lateral diffusive trends
107      !!      is provided by a 2nd order operator rotated along neural or geo-
108      !!      potential surfaces to which an eddy induced advection can be
109      !!      added. It is computed using before fields (forward in time) and
110      !!      isopycnal or geopotential slopes computed in routine ldfslp.
111      !!
112      !!      Second part: vertical trend associated with the vertical physics
113      !!      ===========  (including the vertical flux proportional to dk[t]
114      !!                  associated with the lateral mixing, through the
115      !!                  update of avt)
116      !!      The vertical diffusion of tracers  is given by:
117      !!             difft = dz( avt dz(t) ) = 1/e3t dk+1( avt/e3w dk(t) )
118      !!      It is computed using a backward time scheme (t=tra).
119      !!      Surface and bottom boundary conditions: no diffusive flux on
120      !!      both tracers (bottom, applied through the masked field avt).
121      !!      Add this trend to the general trend tra :
122      !!         tra = tra + dz( avt dz(t) )
123      !!         (tra = tra + dz( avs dz(t) ) if lk_trc_zdfddm=T )
124      !!
125      !!      Third part: recover avt resulting from the vertical physics
126      !!      ==========  alone, for further diagnostics (for example to
127      !!                  compute the turbocline depth in diamld).
128      !!         avt = zavt
129      !!         (avs = zavs if lk_trc_zdfddm=T )
130      !!
131      !!      'key_trdtra' defined: trend saved for futher diagnostics.
132      !!
133      !!      macro-tasked on vertical slab (jj-loop)
134      !!
135      !! ** Action  : - Update tra with before vertical diffusion trend
136      !!              - Save the trend in trtrd  ('key_trc_diatrd')
137      !!
138      !! History :
139      !!   6.0  !  90-10  (B. Blanke)  Original code
140      !!   7.0  !  91-11 (G. Madec)
141      !!        !  92-06 (M. Imbard) correction on tracer trend loops
142      !!        !  96-01 (G. Madec) statement function for e3
143      !!        !  97-05 (G. Madec) vertical component of isopycnal
144      !!        !  97-07 (G. Madec) geopotential diffusion in s-coord
145      !!        !  98-03  (L. Bopp MA Foujols) passive tracer generalisation
146      !!        !  00-05  (MA Foujols) add lbc for tracer trends
147      !!        !  00-06  (O Aumont)  correct isopycnal scheme suppress
148      !!        !                     avt multiple correction
149      !!        !  00-08  (G. Madec)  double diffusive mixing
150      !!   8.5  !  02-08  (G. Madec)  F90: Free form and module
151      !!   9.0  !  04-03  (C. Ethe )  adapted for passive tracers
152      !!---------------------------------------------------------------------
153      !! * Modules used
154      USE oce_trc, ONLY :   zwd   => ua,  &  ! ua, va used as
155                            zws   => va      ! workspace
156      !! * Arguments
157      INTEGER, INTENT( in ) ::   kt          ! ocean time-step index
158
159      !! * Local declarations
160      INTEGER ::   ji, jj, jk,jn                ! dummy loop indices
161      REAL(wp) ::   &
162         zavi, zrhs                          ! temporary scalars
163      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   &
164         zwi, zwt, zavsi                     ! temporary workspace arrays
165      REAL(wp) ::    ztra              !temporary scalars
166      !!---------------------------------------------------------------------
167
168
169      ! I. Local constant initialization
170      ! --------------------------------
171      ! ... time step = 2 rdttra ex
172      IF( ln_trcadv_cen2 .OR. ln_trcadv_tvd ) THEN
173         ! time step = 2 rdttra with Arakawa or TVD advection scheme
174         IF( neuler == 0 .AND. kt == nittrc000 ) THEN
175            rdttrc(:) =  rdttra(:) * FLOAT(ndttrc)             ! restarting with Euler time stepping
176         ELSEIF( kt <= nittrc000 + 1 ) THEN
177            rdttrc(:) = 2. * rdttra(:) * FLOAT(ndttrc)         ! leapfrog
178         ENDIF
179      ELSE
180         rdttrc(:) =  rdttra(:) * FLOAT(ndttrc)     
181      ENDIF
182
183      DO jn = 1, jptra
184         
185         zwd( 1 ,:,:)=0.e0     ;     zwd(jpi,:,:)=0.e0
186         zws( 1 ,:,:)=0.e0     ;     zws(jpi,:,:)=0.e0
187         zwi( 1 ,:,:)=0.e0     ;     zwi(jpi,:,:)=0.e0
188
189         zwt( 1 ,:,:)=0.e0     ;     zwt(jpi,:,:)=0.e0
190         zwt(  :,:,1)=0.e0     ;     zwt(:,:,jpk)= 0.e0
191         zavsi( 1 ,:,:)=0.e0   ;     zavsi(jpi,:,:)=0.e0 
192         zavsi(  :,:,1)=0.e0   ;     zavsi(:,:,jpk)=0.e0
193
194
195         ! II. Vertical trend associated with the vertical physics
196         ! =======================================================
197         !     (including the vertical flux proportional to dk[t] associated
198         !      with the lateral mixing, through the avt update)
199         !     dk[ avt dk[ (t,s) ] ] diffusive trends
200
201
202         ! II.0 Matrix construction
203         ! ------------------------       
204         ! update and save of avt (and avs if double diffusive mixing)
205         DO jk = 2, jpkm1
206            DO jj = 2, jpjm1
207               DO ji = fs_2, fs_jpim1   ! vector opt.
208                  zavi = fsahtw(ji,jj,jk) * (                 &   ! vertical mixing coef. due to lateral mixing
209                     &                           wslpi(ji,jj,jk) * wslpi(ji,jj,jk)      &
210                     &                         + wslpj(ji,jj,jk) * wslpj(ji,jj,jk)  )
211                  zavsi(ji,jj,jk) = fstravs(ji,jj,jk) + zavi        ! dd mixing: zavsi = total vertical mixing coef. on tracer
212
213               END DO
214            END DO
215         END DO
216
217
218         ! II.2 Vertical diffusion on tracer
219         ! ---------------------------========
220
221         ! Rebuild the Matrix as avt /= avs
222
223         ! Diagonal, inferior, superior  (including the bottom boundary condition via avs masked)
224         DO jk = 1, jpkm1
225            DO jj = 2, jpjm1
226               DO ji = fs_2, fs_jpim1   ! vector opt.
227                  zwi(ji,jj,jk) = - rdttrc(jk) * zavsi(ji,jj,jk  ) / ( fse3t(ji,jj,jk) * fse3w(ji,jj,jk  ) )
228                  zws(ji,jj,jk) = - rdttrc(jk) * zavsi(ji,jj,jk+1) / ( fse3t(ji,jj,jk) * fse3w(ji,jj,jk+1) )
229                  zwd(ji,jj,jk) = 1. - zwi(ji,jj,jk) - zws(ji,jj,jk)
230               END DO
231            END DO
232         END DO
233
234         ! Surface boudary conditions
235         DO jj = 2, jpjm1
236            DO ji = fs_2, fs_jpim1   ! vector opt.
237               zwi(ji,jj,1) = 0.e0
238               zwd(ji,jj,1) = 1. - zws(ji,jj,1)
239            END DO
240         END DO
241
242
243         !! Matrix inversion from the first level
244         !!----------------------------------------------------------------------
245         !   solve m.x = y  where m is a tri diagonal matrix ( jpk*jpk )
246         !
247         !        ( zwd1 zws1   0    0    0  )( zwx1 ) ( zwy1 )
248         !        ( zwi2 zwd2 zws2   0    0  )( zwx2 ) ( zwy2 )
249         !        (  0   zwi3 zwd3 zws3   0  )( zwx3 )=( zwy3 )
250         !        (        ...               )( ...  ) ( ...  )
251         !        (  0    0    0   zwik zwdk )( zwxk ) ( zwyk )
252         !
253         !   m is decomposed in the product of an upper and lower triangular
254         !   matrix
255         !   The 3 diagonal terms are in 2d arrays: zwd, zws, zwi
256         !   The second member is in 2d array zwy
257         !   The solution is in 2d array zwx
258         !   The 3d arry zwt is a work space array
259         !   zwy is used and then used as a work space array : its value is modified!
260
261         ! first recurrence:   Tk = Dk - Ik Sk-1 / Tk-1   (increasing k)
262         DO jj = 2, jpjm1
263            DO ji = fs_2, fs_jpim1
264               zwt(ji,jj,1) = zwd(ji,jj,1)
265            END DO
266         END DO
267         DO jk = 2, jpkm1
268            DO jj = 2, jpjm1
269               DO ji = fs_2, fs_jpim1
270                  zwt(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1)  /zwt(ji,jj,jk-1)
271               END DO
272            END DO
273         END DO
274
275         ! second recurrence:    Zk = Yk - Ik / Tk-1  Zk-1
276         DO jj = 2, jpjm1
277            DO ji = fs_2, fs_jpim1
278               tra(ji,jj,1,jn) = trb(ji,jj,1,jn) + rdttrc(1) * tra(ji,jj,1,jn)
279            END DO
280         END DO
281         DO jk = 2, jpkm1
282            DO jj = 2, jpjm1
283               DO ji = fs_2, fs_jpim1
284                  zrhs = trb(ji,jj,jk,jn) + rdttrc(jk) * tra(ji,jj,jk,jn)   ! zrhs=right hand side
285                  tra(ji,jj,jk,jn) = zrhs - zwi(ji,jj,jk) / zwt(ji,jj,jk-1) * tra(ji,jj,jk-1,jn)
286               END DO
287            END DO
288         END DO
289
290         ! third recurrence: Xk = (Zk - Sk Xk+1 ) / Tk
291         ! Save the masked passive tracer after in tra
292         ! (c a u t i o n: passive tracer not its trend, Leap-frog scheme done it will not be done in tranxt)
293         DO jj = 2, jpjm1
294            DO ji = fs_2, fs_jpim1
295               tra(ji,jj,jpkm1,jn) = tra(ji,jj,jpkm1,jn) / zwt(ji,jj,jpkm1) * tmask(ji,jj,jpkm1)
296            END DO
297         END DO
298         DO jk = jpk-2, 1, -1
299            DO jj = 2, jpjm1
300               DO ji = fs_2, fs_jpim1
301                  tra(ji,jj,jk,jn) = ( tra(ji,jj,jk,jn) - zws(ji,jj,jk) * tra(ji,jj,jk+1,jn) ) / zwt(ji,jj,jk) * tmask(ji,jj,jk)
302               END DO
303            END DO
304         END DO
305
306#if defined key_trc_diatrd
307         ! Compute and save the vertical diffusive passive tracer trends
308#  if defined key_trc_ldfiso 
309         DO jk = 1, jpkm1
310            DO jj = 2, jpjm1
311               DO ji = fs_2, fs_jpim1   ! vector opt.
312                  ztra = ( tra(ji,jj,jk,jn) - trb(ji,jj,jk,jn) ) / rdttrc(jk)
313                  trtrd(ji,jj,jk,jn,6) = ztra - tra(ji,jj,jk,jn) + trtrd(ji,jj,jk,jn,6)
314               END DO
315            END DO
316         END DO
317#  else
318         DO jk = 1, jpkm1
319            DO jj = 2, jpjm1
320               DO ji = fs_2, fs_jpim1   ! vector opt.
321                  ztra = ( tra(ji,jj,jk,jn) - trb(ji,jj,jk,jn) ) / rdttrc(jk)
322                  trtrd(ji,jj,jk,jn,6) = ztra - tra(ji,jj,jk,jn)
323               END DO
324            END DO
325         END DO
326#  endif
327#endif
328
329      END DO
330
331   END SUBROUTINE trc_zdf_zdf
332
333
334   SUBROUTINE trc_zdf_iso
335      !!----------------------------------------------------------------------
336      !!                  ***  ROUTINE trc_zdf_iso  ***
337      !!
338      !! ** Purpose :
339      !!     Compute the trend due to the vertical tracer diffusion inclu-
340      !!     ding the vertical component of lateral mixing (only for second
341      !!     order operator, for fourth order it is already computed and
342      !!     add to the general trend in traldf.F) and add it to the general
343      !!     trend of the tracer equations.
344      !!
345      !! ** Method :
346      !!         The vertical component of the lateral diffusive trends is
347      !!      provided by a 2nd order operator rotated along neural or geopo-
348      !!      tential surfaces to which an eddy induced advection can be added
349      !!      It is computed using before fields (forward in time) and isopyc-
350      !!      nal or geopotential slopes computed in routine ldfslp.
351      !!
352      !!      First part: vertical trends associated with the lateral mixing
353      !!      ==========  (excluding the vertical flux proportional to dk[t] )
354      !!      vertical fluxes associated with the rotated lateral mixing:
355      !!         zftw =-aht {  e2t*wslpi di[ mi(mk(trb)) ]
356      !!                     + e1t*wslpj dj[ mj(mk(trb)) ]  }
357      !!      save avt coef. resulting from vertical physics alone in zavt:
358      !!         zavt = avt
359      !!      update and save in zavt the vertical eddy viscosity coefficient:
360      !!         avt = avt + wslpi^2+wslj^2
361      !!      add vertical Eddy Induced advective fluxes (lk_traldf_eiv=T):
362      !!         zftw = zftw + { di[aht e2u mi(wslpi)]
363      !!                    +dj[aht e1v mj(wslpj)] } mk(trb)
364      !!      take the horizontal divergence of the fluxes:
365      !!         difft = 1/(e1t*e2t*e3t) dk[ zftw ]
366      !!      Add this trend to the general trend tra :
367      !!         tra = tra + difft
368      !!
369      !! ** Action :
370      !!         Update tra arrays with the before vertical diffusion trend
371      !!         Save in trtrd arrays the trends if 'key_trc_diatrd' defined
372      !!
373      !! History :
374      !!   6.0  !  90-10  (B. Blanke)  Original code
375      !!   7.0  !  91-11  (G. Madec)
376      !!        !  92-06  (M. Imbard) correction on tracer trend loops
377      !!        !  96-01  (G. Madec) statement function for e3
378      !!        !  97-05  (G. Madec) vertical component of isopycnal
379      !!        !  97-07  (G. Madec) geopotential diffusion in s-coord
380      !!        !  98-03  (L. Bopp MA Foujols) passive tracer generalisation
381      !!        !  00-05  (MA Foujols) add lbc for tracer trends
382      !!        !  00-06  (O Aumont)  correct isopycnal scheme suppress
383      !!        !                     avt multiple correction
384      !!        !  00-08  (G. Madec)  double diffusive mixing
385      !!   8.5  !  02-08  (G. Madec)  F90: Free form and module
386      !!   9.0  !  04-03  (C. Ethe )  adapted for passive tracers
387      !!---------------------------------------------------------------------
388      !! * Modules used
389      USE oce_trc, ONLY :   zwx => ua,  &  ! use ua, va as
390                            zwy => va      ! workspace arrays
391
392      !! * Local declarations
393      INTEGER ::   ji, jj, jk,jn       ! dummy loop indices
394#if defined key_partial_steps
395      INTEGER ::   iku, ikv
396#endif
397      REAL(wp) ::   &
398         ztavg,                  &  ! temporary scalars
399         zcoef0, zcoef3,         &  !    "         "
400         zcoef4,                 &  !    "         "
401         zbtr, zmku, zmkv,       &  !    "         "
402#if defined key_trcldf_eiv
403         zcoeg3,                 &  !    "         "
404         zuwki, zvwki,           &  !    "         "
405         zuwk, zvwk,             &  !    "         "
406#endif
407         ztav
408      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   &
409         zwz, zwt, ztfw             ! temporary workspace arrays
410      !!---------------------------------------------------------------------
411
412      DO jn = 1, jptra
413
414         ! 0. Local constant initialization
415         ! --------------------------------
416         ztavg = 0.e0
417
418         zwx( 1 ,:,:)=0.e0     ;     zwx(jpi,:,:)=0.e0
419         zwy( 1 ,:,:)=0.e0     ;     zwy(jpi,:,:)=0.e0
420         zwz( 1 ,:,:)=0.e0     ;     zwz(jpi,:,:)=0.e0
421         zwt( 1 ,:,:)=0.e0     ;     zwt(jpi,:,:)=0.e0
422         ztfw( 1 ,:,:)=0.e0    ;     ztfw(jpi,:,:)=0.e0
423
424         ! I. Vertical trends associated with lateral mixing
425         ! -------------------------------------------------
426         !    (excluding the vertical flux proportional to dk[t] )
427
428
429         ! I.1 horizontal tracer gradient
430         ! ------------------------------
431
432         DO jk = 1, jpkm1
433            DO jj = 1, jpjm1
434               DO ji = 1, fs_jpim1   ! vector opt.
435                  ! i-gradient of passive tracer at ji
436                  zwx (ji,jj,jk) = ( trb(ji+1,jj,jk,jn)-trb(ji,jj,jk,jn) ) * umask(ji,jj,jk)
437                  ! j-gradient of passive tracer at jj
438                  zwy (ji,jj,jk) = ( trb(ji,jj+1,jk,jn)-trb(ji,jj,jk,jn) ) * vmask(ji,jj,jk)
439               END DO
440            END DO
441         END DO
442#  if defined key_partial_steps
443         ! partial steps correction at the bottom ocean level
444         DO jj = 1, jpjm1
445            DO ji = 1, fs_jpim1   ! vector opt.
446               ! last ocean level
447               iku  = MIN( mbathy(ji,jj), mbathy(ji+1,jj  ) ) - 1
448               ikv  = MIN( mbathy(ji,jj), mbathy(ji  ,jj+1) ) - 1
449               ! i-gradient of passive tracer
450               zwx (ji,jj,iku) = gtru(ji,jj,jn)
451               ! j-gradient of passive tracer
452               zwy (ji,jj,ikv) = gtrv(ji,jj,jn) 
453            END DO
454         END DO
455#endif
456
457
458         ! I.2 Vertical fluxes
459         ! -------------------
460
461         ! Surface and bottom vertical fluxes set to zero
462         ztfw(:,:, 1 ) = 0.e0
463         ztfw(:,:,jpk) = 0.e0
464
465         ! interior (2=<jk=<jpk-1)
466         DO jk = 2, jpkm1
467            DO jj = 2, jpjm1
468               DO ji = fs_2, fs_jpim1   ! vector opt.
469                  zcoef0 = - fsahtw(ji,jj,jk) * tmask(ji,jj,jk)
470
471                  zmku = 1./MAX(   umask(ji  ,jj,jk-1) + umask(ji-1,jj,jk)      &
472                     &           + umask(ji-1,jj,jk-1) + umask(ji  ,jj,jk), 1.  )
473
474                  zmkv = 1./MAX(   vmask(ji,jj  ,jk-1) + vmask(ji,jj-1,jk)      &
475                     &           + vmask(ji,jj-1,jk-1) + vmask(ji,jj  ,jk), 1.  )
476
477                  zcoef3 = zcoef0 * e2t(ji,jj) * zmku * wslpi (ji,jj,jk)
478                  zcoef4 = zcoef0 * e1t(ji,jj) * zmkv * wslpj (ji,jj,jk)
479
480                  ztfw(ji,jj,jk) = zcoef3 * (   zwx(ji  ,jj  ,jk-1) + zwx(ji-1,jj  ,jk)      &
481                     &                        + zwx(ji-1,jj  ,jk-1) + zwx(ji  ,jj  ,jk)  )   &
482                     &           + zcoef4 * (   zwy(ji  ,jj  ,jk-1) + zwy(ji  ,jj-1,jk)      &
483                     &                        + zwy(ji  ,jj-1,jk-1) + zwy(ji  ,jj  ,jk)  )
484               END DO
485            END DO
486         END DO
487
488#if defined key_trcldf_eiv
489         !                              ! ---------------------------------------!
490         !                              ! Eddy induced vertical advective fluxes !
491         !                              ! ---------------------------------------!
492         zwx(:,:, 1 ) = 0.e0
493         zwx(:,:,jpk) = 0.e0
494
495         DO jk = 2, jpkm1
496            DO jj = 2, jpjm1
497               DO ji = fs_2, fs_jpim1   ! vector opt.
498#   if defined key_traldf_c2d || defined key_traldf_c3d
499                  zuwki = ( wslpi(ji,jj,jk) + wslpi(ji-1,jj,jk) )   &
500                     &  * fsaeitru(ji-1,jj,jk) * e2u(ji-1,jj) * umask(ji-1,jj,jk)
501                  zuwk  = ( wslpi(ji,jj,jk) + wslpi(ji+1,jj,jk) )   &
502                     &  * fsaeitru(ji  ,jj,jk) * e2u(ji  ,jj) * umask(ji  ,jj,jk)
503                  zvwki = ( wslpj(ji,jj,jk) + wslpj(ji,jj-1,jk) )   &
504                     &  * fsaeitrv(ji,jj-1,jk) * e1v(ji,jj-1) * vmask(ji,jj-1,jk)
505                  zvwk  = ( wslpj(ji,jj,jk) + wslpj(ji,jj+1,jk) )   &
506                     &  * fsaeitrv(ji,jj  ,jk) * e1v(ji  ,jj) * vmask(ji  ,jj,jk)
507
508                  zcoeg3 = + 0.25 * tmask(ji,jj,jk) * ( zuwk - zuwki + zvwk - zvwki )
509#   else
510                  zuwki = ( wslpi(ji,jj,jk) + wslpi(ji-1,jj,jk) )   &
511                     &  * e2u(ji-1,jj) * umask(ji-1,jj,jk)
512                  zuwk  = ( wslpi(ji,jj,jk) + wslpi(ji+1,jj,jk) )   &
513                     &  * e2u(ji  ,jj) * umask(ji  ,jj,jk)
514                  zvwki = ( wslpj(ji,jj,jk) + wslpj(ji,jj-1,jk) )   &
515                     &  * e1v(ji,jj-1) * vmask(ji,jj-1,jk)
516                  zvwk  = ( wslpj(ji,jj,jk) + wslpj(ji,jj+1,jk) )   &
517                     &  * e1v(ji  ,jj) * vmask(ji  ,jj,jk)
518
519                  zcoeg3 = + 0.25 * tmask(ji,jj,jk) * fsaeiw(ji,jj,jk)   &
520                     &            * ( zuwk - zuwki + zvwk - zvwki )
521#   endif
522                  zwx(ji,jj,jk) = + zcoeg3 * ( trb(ji,jj,jk,jn) + trb(ji,jj,jk-1,jn) )
523
524                  ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + zwx(ji,jj,jk)
525#   if defined key_diaeiv
526                  w_trc_eiv(ji,jj,jk) = -2. * zcoeg3 / ( e1t(ji,jj)*e2t(ji,jj) )
527#   endif
528               END DO
529            END DO
530         END DO
531#endif
532
533         ! I.5 Divergence of vertical fluxes added to the general tracer trend
534         ! -------------------------------------------------------------------
535
536         DO jk = 1, jpkm1
537            DO jj = 2, jpjm1
538               DO ji = fs_2, fs_jpim1   ! vector opt.
539                  zbtr =  1. / ( e1t(ji,jj)*e2t(ji,jj)*fse3t(ji,jj,jk) )
540                  ztav = (  ztfw(ji,jj,jk) - ztfw(ji,jj,jk+1)  ) * zbtr
541                  tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + ztav
542#if defined key_trc_diatrd
543#   if defined key_trcldf_eiv
544                  ztavg = ( zwx(ji,jj,jk) - zwx(ji,jj,jk+1) ) * zbtr
545                  !  WARNING trtrd(ji,jj,jk,7) used for vertical gent velocity trend  not for damping !!!
546                  trtrd(ji,jj,jk,jn,7) = ztavg
547#   endif
548                  trtrd(ji,jj,jk,jn,6) = ztav - ztavg
549#endif
550               END DO
551            END DO
552         END DO
553
554      END DO
555
556   END SUBROUTINE trc_zdf_iso
557
558#else
559   !!----------------------------------------------------------------------
560   !!   Dummy module :             NO rotation of the lateral mixing tensor
561   !!----------------------------------------------------------------------
562CONTAINS
563   SUBROUTINE trc_zdf_iso_vopt( kt )              ! empty routine
564      WRITE(*,*) 'trc_zdf_iso_vopt: You should not have seen this print! error?', kt
565   END SUBROUTINE trc_zdf_iso_vopt
566#endif
567
568   !!==============================================================================
569END MODULE trczdf_iso_vopt
Note: See TracBrowser for help on using the repository browser.