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.
trazdf_iso.F90 in trunk/NEMO/OPA_SRC/TRA – NEMO

source: trunk/NEMO/OPA_SRC/TRA/trazdf_iso.F90 @ 3

Last change on this file since 3 was 3, checked in by opalod, 20 years ago

Initial revision
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 22.1 KB
Line 
1MODULE trazdf_iso
2   !!==============================================================================
3   !!                    ***  MODULE  trazdf_iso  ***
4   !! Ocean active tracers:  vertical component of the tracer mixing trend
5   !!==============================================================================
6#if defined key_ldfslp   ||   defined key_esopa
7   !!----------------------------------------------------------------------
8   !!   'key_ldfslp'                  rotation of the lateral mixing tensor
9   !!----------------------------------------------------------------------
10   !!   tra_zdf_iso  : update the tracer trend with the vertical part of
11   !!                  the isopycnal or geopotential s-coord. operator and
12   !!                  the vertical diffusion
13   !!----------------------------------------------------------------------
14   !! * Modules used
15   USE oce             ! ocean dynamics and tracers variables
16   USE dom_oce         ! ocean space and time domain variables
17   USE ldfslp          ! Make iso-neutral slopes available
18   USE ldftra_oce     ! ocean active tracers: lateral physics
19   USE zdf_oce         ! ocean vertical physics
20   USE zdfddm          ! ocean vertical physics: double diffusion
21   USE trdtra_oce     ! tracers trends diagnostics variables
22   USE in_out_manager  ! I/O manager
23   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
24
25   IMPLICIT NONE
26   PRIVATE
27
28   !! * Accessibility
29   PUBLIC tra_zdf_iso    ! routine called by step.F90
30
31   !! * Substitutions
32#  include "domzgr_substitute.h90"
33#  include "ldftra_substitute.h90"
34#  include "ldfeiv_substitute.h90"
35#  include "zdfddm_substitute.h90"
36   !!----------------------------------------------------------------------
37
38CONTAINS
39
40   SUBROUTINE tra_zdf_iso( kt )
41      !!----------------------------------------------------------------------
42      !!                  ***  ROUTINE tra_zdf_iso  ***
43      !!
44      !! ** Purpose :
45      !!     Compute the trend due to the vertical tracer diffusion inclu-
46      !!     ding the vertical component of lateral mixing (only for second
47      !!     order operator, for fourth order it is already computed and
48      !!     add to the general trend in traldf.F) and add it to the general
49      !!     trend of the tracer equations.
50      !!
51      !! ** Method :
52      !!         The vertical component of the lateral diffusive trends is
53      !!      provided by a 2nd order operator rotated along neural or geopo-
54      !!      tential surfaces to which an eddy induced advection can be added
55      !!      It is computed using before fields (forward in time) and isopyc-
56      !!      nal or geopotential slopes computed in routine ldfslp.
57      !!
58      !!      First part: vertical trends associated with the lateral mixing
59      !!      ==========  (excluding the vertical flux proportional to dk[t] )
60      !!      vertical fluxes associated with the rotated lateral mixing:
61      !!         zftw =-aht {  e2t*wslpi di[ mi(mk(tb)) ]
62      !!                     + e1t*wslpj dj[ mj(mk(tb)) ]  }
63      !!      save avt coef. resulting from vertical physics alone in zavt:
64      !!         zavt = avt
65      !!      update and save in zavt the vertical eddy viscosity coefficient:
66      !!         avt = avt + wslpi^2+wslj^2
67      !!      add vertical Eddy Induced advective fluxes ('lk_traldf_eiv=T):
68      !!         zftw = zftw + { di[aht e2u mi(wslpi)]
69      !!                    +dj[aht e1v mj(wslpj)] } mk(tb)
70      !!      take the horizontal divergence of the fluxes:
71      !!         difft = 1/(e1t*e2t*e3t) dk[ zftw ]
72      !!      Add this trend to the general trend (ta,sa):
73      !!         ta = ta + difft
74      !!
75      !!      Second part: vertical trend associated with the vertical physics
76      !!      ===========  (including the vertical flux proportional to dk[t]
77      !!                  associated with the lateral mixing, through the
78      !!                  update of avt)
79      !!      The vertical diffusion of tracers (t & s) is given by:
80      !!             difft = dz( avt dz(t) ) = 1/e3t dk+1( avt/e3w dk(t) )
81      !!      It is computed using a backward time scheme, t=ta.
82      !!      Surface and bottom boundary conditions: no diffusive flux on
83      !!      both tracers (bottom, applied through the masked field avt).
84      !!      Add this trend to the general trend ta,sa :
85      !!         ta = ta + dz( avt dz(t) )
86      !!         (sa = sa + dz( avs dz(t) ) if lk_zdfddm=T )
87      !!
88      !!      Third part: recover avt resulting from the vertical physics
89      !!      ==========  alone, for further diagnostics (for example to
90      !!                  compute the turbocline depth in diamld).
91      !!         avt = zavt
92      !!         (avs = zavs if lk_zdfddm=T )
93      !!
94      !!      'key_trdtra' defined: trend saved for futher diagnostics.
95      !!
96      !!      macro-tasked on vertical slab (jj-loop)
97      !!
98      !! ** Action :
99      !!         Update (ta,sa) arrays with the before vertical diffusion trend
100      !!         Save in (ttrd,strd) arrays the trends if 'key_diatrends' defined
101      !!
102      !! History :
103      !!   7.0  !  91-11  (G. Madec)  Original code
104      !!        !  92-06  (M. Imbard)  correction on tracer trend loops
105      !!        !  96-01  (G. Madec)  statement function for e3
106      !!        !  97-05  (G. Madec)  vertical component of isopycnal
107      !!        !  97-07  (G. Madec)  geopotential diffusion in s-coord
108      !!        !  00-08  (G. Madec)  double diffusive mixing
109      !!   8.5  !  02-08  (G. Madec)  F90: Free form and module
110      !!---------------------------------------------------------------------
111      !! * Arguments
112      INTEGER, INTENT( in ) ::   kt           ! ocean time-step index
113
114      !! * Local save
115      REAL(wp), DIMENSION(jpk), SAVE ::  &
116         z2dt
117
118      !! * Local declarations
119      INTEGER ::   ji, jj, jk                 ! dummy loop indices
120      INTEGER ::   ikst, ikenm2, ikstp1       ! temporary integers
121#if defined key_partial_steps
122      INTEGER ::   iku, ikv, ikv1             ! temporary integers
123#endif
124      REAL(wp) ::   zta, zsa
125      REAL(wp) ::   &
126         ztavg, zsavg,          &  ! ???
127         zcoef0, zcoef3,        &  ! ???
128         zcoef4, zcoeg3, zavi,  &  ! ???
129         zbtr, zmku, zmkv,      &  !
130         zuwki, zvwki, zuwk,    &  !
131         zvwk, ztav, zsav
132      REAL(wp), DIMENSION(jpi,jpk) ::   &
133         zwd, zws, zwi,         &  ! ???
134         zwx, zwy, zwz, zwt        ! ???
135      REAL(wp), DIMENSION(jpi,jpk) ::   &
136         ztfw, zdit, zdjt, zdj1t,   &
137         zsfw, zdis, zdjs, zdj1s,   &
138         zavt, zavs
139#if defined key_traldf_eiv
140      REAL(wp), DIMENSION(jpi,jpk) ::   &
141         ztfwg, zsfwg
142#endif
143      !!---------------------------------------------------------------------
144      !!  OPA 8.5, LODYC-IPSL (2002)
145      !!---------------------------------------------------------------------
146
147      IF( kt == nit000 ) THEN
148         IF(lwp) WRITE(numout,*)
149         IF(lwp) WRITE(numout,*) 'tra_zdf_iso : vertical mixing (including isopycnal component)'
150         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
151#if defined key_diaeiv
152         w_eiv(:,:,:) = 0.e0
153#endif
154      ENDIF
155
156      ! 0. Local constant initialization
157      ! --------------------------------
158
159      ztavg = 0.e0
160      zsavg = 0.e0
161      ! time step = 2 rdttra ex
162      IF( neuler == 0 .AND. kt == nit000 ) THEN
163         z2dt(:) =  rdttra(:)              ! restarting with Euler time stepping
164      ELSEIF( kt <= nit000 + 1) THEN
165         z2dt(:) = 2. * rdttra(:)          ! leapfrog
166      ENDIF
167
168      !                                                ! ===============
169      DO jj = 2, jpjm1                                 !  Vertical slab
170         !                                             ! ===============
171
172         ! I. vertical trends associated with the lateral mixing
173         ! =====================================================
174         !  (excluding the vertical flux proportional to dk[t]
175
176
177         ! I.1 horizontal tracer gradient
178         ! ------------------------------
179
180         DO jk = 1, jpkm1
181            DO ji = 1, jpim1
182               ! i-gradient of T and S at jj
183               zdit (ji,jk) = ( tb(ji+1,jj,jk)-tb(ji,jj,jk) ) * umask(ji,jj,jk)
184               zdis (ji,jk) = ( sb(ji+1,jj,jk)-sb(ji,jj,jk) ) * umask(ji,jj,jk)
185               ! j-gradient of T and S at jj
186               zdjt (ji,jk) = ( tb(ji,jj+1,jk)-tb(ji,jj,jk) ) * vmask(ji,jj,jk)
187               zdjs (ji,jk) = ( sb(ji,jj+1,jk)-sb(ji,jj,jk) ) * vmask(ji,jj,jk)
188               ! j-gradient of T and S at jj+1
189               zdj1t(ji,jk) = ( tb(ji,jj,jk)-tb(ji,jj-1,jk) ) * vmask(ji,jj-1,jk)
190               zdj1s(ji,jk) = ( sb(ji,jj,jk)-sb(ji,jj-1,jk) ) * vmask(ji,jj-1,jk)
191            END DO
192         END DO
193#  if defined key_partial_steps
194         ! partial steps correction at the bottom ocean level
195         DO ji = 1, jpim1
196            ! last ocean level
197            iku  = MIN( mbathy(ji,jj), mbathy(ji+1,jj  ) ) - 1
198            ikv  = MIN( mbathy(ji,jj), mbathy(ji  ,jj+1) ) - 1
199            ikv1 = MIN( mbathy(ji,jj), mbathy(ji  ,jj-1) ) - 1
200            ! i-gradient of T and S at jj
201            zdit (ji,iku) = gtu(ji,jj)
202            zdis (ji,iku) = gsu(ji,jj)
203            ! j-gradient of T and S at jj
204            zdjt (ji,ikv) = gtv(ji,jj) 
205            zdjs (ji,ikv) = gsv(ji,jj) 
206            ! j-gradient of T and S at jj+1
207            zdj1t(ji,ikv1)= gtv(ji,jj-1)
208            zdj1s(ji,ikv1)= gsv(ji,jj-1)
209         END DO
210#endif
211
212
213         ! I.2 Vertical fluxes
214         ! -------------------
215
216         ! Surface and bottom vertical fluxes set to zero
217         ztfw(:, 1 ) = 0.e0
218         zsfw(:, 1 ) = 0.e0
219         ztfw(:,jpk) = 0.e0
220         zsfw(:,jpk) = 0.e0
221         IF( lk_traldf_eiv ) THEN
222            ztfwg(:, 1 ) = 0.e0
223            zsfwg(:, 1 ) = 0.e0
224            ztfwg(:,jpk) = 0.e0
225            zsfwg(:,jpk) = 0.e0
226         ENDIF
227
228         ! interior (2=<jk=<jpk-1)
229         DO jk = 2, jpkm1
230            DO ji = 2, jpim1
231               zcoef0 = - fsahtw(ji,jj,jk) * tmask(ji,jj,jk)
232
233               zmku = 1./MAX( umask(ji  ,jj,jk-1) + umask(ji-1,jj,jk)   &
234                  &          +umask(ji-1,jj,jk-1) + umask(ji  ,jj,jk), 1. )
235
236               zmkv = 1./MAX( vmask(ji,jj  ,jk-1) + vmask(ji,jj-1,jk)   &
237                  &          +vmask(ji,jj-1,jk-1) + vmask(ji,jj  ,jk), 1. )
238
239               zcoef3 = zcoef0 * e2t(ji,jj) * zmku * wslpi (ji,jj,jk)
240               zcoef4 = zcoef0 * e1t(ji,jj) * zmkv * wslpj (ji,jj,jk)
241
242               ztfw(ji,jk) = zcoef3 * ( zdit (ji  ,jk-1) + zdit (ji-1,jk)     &
243                  &                    +zdit (ji-1,jk-1) + zdit (ji  ,jk) )   &
244                  &        + zcoef4 * ( zdjt (ji  ,jk-1) + zdj1t(ji  ,jk)     &
245                  &                    +zdj1t(ji  ,jk-1) + zdjt (ji  ,jk) )
246
247               zsfw(ji,jk) = zcoef3 * ( zdis (ji  ,jk-1) + zdis (ji-1,jk)     &
248                  &                    +zdis (ji-1,jk-1) + zdis (ji  ,jk) )   &
249                  &        + zcoef4 * ( zdjs (ji  ,jk-1) + zdj1s(ji  ,jk)     &
250                  &                    +zdj1s(ji  ,jk-1) + zdjs (ji  ,jk) )
251            END DO
252         END DO
253
254
255         ! I.3  update and save of avt (and avs if double diffusive mixing)
256         ! ---------------------------
257
258         DO jk = 2, jpkm1
259            DO ji = 2, jpim1
260
261               zavi = fsahtw(ji,jj,jk)*( wslpi(ji,jj,jk)*wslpi(ji,jj,jk)   &
262                  &                     +wslpj(ji,jj,jk)*wslpj(ji,jj,jk) )
263
264               ! save avt in zavt to recover avt for mixed layer depth diag.
265               zavt(ji,jk) = avt(ji,jj,jk)
266               ! add isopycnal vertical coeff. to avt
267               avt(ji,jj,jk) = avt(ji,jj,jk) + zavi
268               ! same procedure on avs if necessary
269#if defined key_zdfddm
270               ! save avs in zavs to recover avs in output files
271               zavs(ji,jk) = fsavs(ji,jj,jk)
272               ! add isopycnal vertical coeff. to avs
273               fsavs(ji,jj,jk) = fsavs(ji,jj,jk) + zavi
274#endif
275            END DO
276         END DO
277
278         !                              ! ---------------------------------------!
279         IF( lk_traldf_eiv ) THEN       ! Eddy induced vertical advective fluxes !
280            !                           ! ---------------------------------------!
281#if defined key_traldf_c2d || defined key_traldf_c3d
282            DO jk = 2, jpkm1
283               DO ji = 2, jpim1
284                  zuwki = ( wslpi(ji,jj,jk) + wslpi(ji-1,jj,jk) )   &
285                     &  * fsaeiu(ji-1,jj,jk) * e2u(ji-1,jj)*umask(ji-1,jj,jk)
286                  zuwk  = ( wslpi(ji,jj,jk) + wslpi(ji+1,jj,jk) )   &
287                     &  * fsaeiu(ji  ,jj,jk) * e2u(ji  ,jj)*umask(ji  ,jj,jk)
288                  zvwki = ( wslpj(ji,jj,jk) + wslpj(ji,jj-1,jk) )   &
289                     &  * fsaeiv(ji,jj-1,jk) * e1v(ji,jj-1)*vmask(ji,jj-1,jk)
290                  zvwk  = ( wslpj(ji,jj,jk) + wslpj(ji,jj+1,jk) )   &
291                     &  * fsaeiv(ji,jj  ,jk) * e1v(ji  ,jj)*vmask(ji  ,jj,jk)
292
293                  zcoeg3 = + 0.25 * tmask(ji,jj,jk) * ( zuwk - zuwki + zvwk - zvwki )
294   
295                  ztfwg(ji,jk) = + zcoeg3 * ( tb(ji,jj,jk) + tb(ji,jj,jk-1) )
296                  zsfwg(ji,jk) = + zcoeg3 * ( sb(ji,jj,jk) + sb(ji,jj,jk-1) )
297   
298                  ztfw(ji,jk) = ztfw(ji,jk) + ztfwg(ji,jk)
299                  zsfw(ji,jk) = zsfw(ji,jk) + zsfwg(ji,jk)
300# if defined key_diaeiv
301                  w_eiv(ji,jj,jk) = -2. *  zcoeg3 / ( e1t(ji,jj)*e2t(ji,jj) )
302# endif
303               END DO
304            END DO
305
306#else
307            DO jk = 2, jpkm1
308               DO ji = 2, jpim1
309                  zuwki = ( wslpi(ji,jj,jk) + wslpi(ji-1,jj,jk) )   &
310                     &  * e2u(ji-1,jj)*umask(ji-1,jj,jk)
311                  zuwk  = ( wslpi(ji,jj,jk) + wslpi(ji+1,jj,jk) )   &
312                     &  * e2u(ji  ,jj)*umask(ji  ,jj,jk)
313                  zvwki = ( wslpj(ji,jj,jk) + wslpj(ji,jj-1,jk) )   &
314                     &  * e1v(ji,jj-1)*vmask(ji,jj-1,jk)
315                  zvwk  = ( wslpj(ji,jj,jk) + wslpj(ji,jj+1,jk) )   &
316                     &  * e1v(ji  ,jj)*vmask(ji  ,jj,jk)
317   
318                  zcoeg3 = + 0.25 * tmask(ji,jj,jk) * fsaeiw(ji,jj,jk)   &
319                     &            * ( zuwk - zuwki + zvwk - zvwki )
320
321                  ztfwg(ji,jk) = + zcoeg3 * ( tb(ji,jj,jk) + tb(ji,jj,jk-1) )
322                  zsfwg(ji,jk) = + zcoeg3 * ( sb(ji,jj,jk) + sb(ji,jj,jk-1) )
323
324                  ztfw(ji,jk) = ztfw(ji,jk) + ztfwg(ji,jk)
325                  zsfw(ji,jk) = zsfw(ji,jk) + zsfwg(ji,jk)
326# if defined key_diaeiv
327                  w_eiv(ji,jj,jk) = -2. *  zcoeg3 / ( e1t(ji,jj)*e2t(ji,jj) )
328# endif
329               END DO
330            END DO
331#endif
332         ENDIF
333
334         ! I.5 Divergence of vertical fluxes added to the general tracer trend
335         ! -------------------------------------------------------------------
336
337         DO jk = 1, jpkm1
338            DO ji = 2, jpim1
339               zbtr =  1. / ( e1t(ji,jj)*e2t(ji,jj)*fse3t(ji,jj,jk) )
340               ztav = (  ztfw(ji,jk) - ztfw(ji,jk+1)  ) * zbtr
341               zsav = (  zsfw(ji,jk) - zsfw(ji,jk+1)  ) * zbtr
342               ta(ji,jj,jk) = ta(ji,jj,jk) + ztav
343               sa(ji,jj,jk) = sa(ji,jj,jk) + zsav
344#if defined key_trdtra || defined key_trdmld
345#   if defined key_traldf_eiv
346               ztavg = ( ztfwg(ji,jk) - ztfwg(ji,jk+1) ) * zbtr
347               zsavg = ( zsfwg(ji,jk) - zsfwg(ji,jk+1) ) * zbtr
348               !  WARNING ttrd(ji,jj,jk,6) used for vertical gent velocity trend
349               !                           not for damping !!!
350               ttrd(ji,jj,jk,6) = ztavg
351               strd(ji,jj,jk,6) = zsavg
352#   endif
353               ttrd(ji,jj,jk,4) = ztav - ztavg
354               strd(ji,jj,jk,4) = zsav - zsavg
355#endif
356            END DO
357         END DO
358         !                                             ! ===============
359      END DO                                           !   End of slab
360      !                                                ! ===============
361
362
363      IF( l_ctl .AND. lwp ) THEN         ! print mean trends (used for debugging)
364         zta = SUM( ta(2:jpim1,2:jpjm1,1:jpkm1) * tmask(2:jpim1,2:jpjm1,1:jpkm1) )
365         zsa = SUM( sa(2:jpim1,2:jpjm1,1:jpkm1) * tmask(2:jpim1,2:jpjm1,1:jpkm1) )
366         WRITE(numout,*) ' zdf 1- Ta: ', zta-t_ctl, ' Sa: ', zsa-s_ctl
367         t_ctl = zta   ;   s_ctl = zsa
368      ENDIF
369
370      !                                                ! ===============
371      DO jj = 2, jpjm1                                 !  Vertical slab
372         !                                             ! ===============
373
374         ! II. Vertical trend associated with the vertical physics
375         ! =======================================================
376         !     (including the vertical flux proportional to dk[t] associated
377         !      with the lateral mixing, through the avt update)
378         !     dk[ avt dk[ (t,s) ] ] diffusive trends
379
380
381         ! II.0 Matrix construction
382         ! ------------------------
383
384         ! Diagonal, inferior, superior
385         ! (including the bottom boundary condition via avt masked)
386         DO jk = 1, jpkm1
387            DO ji = 2, jpim1
388               zwi(ji,jk) = - z2dt(jk) * avt(ji,jj,jk  )   &
389                                       / ( fse3t(ji,jj,jk) * fse3w(ji,jj,jk  ) )
390               zws(ji,jk) = - z2dt(jk) * avt(ji,jj,jk+1)   &
391                                       / ( fse3t(ji,jj,jk) * fse3w(ji,jj,jk+1) )
392               zwd(ji,jk) = 1. - zwi(ji,jk) - zws(ji,jk)
393            END DO
394         END DO
395
396         ! Surface boudary conditions
397         DO ji = 2, jpim1
398            zwi(ji,1) = 0.e0
399            zwd(ji,1) = 1. - zws(ji,1)
400         END DO
401
402
403         ! II.1. Vertical diffusion on t
404         ! ---------------------------
405
406         ! Second member construction
407         DO jk = 1, jpkm1
408            DO ji = 2, jpim1
409               zwy(ji,jk) = tb(ji,jj,jk) + z2dt(jk) * ta(ji,jj,jk)
410            END DO
411         END DO
412
413         ! Matrix inversion from the first level
414         ikst = 1
415#   include "zdf.matrixsolver.h90"
416
417#if defined key_trdtra   ||   defined key_trdmld
418         ! Compute and save the vertical diffusive temperature trends
419         IF( l_traldf_iso ) THEN
420            DO jk = 1, jpkm1
421               DO ji = 2, jpim1
422                  zta = ( zwx(ji,jk) - tb(ji,jj,jk) ) / z2dt(jk)
423                  ttrd(ji,jj,jk,4) = zta - ta(ji,jj,jk) + ttrd(ji,jj,jk,4)
424               END DO
425            END DO
426         ELSE
427            DO jk = 1, jpkm1
428               DO ji = 2, jpim1
429                  zta = ( zwx(ji,jk) - tb(ji,jj,jk) ) / z2dt(jk)
430                  ttrd(ji,jj,jk,4) = zta - ta(ji,jj,jk)
431               END DO
432            END DO
433         ENDIF
434#endif
435
436         ! Save the masked temperature after in ta
437         ! (c a u t i o n:  temperature not its trend, Leap-frog scheme done
438         !                  it will not be done in tranxt)
439         DO jk = 1, jpkm1
440            DO ji = 2, jpim1
441               ta(ji,jj,jk) = zwx(ji,jk) * tmask(ji,jj,jk)
442            END DO
443         END DO
444
445
446         ! II.2 Vertical diffusion on s
447         ! ---------------------------
448
449#if defined key_zdfddm
450         ! Rebuild the Matrix as avt /= avs
451
452         ! Diagonal, inferior, superior
453         ! (including the bottom boundary condition via avs masked)
454         DO jk = 1, jpkm1
455            DO ji = 2, jpim1
456               zwi(ji,jk) = - z2dt(jk) * fsavs(ji,jj,jk  )   &
457                                       /( fse3t(ji,jj,jk) * fse3w(ji,jj,jk  ) )
458               zws(ji,jk) = - z2dt(jk) * fsavs(ji,jj,jk+1)   &
459                                       /( fse3t(ji,jj,jk) * fse3w(ji,jj,jk+1) )
460               zwd(ji,jk) = 1. - zwi(ji,jk) - zws(ji,jk)
461            END DO
462         END DO
463
464         ! Surface boudary conditions
465         DO ji = 2, jpim1
466            zwi(ji,1) = 0.e0
467            zwd(ji,1) = 1. - zws(ji,1)
468         END DO
469#endif
470         ! Second member construction
471         DO jk = 1, jpkm1
472            DO ji = 2, jpim1
473               zwy(ji,jk) = sb(ji,jj,jk) + z2dt(jk) * sa(ji,jj,jk)
474            END DO
475         END DO
476
477         ! Matrix inversion from the first level
478         ikst = 1
479#   include "zdf.matrixsolver.h90"
480
481#if defined key_trdtra   ||   defined key_trdmld
482         ! Compute and save the vertical diffusive salinity trends
483         IF( l_traldf_iso ) THEN
484            DO jk = 1, jpkm1
485               DO ji = 2, jpim1
486                  zsa = ( zwx(ji,jk) - sb(ji,jj,jk) ) / z2dt(jk)
487                  strd(ji,jj,jk,4) = zsa - sa(ji,jj,jk) + strd(ji,jj,jk,4)
488               END DO
489            END DO
490         ELSE
491            DO jk = 1, jpkm1
492               DO ji = 2, jpim1
493                  zsa = ( zwx(ji,jk) - sb(ji,jj,jk) ) / z2dt(jk)
494                  strd(ji,jj,jk,4) = zsa - sa(ji,jj,jk)
495               END DO
496            END DO
497         ENDIF
498#endif
499
500         ! Save the masked salinity after in sa
501         ! (c a u t i o n:  salinity not its trend, Leap-frog scheme done
502         !                  it will not be done in tranxt)
503         DO jk = 1, jpkm1
504            DO ji = 2, jpim1
505               sa(ji,jj,jk) = zwx(ji,jk)  * tmask(ji,jj,jk)
506            END DO
507         END DO
508
509
510         ! III. recover the avt (avs) resulting from vertical physics only
511         ! ===============================================================
512
513         DO jk = 2, jpkm1
514            DO ji = 2, jpim1
515               avt(ji,jj,jk) = zavt(ji,jk)
516#if defined key_zdfddm
517               fsavs(ji,jj,jk) = zavs(ji,jk)
518#endif
519            END DO
520         END DO
521
522         !                                             ! ===============
523      END DO                                           !   End of slab
524      !                                                ! ===============
525
526      IF( l_ctl .AND. lwp ) THEN         ! print mean trends (used for debugging)
527         zta = SUM( ta(2:jpim1,2:jpjm1,1:jpkm1) * tmask(2:jpim1,2:jpjm1,1:jpkm1) )
528         zsa = SUM( sa(2:jpim1,2:jpjm1,1:jpkm1) * tmask(2:jpim1,2:jpjm1,1:jpkm1) )
529         WRITE(numout,*) ' zdf 2- Ta: ', zta, ' Sa: ', zsa
530      ENDIF
531
532   END SUBROUTINE tra_zdf_iso
533
534#else
535   !!----------------------------------------------------------------------
536   !!   Dummy module               NO rotation of the lateral mixing tensor
537   !!----------------------------------------------------------------------
538CONTAINS
539   SUBROUTINE tra_zdf_iso( kt )              ! empty routine
540      WRITE(*,*) kt
541   END SUBROUTINE tra_zdf_iso
542#endif
543
544   !!==============================================================================
545END MODULE trazdf_iso
Note: See TracBrowser for help on using the repository browser.