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

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

style of all top namelist has been modified ; update modules to take it into account, see ticket:196

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