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 @ 247

Last change on this file since 247 was 247, checked in by opalod, 19 years ago

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