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.
traldf_iso.F90 in branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/OPA_SRC/TRA – NEMO

source: branches/2014/dev_CNRS0_NOC1_LDF/NEMOGCM/NEMO/OPA_SRC/TRA/traldf_iso.F90 @ 4616

Last change on this file since 4616 was 4616, checked in by gm, 10 years ago

#1260 : see the associated wiki page for explanation
  • Property svn:keywords set to Id
File size: 22.2 KB
Line 
1MODULE traldf_iso
2   !!======================================================================
3   !!                   ***  MODULE  traldf_iso  ***
4   !! Ocean  tracers:  horizontal component of the lateral tracer mixing trend
5   !!======================================================================
6   !! History :  OPA  ! 1994-08  (G. Madec, M. Imbard)
7   !!            8.0  ! 1997-05  (G. Madec)  split into traldf and trazdf
8   !!            NEMO ! 2002-08  (G. Madec)  Free form, F90
9   !!            1.0  ! 2005-11  (G. Madec)  merge traldf and trazdf :-)
10   !!            3.3  ! 2010-09  (C. Ethe, G. Madec) Merge TRA-TRC
11   !!            3.7  ! 2014-01  (G. Madec)  restructuration/simplification of aht/aeiv specification
12   !!            3.7  ! 2014-02  (F. Lemarie, G. Madec)  triad operator (Griffies) + Method of Stabilizing Correction
13   !!----------------------------------------------------------------------
14
15   !!----------------------------------------------------------------------
16   !!   tra_ldf_iso  : update the tracer trend with the horizontal component of a iso-neutral laplacian operator
17   !!                  and with the vertical part of the isopycnal or geopotential s-coord. operator
18   !!----------------------------------------------------------------------
19   USE oce             ! ocean dynamics and active tracers
20   USE dom_oce         ! ocean space and time domain
21   USE trc_oce         ! share passive tracers/Ocean variables
22   USE zdf_oce         ! ocean vertical physics
23   USE ldftra          ! lateral diffusion: tracer eddy coefficients
24   USE ldfslp          ! iso-neutral slopes
25   USE diaptr          ! poleward transport diagnostics
26   !
27   USE in_out_manager  ! I/O manager
28   USE iom             ! I/O library
29   USE phycst          ! physical constants
30   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
31   USE wrk_nemo        ! Memory Allocation
32   USE timing          ! Timing
33
34   IMPLICIT NONE
35   PRIVATE
36
37   PUBLIC   tra_ldf_iso   ! routine called by step.F90
38
39   !! * Substitutions
40#  include "domzgr_substitute.h90"
41#  include "vectopt_loop_substitute.h90"
42   !!----------------------------------------------------------------------
43   !! NEMO/OPA 3.7 , NEMO Consortium (2014)
44   !! $Id$
45   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
46   !!----------------------------------------------------------------------
47CONTAINS
48
49  SUBROUTINE tra_ldf_iso( kt, kit000, cdtype, pgu, pgv , pahu, pahv ,       &
50       &                                      ptb, ptbb, pta , kjpt , kpass )
51      !!----------------------------------------------------------------------
52      !!                  ***  ROUTINE tra_ldf_iso  ***
53      !!
54      !! ** Purpose :   Compute the before horizontal tracer (t & s) diffusive
55      !!      trend for a laplacian tensor (ezxcept the dz[ dz[.] ] term) and
56      !!      add it to the general trend of tracer equation.
57      !!
58      !! ** Method  :   The horizontal component of the lateral diffusive trends
59      !!      is provided by a 2nd order operator rotated along neural or geopo-
60      !!      tential surfaces to which an eddy induced advection can be added
61      !!      It is computed using before fields (forward in time) and isopyc-
62      !!      nal or geopotential slopes computed in routine ldfslp.
63      !!
64      !!      1st part :  masked horizontal derivative of T  ( di[ t ] )
65      !!      ========    with partial cell update if ln_zps=T.
66      !!
67      !!      2nd part :  horizontal fluxes of the lateral mixing operator
68      !!      ========   
69      !!         zftu =  pahu e2u*e3u/e1u di[ tb ]
70      !!               - pahu e2u*uslp    dk[ mi(mk(tb)) ]
71      !!         zftv =  pahv e1v*e3v/e2v dj[ tb ]
72      !!               - pahv e2u*vslp    dk[ mj(mk(tb)) ]
73      !!      take the horizontal divergence of the fluxes:
74      !!         difft = 1/(e1e2t*e3t) {  di-1[ zftu ] +  dj-1[ zftv ]  }
75      !!      Add this trend to the general trend (ta,sa):
76      !!         ta = ta + difft
77      !!
78      !!      3rd part: vertical trends of the lateral mixing operator
79      !!      ========  (excluding the vertical flux proportional to dk[t] )
80      !!      vertical fluxes associated with the rotated lateral mixing:
81      !!         zftw = - {  mi(mk(pahu)) * e2t*wslpi di[ mi(mk(tb)) ]
82      !!                   + mj(mk(pahv)) * e1t*wslpj dj[ mj(mk(tb)) ]  }
83      !!      take the horizontal divergence of the fluxes:
84      !!         difft = 1/(e1e2t*e3t) dk[ zftw ]
85      !!      Add this trend to the general trend (ta,sa):
86      !!         pta = pta + difft
87      !!
88      !! ** Action :   Update pta arrays with the before rotated diffusion
89      !!----------------------------------------------------------------------
90      INTEGER                              , INTENT(in   ) ::   kt         ! ocean time-step index
91      INTEGER                              , INTENT(in   ) ::   kit000     ! first time step index
92      CHARACTER(len=3)                     , INTENT(in   ) ::   cdtype     ! =TRA or TRC (tracer indicator)
93      INTEGER                              , INTENT(in   ) ::   kjpt       ! number of tracers
94      INTEGER                              , INTENT(in   ) ::   kpass      ! =1/2 first or second passage
95      REAL(wp), DIMENSION(jpi,jpj    ,kjpt), INTENT(in   ) ::   pgu , pgv  ! tracer gradient at pstep levels
96      REAL(wp), DIMENSION(jpi,jpj,jpk)     , INTENT(in   ) ::   pahu, pahv ! eddy diffusivity at u- and v-points  [m2/s]
97      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb        ! tracer (kpass=1) or laplacian of tracer (kpass=2)
98      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptbb       ! tracer (only used in kpass=2)
99      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta        ! tracer trend
100      !
101      INTEGER  ::  ji, jj, jk, jn   ! dummy loop indices
102      INTEGER  ::  ierr             ! local integer
103      REAL(wp) ::  zmsku, zahu_w, zabe1, zcof1, zcoef3   ! local scalars
104      REAL(wp) ::  zmskv, zahv_w, zabe2, zcof2, zcoef4   !   -      -
105      REAL(wp) ::  zcoef0, ze3w_2, zsign, z2dt, z1_2dt   !   -      -
106#if defined key_diaar5
107      REAL(wp)                         ::   zztmp               ! local scalar
108#endif
109      REAL(wp), POINTER, DIMENSION(:,:)   ::   zdkt, zdk1t, z2d
110      REAL(wp), POINTER, DIMENSION(:,:,:) ::   zdit, zdjt, zftu, zftv, ztfw 
111      !!----------------------------------------------------------------------
112      !
113      IF( nn_timing == 1 )  CALL timing_start('tra_ldf_iso')
114      !
115      CALL wrk_alloc( jpi, jpj,      zdkt, zdk1t, z2d ) 
116      CALL wrk_alloc( jpi, jpj, jpk, zdit, zdjt, zftu, zftv, ztfw  ) 
117      !
118      IF( kt == kit000 )  THEN
119         IF(lwp) WRITE(numout,*)
120         IF(lwp) WRITE(numout,*) 'tra_ldf_iso : rotated laplacian diffusion operator on ', cdtype
121         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
122         !
123         akz     (:,:,:) = 0._wp     
124         ah_wslp2(:,:,:) = 0._wp
125      ENDIF
126      !
127      !                                               ! set time step size (Euler/Leapfrog)
128      IF( neuler == 0 .AND. kt == nit000 ) THEN   ;   z2dt =     rdttra(1)      ! at nit000   (Euler)
129      ELSE                                        ;   z2dt = 2.* rdttra(1)      !             (Leapfrog)
130      ENDIF
131      z1_2dt = 1._wp / z2dt
132      !
133      IF( kpass == 1 ) THEN   ;   zsign =  1._wp      ! bilaplacian operator require a minus sign (eddy diffusivity >0)
134      ELSE                    ;   zsign = -1._wp
135      ENDIF
136         
137         
138      !!----------------------------------------------------------------------
139      !!   0 - calculate  ah_wslp2 and akz
140      !!----------------------------------------------------------------------
141      !
142      IF( kpass == 1 ) THEN                  !==  first pass only  ==!
143         !
144         DO jk = 2, jpkm1
145            DO jj = 2, jpjm1
146               DO ji = fs_2, fs_jpim1   ! vector opt.
147                  !
148                  zmsku = tmask(ji,jj,jk) / MAX(   umask(ji  ,jj,jk-1) + umask(ji-1,jj,jk)          &
149                     &                           + umask(ji-1,jj,jk-1) + umask(ji  ,jj,jk) , 1._wp  )
150                  zmskv = tmask(ji,jj,jk) / MAX(   vmask(ji,jj  ,jk-1) + vmask(ji,jj-1,jk)          &
151                     &                           + vmask(ji,jj-1,jk-1) + vmask(ji,jj  ,jk) , 1._wp  )
152                     !
153                  zahu_w = (   pahu(ji  ,jj,jk-1) + pahu(ji-1,jj,jk)    &
154                     &       + pahu(ji-1,jj,jk-1) + pahu(ji  ,jj,jk)  ) * zmsku
155                  zahv_w = (   pahv(ji,jj  ,jk-1) + pahv(ji,jj-1,jk)    &
156                     &       + pahv(ji,jj-1,jk-1) + pahv(ji,jj  ,jk)  ) * zmskv
157                     !
158                  ah_wslp2(ji,jj,jk) = zahu_w * wslpi(ji,jj,jk) * wslpi(ji,jj,jk)   &
159                     &               + zahv_w * wslpj(ji,jj,jk) * wslpj(ji,jj,jk)
160               END DO
161            END DO
162         END DO
163         !
164         IF( ln_traldf_msc ) THEN                ! stabilizing vertical diffusivity coefficient
165            DO jk = 2, jpkm1
166               DO jj = 2, jpjm1
167                  DO ji = fs_2, fs_jpim1
168                     akz(ji,jj,jk) = 0.25_wp * (                                                                     &
169                        &              ( pahu(ji  ,jj,jk) + pahu(ji  ,jj,jk-1) ) / ( e1u(ji  ,jj) * e1u(ji  ,jj) )   &
170                        &            + ( pahu(ji-1,jj,jk) + pahu(ji-1,jj,jk-1) ) / ( e1u(ji-1,jj) * e1u(ji-1,jj) )   &
171                        &            + ( pahv(ji,jj  ,jk) + pahv(ji,jj  ,jk-1) ) / ( e2v(ji,jj  ) * e2v(ji,jj  ) )   &
172                        &            + ( pahv(ji,jj-1,jk) + pahv(ji,jj-1,jk-1) ) / ( e2v(ji,jj-1) * e2v(ji,jj-1) )   )
173                  END DO
174               END DO
175            END DO
176            !
177            IF( ln_traldf_blp ) THEN                ! bilaplacian operator
178               DO jk = 2, jpkm1
179                  DO jj = 1, jpjm1
180                     DO ji = 1, fs_jpim1
181                        akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk)   &
182                           &          * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( fse3w(ji,jj,jk) * fse3w(ji,jj,jk) )  )
183                     END DO
184                  END DO
185               END DO
186            ELSEIF( ln_traldf_lap ) THEN              ! laplacian operator
187               DO jk = 2, jpkm1
188                  DO jj = 1, jpjm1
189                     DO ji = 1, fs_jpim1
190                        ze3w_2 = fse3w(ji,jj,jk) * fse3w(ji,jj,jk)
191                        zcoef0 = z2dt * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2  )
192                        akz(ji,jj,jk) = MAX( zcoef0 - 0.5_wp , 0._wp ) * ze3w_2 * z1_2dt
193                     END DO
194                  END DO
195               END DO
196           ENDIF
197           !
198         ELSE                                    ! 33 flux set to zero with akz=ah_wslp2 ==>> computed in full implicit
199            akz(:,:,:) = ah_wslp2(:,:,:)     
200         ENDIF
201      ENDIF
202      !
203      !                                                          ! ===========
204      DO jn = 1, kjpt                                            ! tracer loop
205         !                                                       ! ===========
206         !                                               
207         !!----------------------------------------------------------------------
208         !!   I - masked horizontal derivative
209         !!----------------------------------------------------------------------
210!!gm : bug.... why (x,:,:)?   (1,jpj,:) and (jpi,1,:) should be sufficient....
211         zdit (1,:,:) = 0._wp     ;     zdit (jpi,:,:) = 0._wp
212         zdjt (1,:,:) = 0._wp     ;     zdjt (jpi,:,:) = 0._wp
213         !!end
214
215         ! Horizontal tracer gradient
216         DO jk = 1, jpkm1
217            DO jj = 1, jpjm1
218               DO ji = 1, fs_jpim1   ! vector opt.
219                  zdit(ji,jj,jk) = ( ptb(ji+1,jj  ,jk,jn) - ptb(ji,jj,jk,jn) ) * umask(ji,jj,jk)
220                  zdjt(ji,jj,jk) = ( ptb(ji  ,jj+1,jk,jn) - ptb(ji,jj,jk,jn) ) * vmask(ji,jj,jk)
221               END DO
222            END DO
223         END DO
224         IF( ln_zps ) THEN      ! partial steps correction at the last ocean level
225            DO jj = 1, jpjm1
226               DO ji = 1, fs_jpim1   ! vector opt.
227                  zdit(ji,jj,mbku(ji,jj)) = pgu(ji,jj,jn)         
228                  zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn)     
229               END DO
230            END DO
231         ENDIF
232
233         !!----------------------------------------------------------------------
234         !!   II - horizontal trend  (full)
235         !!----------------------------------------------------------------------
236         !
237         DO jk = 1, jpkm1                                 ! Horizontal slab
238            !
239            !                             !== Vertical tracer gradient
240            zdk1t(:,:) = ( ptb(:,:,jk,jn) - ptb(:,:,jk+1,jn) ) * tmask(:,:,jk+1)     ! level jk+1
241            !
242            IF( jk == 1 ) THEN   ;   zdkt(:,:) = zdk1t(:,:)                          ! surface: zdkt(jk=1)=zdkt(jk=2)
243            ELSE                 ;   zdkt(:,:) = ( ptb(:,:,jk-1,jn) - ptb(:,:,jk,jn) ) * tmask(:,:,jk)
244            ENDIF
245
246            DO jj = 1 , jpjm1            !==  Horizontal fluxes
247               DO ji = 1, fs_jpim1   ! vector opt.
248                  zabe1 = pahu(ji,jj,jk) * e2_e1u(ji,jj) * fse3u(ji,jj,jk)
249                  zabe2 = pahv(ji,jj,jk) * e1_e2v(ji,jj) * fse3v(ji,jj,jk)
250                  !
251                  zmsku = 1. / MAX(  tmask(ji+1,jj,jk  ) + tmask(ji,jj,jk+1)   &
252                     &             + tmask(ji+1,jj,jk+1) + tmask(ji,jj,jk  ), 1. )
253                  !
254                  zmskv = 1. / MAX(  tmask(ji,jj+1,jk  ) + tmask(ji,jj,jk+1)   &
255                     &             + tmask(ji,jj+1,jk+1) + tmask(ji,jj,jk  ), 1. )
256                  !
257                  zcof1 = - pahu(ji,jj,jk) * e2u(ji,jj) * uslp(ji,jj,jk) * zmsku
258                  zcof2 = - pahv(ji,jj,jk) * e1v(ji,jj) * vslp(ji,jj,jk) * zmskv
259                  !
260                  zftu(ji,jj,jk ) = (  zabe1 * zdit(ji,jj,jk)   &
261                     &               + zcof1 * (  zdkt (ji+1,jj) + zdk1t(ji,jj)      &
262                     &                          + zdk1t(ji+1,jj) + zdkt (ji,jj)  )  ) * umask(ji,jj,jk)
263                  zftv(ji,jj,jk) = (  zabe2 * zdjt(ji,jj,jk)   &
264                     &               + zcof2 * (  zdkt (ji,jj+1) + zdk1t(ji,jj)      &
265                     &                          + zdk1t(ji,jj+1) + zdkt (ji,jj)  )  ) * vmask(ji,jj,jk)                 
266               END DO
267            END DO
268            !
269            DO jj = 2 , jpjm1          !== horizontal divergence and add to pta
270               DO ji = fs_2, fs_jpim1   ! vector opt.
271                  pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + zsign * (  zftu(ji,jj,jk) - zftu(ji-1,jj,jk)      &
272                     &                                           + zftv(ji,jj,jk) - zftv(ji,jj-1,jk)  )   &
273                     &                                        / (  e1e2t(ji,jj) * fse3t(ji,jj,jk)  )
274               END DO
275            END DO
276         END DO                                        !   End of slab 
277
278
279         !!----------------------------------------------------------------------
280         !!   III - vertical trend (full)
281         !!----------------------------------------------------------------------
282         
283         ztfw(1,:,:) = 0._wp     ;     ztfw(jpi,:,:) = 0._wp
284         
285         ! Vertical fluxes
286         ! ---------------
287         
288         ! Surface and bottom vertical fluxes set to zero
289         ztfw(:,:, 1 ) = 0._wp      ;      ztfw(:,:,jpk) = 0._wp
290         
291         ! interior (2=<jk=<jpk-1)
292         DO jk = 2, jpkm1
293            DO jj = 2, jpjm1
294               DO ji = fs_2, fs_jpim1   ! vector opt.
295                  !
296                  zmsku = tmask(ji,jj,jk) / MAX(   umask(ji  ,jj,jk-1) + umask(ji-1,jj,jk)          &
297                     &                           + umask(ji-1,jj,jk-1) + umask(ji  ,jj,jk) , 1._wp  )
298                  zmskv = tmask(ji,jj,jk) / MAX(   vmask(ji,jj  ,jk-1) + vmask(ji,jj-1,jk)          &
299                     &                           + vmask(ji,jj-1,jk-1) + vmask(ji,jj  ,jk) , 1._wp  )
300                     !
301                  zahu_w = (   pahu(ji  ,jj,jk-1) + pahu(ji-1,jj,jk)    &
302                     &       + pahu(ji-1,jj,jk-1) + pahu(ji  ,jj,jk)  ) * zmsku
303                  zahv_w = (   pahv(ji,jj  ,jk-1) + pahv(ji,jj-1,jk)    &
304                     &       + pahv(ji,jj-1,jk-1) + pahv(ji,jj  ,jk)  ) * zmskv
305                     !
306                  ah_wslp2(ji,jj,jk) = zahu_w * wslpi(ji,jj,jk) * wslpi(ji,jj,jk)   &
307                     &               + zahv_w * wslpj(ji,jj,jk) * wslpj(ji,jj,jk)
308                     !
309                  zcoef3 = - zahu_w * e2t(ji,jj) * zmsku * wslpi (ji,jj,jk)
310                  zcoef4 = - zahv_w * e1t(ji,jj) * zmskv * wslpj (ji,jj,jk)
311                  !
312                  ztfw(ji,jj,jk) = zcoef3 * (   zdit(ji  ,jj  ,jk-1) + zdit(ji-1,jj  ,jk)      &
313                     &                        + zdit(ji-1,jj  ,jk-1) + zdit(ji  ,jj  ,jk)  )   &
314                     &           + zcoef4 * (   zdjt(ji  ,jj  ,jk-1) + zdjt(ji  ,jj-1,jk)      &
315                     &                        + zdjt(ji  ,jj-1,jk-1) + zdjt(ji  ,jj  ,jk)  )
316               END DO
317            END DO
318         END DO
319         !
320         !                                !==  add the vertical 33 flux  ==!
321         IF( ln_traldf_lap ) THEN               ! laplacian case: eddy coef = ah_wslp2 - akz
322            DO jk = 2, jpkm1       
323               DO jj = 1, jpjm1
324                  DO ji = fs_2, fs_jpim1
325                     ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / fse3w(ji,jj,jk) * tmask(ji,jj,jk)   &
326                        &                            * ( ah_wslp2(ji,jj,jk) - akz(ji,jj,jk) )             &
327                        &                            * ( ptb(ji,jj,jk-1,jn) - ptb(ji,jj,jk,jn) )
328                  END DO
329               END DO
330            END DO
331            !
332         ELSE                                   ! bilaplacian
333            SELECT CASE( kpass )
334            CASE(  1  )                            ! 1st pass : eddy coef = ah_wslp2
335               DO jk = 2, jpkm1 
336                  DO jj = 1, jpjm1
337                     DO ji = fs_2, fs_jpim1
338                        ztfw(ji,jj,jk) = ztfw(ji,jj,jk)    &
339                           &           + ah_wslp2(ji,jj,jk) * e1e2t(ji,jj)   &
340                           &           * ( ptb(ji,jj,jk-1,jn) - ptb(ji,jj,jk,jn) ) * tmask(ji,jj,jk) / fse3w(ji,jj,jk)
341                     END DO
342                  END DO
343               END DO
344            CASE(  2  )                         ! 2nd pass : eddy flux = ah_wslp2 and akz applied on ptb  and ptbb gradients, resp.
345               DO jk = 2, jpkm1 
346                  DO jj = 1, jpjm1
347                     DO ji = fs_2, fs_jpim1
348                        ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / fse3w(ji,jj,jk) * tmask(ji,jj,jk)                      &
349                           &                            * (  ah_wslp2(ji,jj,jk) * ( ptb (ji,jj,jk-1,jn) - ptb (ji,jj,jk,jn) )   &
350                           &                               + akz     (ji,jj,jk) * ( ptbb(ji,jj,jk-1,jn) - ptbb(ji,jj,jk,jn) )   )
351                     END DO
352                  END DO
353               END DO
354            END SELECT
355         ENDIF
356         !         
357         DO jk = 1, jpkm1                 !==  Divergence of vertical fluxes added to pta  ==!
358            DO jj = 2, jpjm1
359               DO ji = fs_2, fs_jpim1   ! vector opt.
360                  pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + zsign * (  ztfw (ji,jj,jk) - ztfw(ji,jj,jk+1)  )   &
361                     &                                        / (  e1e2t(ji,jj) * fse3t(ji,jj,jk)  )
362               END DO
363            END DO
364         END DO
365         !
366         IF( ( kpass == 1 .AND. ln_traldf_lap ) .OR.  &     !==  first pass only (  laplacian)  ==!
367             ( kpass == 2 .AND. ln_traldf_blp ) ) THEN      !==  2nd   pass      (bilaplacian)  ==!
368            !
369            !                             ! "Poleward" diffusive heat or salt transports (T-S case only)
370            IF( cdtype == 'TRA' .AND. ln_diaptr .AND. ( MOD( kt, nn_fptr ) == 0 ) ) THEN
371               ! note sign is reversed to give down-gradient diffusive transports (#1043)
372               IF( jn == jp_tem)   htr_ldf(:) = ptr_vj( -zftv(:,:,:) )
373               IF( jn == jp_sal)   str_ldf(:) = ptr_vj( -zftv(:,:,:) )
374            ENDIF
375#if defined key_diaar5
376            !                             ! AR5 diagnostics:  vertical integrated heat transport
377            IF( cdtype == 'TRA' .AND. jn == jp_tem  ) THEN
378               z2d(:,:) = 0._wp 
379               ! note sign is reversed to give down-gradient diffusive transports (#1043)
380               zztmp = -1.0_wp * rau0 * rcp
381               DO jk = 1, jpkm1
382                  DO jj = 2, jpjm1
383                     DO ji = fs_2, fs_jpim1   ! vector opt.
384                        z2d(ji,jj) = z2d(ji,jj) + zftu(ji,jj,jk) 
385                     END DO
386                  END DO
387               END DO
388               z2d(:,:) = zztmp * z2d(:,:)
389               CALL lbc_lnk( z2d, 'U', -1. )
390               CALL iom_put( "udiff_heattr", z2d )                  ! heat transport in i-direction
391               z2d(:,:) = 0._wp 
392               DO jk = 1, jpkm1
393                  DO jj = 2, jpjm1
394                     DO ji = fs_2, fs_jpim1   ! vector opt.
395                        z2d(ji,jj) = z2d(ji,jj) + zftv(ji,jj,jk) 
396                     END DO
397                  END DO
398               END DO
399               z2d(:,:) = zztmp * z2d(:,:)
400               CALL lbc_lnk( z2d, 'V', -1. )
401               CALL iom_put( "vdiff_heattr", z2d )                  !  heat transport in i-direction
402            END IF
403#endif
404         ENDIF                                                    !== end pass selection  ==!
405         !
406         !                                                        ! ===============
407      END DO                                                      ! end tracer loop
408      !                                                           ! ===============
409      !
410      CALL wrk_dealloc( jpi, jpj,      zdkt, zdk1t, z2d ) 
411      CALL wrk_dealloc( jpi, jpj, jpk, zdit, zdjt, zftu, zftv, ztfw  ) 
412      !
413      IF( nn_timing == 1 )  CALL timing_stop('tra_ldf_iso')
414      !
415   END SUBROUTINE tra_ldf_iso
416
417   !!==============================================================================
418END MODULE traldf_iso
Note: See TracBrowser for help on using the repository browser.