source: NEMO/branches/2019/dev_r12072_TOP-01_ENHANCE-11_CEthe/src/OCE/TRA/traldf_iso.F90 @ 12109

Last change on this file since 12109 was 12109, checked in by cetlod, 10 months ago

check out & merge dev_r11643_ENHANCE-11_CEthe_Shaconemo_diags branch onto dev_r12072_TOP-01_ENHANCE-11_CEthe

  • Property svn:keywords set to Id
File size: 21.3 KB
RevLine 
[3]1MODULE traldf_iso
[503]2   !!======================================================================
[457]3   !!                   ***  MODULE  traldf_iso  ***
[2528]4   !! Ocean  tracers:  horizontal component of the lateral tracer mixing trend
[503]5   !!======================================================================
[5836]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, S. Masson)  restructuration/simplification of aht/aeiv specification
12   !!             -   ! 2014-02  (F. Lemarie, G. Madec)  triad operator (Griffies) + Method of Stabilizing Correction
[503]13   !!----------------------------------------------------------------------
[5836]14
[3]15   !!----------------------------------------------------------------------
[6140]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
[3]18   !!----------------------------------------------------------------------
[6140]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
[7646]26   USE diaar5         ! AR5 diagnostics
[5836]27   !
[6140]28   USE in_out_manager ! I/O manager
29   USE iom            ! I/O library
30   USE phycst         ! physical constants
31   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
[3]32
33   IMPLICIT NONE
34   PRIVATE
35
[503]36   PUBLIC   tra_ldf_iso   ! routine called by step.F90
[3]37
[7646]38   LOGICAL  ::   l_ptr   ! flag to compute poleward transport
39   LOGICAL  ::   l_hst   ! flag to compute heat transport
40
[3]41   !! * Substitutions
42#  include "vectopt_loop_substitute.h90"
43   !!----------------------------------------------------------------------
[9598]44   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
[2528]45   !! $Id$
[10068]46   !! Software governed by the CeCILL license (see ./LICENSE)
[247]47   !!----------------------------------------------------------------------
[3]48CONTAINS
49
[5836]50  SUBROUTINE tra_ldf_iso( kt, kit000, cdtype, pahu, pahv, pgu , pgv ,   &
51      &                                                   pgui, pgvi,   &
52      &                                       ptb , ptbb, pta , kjpt, kpass )
[3]53      !!----------------------------------------------------------------------
54      !!                  ***  ROUTINE tra_ldf_iso  ***
[457]55      !!
[3]56      !! ** Purpose :   Compute the before horizontal tracer (t & s) diffusive
[457]57      !!      trend for a laplacian tensor (ezxcept the dz[ dz[.] ] term) and
58      !!      add it to the general trend of tracer equation.
[3]59      !!
60      !! ** Method  :   The horizontal component of the lateral diffusive trends
61      !!      is provided by a 2nd order operator rotated along neural or geopo-
62      !!      tential surfaces to which an eddy induced advection can be added
63      !!      It is computed using before fields (forward in time) and isopyc-
64      !!      nal or geopotential slopes computed in routine ldfslp.
65      !!
[2528]66      !!      1st part :  masked horizontal derivative of T  ( di[ t ] )
[5836]67      !!      ========    with partial cell update if ln_zps=T
68      !!                  with top     cell update if ln_isfcav
[457]69      !!
70      !!      2nd part :  horizontal fluxes of the lateral mixing operator
71      !!      ========   
[5836]72      !!         zftu =  pahu e2u*e3u/e1u di[ tb ]
73      !!               - pahu e2u*uslp    dk[ mi(mk(tb)) ]
74      !!         zftv =  pahv e1v*e3v/e2v dj[ tb ]
75      !!               - pahv e2u*vslp    dk[ mj(mk(tb)) ]
[3]76      !!      take the horizontal divergence of the fluxes:
[5836]77      !!         difft = 1/(e1e2t*e3t) {  di-1[ zftu ] +  dj-1[ zftv ]  }
[3]78      !!      Add this trend to the general trend (ta,sa):
79      !!         ta = ta + difft
80      !!
[457]81      !!      3rd part: vertical trends of the lateral mixing operator
82      !!      ========  (excluding the vertical flux proportional to dk[t] )
83      !!      vertical fluxes associated with the rotated lateral mixing:
[5836]84      !!         zftw = - {  mi(mk(pahu)) * e2t*wslpi di[ mi(mk(tb)) ]
85      !!                   + mj(mk(pahv)) * e1t*wslpj dj[ mj(mk(tb)) ]  }
[457]86      !!      take the horizontal divergence of the fluxes:
[5836]87      !!         difft = 1/(e1e2t*e3t) dk[ zftw ]
[457]88      !!      Add this trend to the general trend (ta,sa):
[2528]89      !!         pta = pta + difft
[3]90      !!
[2528]91      !! ** Action :   Update pta arrays with the before rotated diffusion
[503]92      !!----------------------------------------------------------------------
[2528]93      INTEGER                              , INTENT(in   ) ::   kt         ! ocean time-step index
[5836]94      INTEGER                              , INTENT(in   ) ::   kit000     ! first time step index
[2528]95      CHARACTER(len=3)                     , INTENT(in   ) ::   cdtype     ! =TRA or TRC (tracer indicator)
96      INTEGER                              , INTENT(in   ) ::   kjpt       ! number of tracers
[5836]97      INTEGER                              , INTENT(in   ) ::   kpass      ! =1/2 first or second passage
98      REAL(wp), DIMENSION(jpi,jpj,jpk)     , INTENT(in   ) ::   pahu, pahv ! eddy diffusivity at u- and v-points  [m2/s]
99      REAL(wp), DIMENSION(jpi,jpj    ,kjpt), INTENT(in   ) ::   pgu, pgv   ! tracer gradient at pstep levels
100      REAL(wp), DIMENSION(jpi,jpj,    kjpt), INTENT(in   ) ::   pgui, pgvi ! tracer gradient at top   levels
101      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb        ! tracer (kpass=1) or laplacian of tracer (kpass=2)
102      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptbb       ! tracer (only used in kpass=2)
103      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta        ! tracer trend
[2715]104      !
[2528]105      INTEGER  ::  ji, jj, jk, jn   ! dummy loop indices
[6140]106      INTEGER  ::  ikt
[5836]107      INTEGER  ::  ierr             ! local integer
108      REAL(wp) ::  zmsku, zahu_w, zabe1, zcof1, zcoef3   ! local scalars
109      REAL(wp) ::  zmskv, zahv_w, zabe2, zcof2, zcoef4   !   -      -
110      REAL(wp) ::  zcoef0, ze3w_2, zsign, z2dt, z1_2dt   !   -      -
[9019]111      REAL(wp), DIMENSION(jpi,jpj)     ::   zdkt, zdk1t, z2d
112      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zdit, zdjt, zftu, zftv, ztfw 
[3]113      !!----------------------------------------------------------------------
[3294]114      !
[9779]115      IF( kpass == 1 .AND. kt == kit000 )  THEN
[3]116         IF(lwp) WRITE(numout,*)
[2528]117         IF(lwp) WRITE(numout,*) 'tra_ldf_iso : rotated laplacian diffusion operator on ', cdtype
[457]118         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
[5836]119         !
[7753]120         akz     (:,:,:) = 0._wp     
121         ah_wslp2(:,:,:) = 0._wp
[3]122      ENDIF
[7646]123      !   
124      l_hst = .FALSE.
125      l_ptr = .FALSE.
[12109]126      IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) )     l_ptr = .TRUE. 
[7646]127      IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. &
128         &                        iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) )   l_hst = .TRUE.
129      !
130      !                                            ! set time step size (Euler/Leapfrog)
[6140]131      IF( neuler == 0 .AND. kt == nit000 ) THEN   ;   z2dt =     rdt      ! at nit000   (Euler)
132      ELSE                                        ;   z2dt = 2.* rdt      !             (Leapfrog)
[5836]133      ENDIF
134      z1_2dt = 1._wp / z2dt
135      !
136      IF( kpass == 1 ) THEN   ;   zsign =  1._wp      ! bilaplacian operator require a minus sign (eddy diffusivity >0)
137      ELSE                    ;   zsign = -1._wp
138      ENDIF
139         
140      !!----------------------------------------------------------------------
141      !!   0 - calculate  ah_wslp2 and akz
142      !!----------------------------------------------------------------------
143      !
144      IF( kpass == 1 ) THEN                  !==  first pass only  ==!
145         !
146         DO jk = 2, jpkm1
147            DO jj = 2, jpjm1
148               DO ji = fs_2, fs_jpim1   ! vector opt.
149                  !
[6140]150                  zmsku = wmask(ji,jj,jk) / MAX(   umask(ji  ,jj,jk-1) + umask(ji-1,jj,jk)          &
[5836]151                     &                           + umask(ji-1,jj,jk-1) + umask(ji  ,jj,jk) , 1._wp  )
[6140]152                  zmskv = wmask(ji,jj,jk) / MAX(   vmask(ji,jj  ,jk-1) + vmask(ji,jj-1,jk)          &
[5836]153                     &                           + vmask(ji,jj-1,jk-1) + vmask(ji,jj  ,jk) , 1._wp  )
154                     !
155                  zahu_w = (   pahu(ji  ,jj,jk-1) + pahu(ji-1,jj,jk)    &
156                     &       + pahu(ji-1,jj,jk-1) + pahu(ji  ,jj,jk)  ) * zmsku
157                  zahv_w = (   pahv(ji,jj  ,jk-1) + pahv(ji,jj-1,jk)    &
158                     &       + pahv(ji,jj-1,jk-1) + pahv(ji,jj  ,jk)  ) * zmskv
159                     !
160                  ah_wslp2(ji,jj,jk) = zahu_w * wslpi(ji,jj,jk) * wslpi(ji,jj,jk)   &
161                     &               + zahv_w * wslpj(ji,jj,jk) * wslpj(ji,jj,jk)
162               END DO
163            END DO
164         END DO
165         !
166         IF( ln_traldf_msc ) THEN                ! stabilizing vertical diffusivity coefficient
167            DO jk = 2, jpkm1
168               DO jj = 2, jpjm1
169                  DO ji = fs_2, fs_jpim1
170                     akz(ji,jj,jk) = 0.25_wp * (                                                                     &
171                        &              ( pahu(ji  ,jj,jk) + pahu(ji  ,jj,jk-1) ) / ( e1u(ji  ,jj) * e1u(ji  ,jj) )   &
172                        &            + ( pahu(ji-1,jj,jk) + pahu(ji-1,jj,jk-1) ) / ( e1u(ji-1,jj) * e1u(ji-1,jj) )   &
173                        &            + ( pahv(ji,jj  ,jk) + pahv(ji,jj  ,jk-1) ) / ( e2v(ji,jj  ) * e2v(ji,jj  ) )   &
174                        &            + ( pahv(ji,jj-1,jk) + pahv(ji,jj-1,jk-1) ) / ( e2v(ji,jj-1) * e2v(ji,jj-1) )   )
175                  END DO
176               END DO
177            END DO
178            !
179            IF( ln_traldf_blp ) THEN                ! bilaplacian operator
180               DO jk = 2, jpkm1
181                  DO jj = 1, jpjm1
182                     DO ji = 1, fs_jpim1
183                        akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk)   &
[6140]184                           &          * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w_n(ji,jj,jk) * e3w_n(ji,jj,jk) )  )
[5836]185                     END DO
186                  END DO
187               END DO
188            ELSEIF( ln_traldf_lap ) THEN              ! laplacian operator
189               DO jk = 2, jpkm1
190                  DO jj = 1, jpjm1
191                     DO ji = 1, fs_jpim1
[6140]192                        ze3w_2 = e3w_n(ji,jj,jk) * e3w_n(ji,jj,jk)
[5836]193                        zcoef0 = z2dt * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2  )
194                        akz(ji,jj,jk) = MAX( zcoef0 - 0.5_wp , 0._wp ) * ze3w_2 * z1_2dt
195                     END DO
196                  END DO
197               END DO
198           ENDIF
199           !
200         ELSE                                    ! 33 flux set to zero with akz=ah_wslp2 ==>> computed in full implicit
[7753]201            akz(:,:,:) = ah_wslp2(:,:,:)     
[5836]202         ENDIF
203      ENDIF
204      !
[2528]205      !                                                          ! ===========
206      DO jn = 1, kjpt                                            ! tracer loop
207         !                                                       ! ===========
208         !                                               
209         !!----------------------------------------------------------------------
210         !!   I - masked horizontal derivative
211         !!----------------------------------------------------------------------
[5836]212!!gm : bug.... why (x,:,:)?   (1,jpj,:) and (jpi,1,:) should be sufficient....
[7753]213         zdit (1,:,:) = 0._wp     ;     zdit (jpi,:,:) = 0._wp
214         zdjt (1,:,:) = 0._wp     ;     zdjt (jpi,:,:) = 0._wp
[2528]215         !!end
[3]216
[2528]217         ! Horizontal tracer gradient
218         DO jk = 1, jpkm1
219            DO jj = 1, jpjm1
220               DO ji = 1, fs_jpim1   ! vector opt.
221                  zdit(ji,jj,jk) = ( ptb(ji+1,jj  ,jk,jn) - ptb(ji,jj,jk,jn) ) * umask(ji,jj,jk)
222                  zdjt(ji,jj,jk) = ( ptb(ji  ,jj+1,jk,jn) - ptb(ji,jj,jk,jn) ) * vmask(ji,jj,jk)
223               END DO
[457]224            END DO
225         END DO
[5836]226         IF( ln_zps ) THEN      ! botton and surface ocean correction of the horizontal gradient
227            DO jj = 1, jpjm1              ! bottom correction (partial bottom cell)
[2528]228               DO ji = 1, fs_jpim1   ! vector opt.
[5120]229                  zdit(ji,jj,mbku(ji,jj)) = pgu(ji,jj,jn)         
230                  zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn)
231               END DO
232            END DO
[5836]233            IF( ln_isfcav ) THEN      ! first wet level beneath a cavity
234               DO jj = 1, jpjm1
235                  DO ji = 1, fs_jpim1   ! vector opt.
236                     IF( miku(ji,jj) > 1 )   zdit(ji,jj,miku(ji,jj)) = pgui(ji,jj,jn)         
237                     IF( mikv(ji,jj) > 1 )   zdjt(ji,jj,mikv(ji,jj)) = pgvi(ji,jj,jn)     
238                  END DO
239               END DO
240            ENDIF
[5120]241         ENDIF
[6140]242         !
[2528]243         !!----------------------------------------------------------------------
244         !!   II - horizontal trend  (full)
245         !!----------------------------------------------------------------------
[5836]246         !
247         DO jk = 1, jpkm1                                 ! Horizontal slab
248            !
249            !                             !== Vertical tracer gradient
[7753]250            zdk1t(:,:) = ( ptb(:,:,jk,jn) - ptb(:,:,jk+1,jn) ) * wmask(:,:,jk+1)     ! level jk+1
[5836]251            !
[7753]252            IF( jk == 1 ) THEN   ;   zdkt(:,:) = zdk1t(:,:)                          ! surface: zdkt(jk=1)=zdkt(jk=2)
253            ELSE                 ;   zdkt(:,:) = ( ptb(:,:,jk-1,jn) - ptb(:,:,jk,jn) ) * wmask(:,:,jk)
[5836]254            ENDIF
255            DO jj = 1 , jpjm1            !==  Horizontal fluxes
[5120]256               DO ji = 1, fs_jpim1   ! vector opt.
[6140]257                  zabe1 = pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u_n(ji,jj,jk)
258                  zabe2 = pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v_n(ji,jj,jk)
[2528]259                  !
[6140]260                  zmsku = 1. / MAX(  wmask(ji+1,jj,jk  ) + wmask(ji,jj,jk+1)   &
261                     &             + wmask(ji+1,jj,jk+1) + wmask(ji,jj,jk  ), 1. )
[2528]262                  !
[6140]263                  zmskv = 1. / MAX(  wmask(ji,jj+1,jk  ) + wmask(ji,jj,jk+1)   &
264                     &             + wmask(ji,jj+1,jk+1) + wmask(ji,jj,jk  ), 1. )
[2528]265                  !
[5836]266                  zcof1 = - pahu(ji,jj,jk) * e2u(ji,jj) * uslp(ji,jj,jk) * zmsku
267                  zcof2 = - pahv(ji,jj,jk) * e1v(ji,jj) * vslp(ji,jj,jk) * zmskv
[2528]268                  !
269                  zftu(ji,jj,jk ) = (  zabe1 * zdit(ji,jj,jk)   &
[5836]270                     &               + zcof1 * (  zdkt (ji+1,jj) + zdk1t(ji,jj)      &
271                     &                          + zdk1t(ji+1,jj) + zdkt (ji,jj)  )  ) * umask(ji,jj,jk)
[2528]272                  zftv(ji,jj,jk) = (  zabe2 * zdjt(ji,jj,jk)   &
[5836]273                     &               + zcof2 * (  zdkt (ji,jj+1) + zdk1t(ji,jj)      &
274                     &                          + zdk1t(ji,jj+1) + zdkt (ji,jj)  )  ) * vmask(ji,jj,jk)                 
[2528]275               END DO
[3]276            END DO
[5836]277            !
278            DO jj = 2 , jpjm1          !== horizontal divergence and add to pta
[5120]279               DO ji = fs_2, fs_jpim1   ! vector opt.
[5836]280                  pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + zsign * (  zftu(ji,jj,jk) - zftu(ji-1,jj,jk)      &
281                     &                                           + zftv(ji,jj,jk) - zftv(ji,jj-1,jk)  )   &
[6140]282                     &                                        * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk)
[2528]283               END DO
[3]284            END DO
[2528]285         END DO                                        !   End of slab 
[3]286
[2528]287         !!----------------------------------------------------------------------
[5836]288         !!   III - vertical trend (full)
[2528]289         !!----------------------------------------------------------------------
[6140]290         !
[11536]291         ztfw(fs_2:1,:,:) = 0._wp     ;     ztfw(jpi:fs_jpim1,:,:) = 0._wp   ! avoid to potentially manipulate NaN values
[6140]292         !
[2528]293         ! Vertical fluxes
294         ! ---------------
[6140]295         !                          ! Surface and bottom vertical fluxes set to zero
[7753]296         ztfw(:,:, 1 ) = 0._wp      ;      ztfw(:,:,jpk) = 0._wp
[2528]297         
[6140]298         DO jk = 2, jpkm1           ! interior (2=<jk=<jpk-1)
[2528]299            DO jj = 2, jpjm1
300               DO ji = fs_2, fs_jpim1   ! vector opt.
301                  !
[6140]302                  zmsku = wmask(ji,jj,jk) / MAX(   umask(ji  ,jj,jk-1) + umask(ji-1,jj,jk)          &
[5836]303                     &                           + umask(ji-1,jj,jk-1) + umask(ji  ,jj,jk) , 1._wp  )
[6140]304                  zmskv = wmask(ji,jj,jk) / MAX(   vmask(ji,jj  ,jk-1) + vmask(ji,jj-1,jk)          &
[5836]305                     &                           + vmask(ji,jj-1,jk-1) + vmask(ji,jj  ,jk) , 1._wp  )
306                     !
307                  zahu_w = (   pahu(ji  ,jj,jk-1) + pahu(ji-1,jj,jk)    &
308                     &       + pahu(ji-1,jj,jk-1) + pahu(ji  ,jj,jk)  ) * zmsku
309                  zahv_w = (   pahv(ji,jj  ,jk-1) + pahv(ji,jj-1,jk)    &
310                     &       + pahv(ji,jj-1,jk-1) + pahv(ji,jj  ,jk)  ) * zmskv
311                     !
312                  zcoef3 = - zahu_w * e2t(ji,jj) * zmsku * wslpi (ji,jj,jk)   !wslpi & j are already w-masked
313                  zcoef4 = - zahv_w * e1t(ji,jj) * zmskv * wslpj (ji,jj,jk)
[2528]314                  !
315                  ztfw(ji,jj,jk) = zcoef3 * (   zdit(ji  ,jj  ,jk-1) + zdit(ji-1,jj  ,jk)      &
316                     &                        + zdit(ji-1,jj  ,jk-1) + zdit(ji  ,jj  ,jk)  )   &
317                     &           + zcoef4 * (   zdjt(ji  ,jj  ,jk-1) + zdjt(ji  ,jj-1,jk)      &
318                     &                        + zdjt(ji  ,jj-1,jk-1) + zdjt(ji  ,jj  ,jk)  )
319               END DO
[457]320            END DO
321         END DO
[5836]322         !                                !==  add the vertical 33 flux  ==!
323         IF( ln_traldf_lap ) THEN               ! laplacian case: eddy coef = ah_wslp2 - akz
324            DO jk = 2, jpkm1       
[11536]325               DO jj = 2, jpjm1
[5836]326                  DO ji = fs_2, fs_jpim1
[6140]327                     ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w_n(ji,jj,jk) * wmask(ji,jj,jk)   &
[5836]328                        &                            * ( ah_wslp2(ji,jj,jk) - akz(ji,jj,jk) )             &
329                        &                            * ( ptb(ji,jj,jk-1,jn) - ptb(ji,jj,jk,jn) )
330                  END DO
331               END DO
332            END DO
333            !
334         ELSE                                   ! bilaplacian
335            SELECT CASE( kpass )
336            CASE(  1  )                            ! 1st pass : eddy coef = ah_wslp2
337               DO jk = 2, jpkm1 
[11536]338                  DO jj = 2, jpjm1
[5836]339                     DO ji = fs_2, fs_jpim1
340                        ztfw(ji,jj,jk) = ztfw(ji,jj,jk)    &
341                           &           + ah_wslp2(ji,jj,jk) * e1e2t(ji,jj)   &
[6140]342                           &           * ( ptb(ji,jj,jk-1,jn) - ptb(ji,jj,jk,jn) ) / e3w_n(ji,jj,jk) * wmask(ji,jj,jk)
[5836]343                     END DO
344                  END DO
345               END DO
346            CASE(  2  )                         ! 2nd pass : eddy flux = ah_wslp2 and akz applied on ptb  and ptbb gradients, resp.
347               DO jk = 2, jpkm1 
[11536]348                  DO jj = 2, jpjm1
[5836]349                     DO ji = fs_2, fs_jpim1
[6140]350                        ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w_n(ji,jj,jk) * wmask(ji,jj,jk)                      &
[5836]351                           &                            * (  ah_wslp2(ji,jj,jk) * ( ptb (ji,jj,jk-1,jn) - ptb (ji,jj,jk,jn) )   &
352                           &                               + akz     (ji,jj,jk) * ( ptbb(ji,jj,jk-1,jn) - ptbb(ji,jj,jk,jn) )   )
353                     END DO
354                  END DO
355               END DO
356            END SELECT
357         ENDIF
358         !         
359         DO jk = 1, jpkm1                 !==  Divergence of vertical fluxes added to pta  ==!
[2528]360            DO jj = 2, jpjm1
361               DO ji = fs_2, fs_jpim1   ! vector opt.
[5836]362                  pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + zsign * (  ztfw (ji,jj,jk) - ztfw(ji,jj,jk+1)  )   &
[6140]363                     &                                        * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk)
[2528]364               END DO
[457]365            END DO
366         END DO
[2528]367         !
[5836]368         IF( ( kpass == 1 .AND. ln_traldf_lap ) .OR.  &     !==  first pass only (  laplacian)  ==!
369             ( kpass == 2 .AND. ln_traldf_blp ) ) THEN      !==  2nd   pass      (bilaplacian)  ==!
370            !
371            !                             ! "Poleward" diffusive heat or salt transports (T-S case only)
[7646]372               ! note sign is reversed to give down-gradient diffusive transports )
373            IF( l_ptr )  CALL dia_ptr_hst( jn, 'ldf', -zftv(:,:,:)  )
374            !                          ! Diffusive heat transports
375            IF( l_hst )  CALL dia_ar5_hst( jn, 'ldf', -zftu(:,:,:), -zftv(:,:,:) )
[5836]376            !
377         ENDIF                                                    !== end pass selection  ==!
378         !
379         !                                                        ! ===============
380      END DO                                                      ! end tracer loop
[503]381      !
[3]382   END SUBROUTINE tra_ldf_iso
383
384   !!==============================================================================
385END MODULE traldf_iso
Note: See TracBrowser for help on using the repository browser.