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.
Changeset 786 for branches/dev_001_GM/NEMO/OPA_SRC/TRA/traadv_tvd.F90 – NEMO

Ignore:
Timestamp:
2008-01-10T18:11:23+01:00 (16 years ago)
Author:
gm
Message:

dev_001_GM - merge TRC-TRA on OPA only, trabbl & zpshde not done and trdmld not OK - compilation OK

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/dev_001_GM/NEMO/OPA_SRC/TRA/traadv_tvd.F90

    r719 r786  
    44   !! Ocean active tracers:  horizontal & vertical advective trend 
    55   !!============================================================================== 
    6    !! History :       !  95-12  (L. Mortier)  Original code 
    7    !!                 !  00-01  (H. Loukos)  adapted to ORCA  
    8    !!                 !  00-10  (MA Foujols E.Kestenare)  include file not routine 
    9    !!                 !  00-12  (E. Kestenare M. Levy)  fix bug in trtrd indexes 
    10    !!                 !  01-07  (E. Durand G. Madec)  adaptation to ORCA config 
     6   !! History :  7.0  !  95-12  (L. Mortier)  Original code 
     7   !!            8.0  !  00-01  (H. Loukos)  adapted to ORCA  
     8   !!             -   !  00-10  (MA Foujols E.Kestenare)  include file not routine 
     9   !!             -   !  00-12  (E. Kestenare M. Levy)  fix bug in trtrd indexes 
     10   !!             -   !  01-07  (E. Durand G. Madec)  adaptation to ORCA config 
    1111   !!            8.5  !  02-06  (G. Madec)  F90: Free form and module 
    12    !!            9.0  !  04-01  (A. de Miranda, G. Madec, J.M. Molines ): advective bbl 
    13    !!            9.0  !  08-04  (S. Cravatte) add the i-, j- & k- trends computation 
    14    !!            " "  !  05-11  (V. Garnier) Surface pressure gradient organization 
    15    !!---------------------------------------------------------------------- 
    16  
    17  
    18    !!---------------------------------------------------------------------- 
    19    !!   tra_adv_tvd  : update the tracer trend with the horizontal 
    20    !!                  and vertical advection trends using a TVD scheme 
    21    !!   nonosc       : compute monotonic tracer fluxes by a nonoscillatory 
    22    !!                  algorithm  
    23    !!---------------------------------------------------------------------- 
    24    USE oce             ! ocean dynamics and active tracers 
     12   !!   NEMO     1.0  !  04-01  (A. de Miranda, G. Madec, J.M. Molines ): advective bbl 
     13   !!             -   !  08-04  (S. Cravatte) add the i-, j- & k- trends computation 
     14   !!             -   !  05-11  (V. Garnier) Surface pressure gradient organization 
     15   !!            2.4  !  08-01  (G. Madec) Merge TRA-TRC 
     16   !!---------------------------------------------------------------------- 
     17 
     18   !!---------------------------------------------------------------------- 
     19   !!   tra_adv_tvd  : update the tracer trend with the horizontal and 
     20   !!                  vertical advection trends using a TVD scheme 
     21   !!   nonosc       : compute monotonic tracer fluxes by a nonoscillatory algorithm 
     22   !!---------------------------------------------------------------------- 
    2523   USE dom_oce         ! ocean space and time domain 
    2624   USE trdmod          ! ocean active tracers trends  
     
    2927   USE dynspg_oce      ! choice/control of key cpp for surface pressure gradient 
    3028   USE trabbl          ! Advective term of BBL 
    31    USE lib_mpp 
     29   USE lib_mpp         ! 
    3230   USE lbclnk          ! ocean lateral boundary condition (or mpp link)  
    3331   USE diaptr          ! poleward transport diagnostics 
    3432   USE prtctl          ! Print control 
    3533 
    36  
    3734   IMPLICIT NONE 
    3835   PRIVATE 
    3936 
    40    PUBLIC   tra_adv_tvd    ! routine called by step.F90 
     37   PUBLIC   tra_adv_tvd    ! routine called by traadv.F90 
    4138 
    4239   !! * Substitutions 
     
    4441#  include "vectopt_loop_substitute.h90" 
    4542   !!---------------------------------------------------------------------- 
    46    !!   OPA 9.0 , LOCEAN-IPSL (2006)  
    47    !! $Header$ 
     43   !! NEMO/OPA 2.4 , LOCEAN-IPSL (2008)  
     44   !! $Id:$  
    4845   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    4946   !!---------------------------------------------------------------------- 
     
    5148CONTAINS 
    5249 
    53    SUBROUTINE tra_adv_tvd( kt, pun, pvn, pwn ) 
     50   SUBROUTINE tra_adv_tvd( kt, cdtype, ktra, pun, pvn, pwn,   & 
     51      &                                      ptb, ptn, pta ) 
    5452      !!---------------------------------------------------------------------- 
    5553      !!                  ***  ROUTINE tra_adv_tvd  *** 
     
    6260      !!       note: - this advection scheme needs a leap-frog time scheme 
    6361      !! 
    64       !! ** Action : - update (ta,sa) with the now advective tracer trends 
     62      !! ** Action : - update pta with the now advective tracer trends 
    6563      !!             - save the trends in (ztrdt,ztrds) ('key_trdtra') 
    6664      !!---------------------------------------------------------------------- 
    67       USE oce              , ztrdt => ua   ! use ua as workspace 
    68       USE oce              , ztrds => va   ! use va as workspace 
    69       !! 
    70       INTEGER , INTENT(in)                         ::   kt    ! ocean time-step index 
    71       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pun   ! ocean velocity u-component 
    72       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pvn   ! ocean velocity v-component 
    73       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pwn   ! ocean velocity w-component 
     65      INTEGER         , INTENT(in   )                         ::   kt              ! ocean time-step index 
     66      CHARACTER(len=3), INTENT(in   )                         ::   cdtype          ! =TRA or TRC (tracer indicator) 
     67      INTEGER         , INTENT(in   )                         ::   ktra            ! tracer index 
     68      REAL(wp)        , INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   pun, pvn, pwn   ! 3 ocean velocity components 
     69      REAL(wp)        , INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   ptb, ptn        ! before and now tracer fields 
     70      REAL(wp)        , INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pta             ! tracer trend  
    7471      !! 
    7572      INTEGER  ::   ji, jj, jk              ! dummy loop indices 
    76       REAL(wp) ::                        &  ! temporary scalar 
    77          ztai, ztaj, ztak,               &  !    "         "    
    78          zsai, zsaj, zsak,               &  !    "         "    
    79          z_hdivn_x, z_hdivn_y, z_hdivn 
    80       REAL(wp) ::   & 
    81          z2dtt, zbtr, zeu, zev,          &  ! temporary scalar 
    82          zew, z2, zbtr1,                 &  ! temporary scalar 
    83          zfp_ui, zfp_vj, zfp_wk,         &  !    "         " 
    84          zfm_ui, zfm_vj, zfm_wk             !    "         " 
     73      REAL(wp) ::   ztai, ztaj, ztak 
     74      REAL(wp) ::   z2dtt, zbtr, zeu, zev   ! temporary scalar 
     75      REAL(wp) ::   zew, z2                          ! temporary scalar 
     76      REAL(wp) ::   zfp_ui, zfp_vj, zfp_wk           !    "         " 
     77      REAL(wp) ::   zfm_ui, zfm_vj, zfm_wk            !    "         " 
    8578      REAL(wp), DIMENSION (jpi,jpj,jpk) ::   zti, ztu, ztv, ztw   ! temporary workspace 
    86       REAL(wp), DIMENSION (jpi,jpj,jpk) ::   zsi, zsu, zsv, zsw   !    "           " 
    87       !!---------------------------------------------------------------------- 
    88  
    89       zti(:,:,:) = 0.e0   ;   zsi(:,:,:) = 0.e0 
     79      !!---------------------------------------------------------------------- 
     80 
     81      zti(:,:,:) = 0.e0  
    9082 
    9183      IF( kt == nit000 .AND. lwp ) THEN 
     
    10193      ! 1. Bottom value : flux set to zero 
    10294      ! --------------- 
    103       ztu(:,:,jpk) = 0.e0   ;   zsu(:,:,jpk) = 0.e0 
    104       ztv(:,:,jpk) = 0.e0   ;   zsv(:,:,jpk) = 0.e0 
    105       ztw(:,:,jpk) = 0.e0   ;   zsw(:,:,jpk) = 0.e0 
    106       zti(:,:,jpk) = 0.e0   ;   zsi(:,:,jpk) = 0.e0 
     95      ztu(:,:,jpk) = 0.e0   ;   ztv(:,:,jpk) = 0.e0 
     96      ztw(:,:,jpk) = 0.e0   ;   zti(:,:,jpk) = 0.e0 
    10797 
    10898 
     
    120110               zfp_vj = zev + ABS( zev ) 
    121111               zfm_vj = zev - ABS( zev ) 
    122                ztu(ji,jj,jk) = zfp_ui * tb(ji,jj,jk) + zfm_ui * tb(ji+1,jj  ,jk) 
    123                ztv(ji,jj,jk) = zfp_vj * tb(ji,jj,jk) + zfm_vj * tb(ji  ,jj+1,jk) 
    124                zsu(ji,jj,jk) = zfp_ui * sb(ji,jj,jk) + zfm_ui * sb(ji+1,jj  ,jk) 
    125                zsv(ji,jj,jk) = zfp_vj * sb(ji,jj,jk) + zfm_vj * sb(ji  ,jj+1,jk) 
     112               ztu(ji,jj,jk) = zfp_ui * ptb(ji,jj,jk) + zfm_ui * ptb(ji+1,jj  ,jk) 
     113               ztv(ji,jj,jk) = zfp_vj * ptb(ji,jj,jk) + zfm_vj * ptb(ji  ,jj+1,jk) 
    126114            END DO 
    127115         END DO 
     
    132120      IF( lk_dynspg_rl .OR. lk_vvl ) THEN               ! rigid lid or variable volume: flux set to zero 
    133121         ztw(:,:,1) = 0.e0 
    134          zsw(:,:,1) = 0.e0 
    135122      ELSE                                              ! free surface 
    136          DO jj = 1, jpj 
    137             DO ji = 1, jpi 
    138                zew = e1t(ji,jj) * e2t(ji,jj) * pwn(ji,jj,1) 
    139                ztw(ji,jj,1) = zew * tb(ji,jj,1) 
    140                zsw(ji,jj,1) = zew * sb(ji,jj,1) 
    141             END DO 
    142          END DO 
     123         ztw(:,:,1) = e1t(:,:) * e2t(:,:) * pwn(:,:,1) * ptb(:,:,1) 
    143124      ENDIF 
    144125 
     
    150131               zfp_wk = zew + ABS( zew ) 
    151132               zfm_wk = zew - ABS( zew ) 
    152                ztw(ji,jj,jk) = zfp_wk * tb(ji,jj,jk) + zfm_wk * tb(ji,jj,jk-1) 
    153                zsw(ji,jj,jk) = zfp_wk * sb(ji,jj,jk) + zfm_wk * sb(ji,jj,jk-1) 
     133               ztw(ji,jj,jk) = zfp_wk * ptb(ji,jj,jk) + zfm_wk * ptb(ji,jj,jk-1) 
    154134            END DO 
    155135         END DO 
     
    165145               ztaj = - ( ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk  ) ) * zbtr 
    166146               ztak = - ( ztw(ji,jj,jk) - ztw(ji  ,jj  ,jk+1) ) * zbtr 
    167                zsai = - ( zsu(ji,jj,jk) - zsu(ji-1,jj  ,jk  ) ) * zbtr 
    168                zsaj = - ( zsv(ji,jj,jk) - zsv(ji  ,jj-1,jk  ) ) * zbtr 
    169                zsak = - ( zsw(ji,jj,jk) - zsw(ji  ,jj  ,jk+1) ) * zbtr 
    170147               ! total intermediate advective trends 
    171148               zti(ji,jj,jk) = ztai + ztaj + ztak 
    172                zsi(ji,jj,jk) = zsai + zsaj + zsak 
    173             END DO 
    174          END DO 
    175       END DO 
    176  
    177  
    178       ! Save the intermediate i / j / k advective trends for diagnostics 
    179       ! ------------------------------------------------------------------- 
    180       ! Warning : We should use zun instead of un in the computations below, but we 
    181       ! also use hdivn which is computed with un, vn (check ???). So we use un, vn 
    182       ! for consistency. Results are therefore approximate with key_trabbl_adv. 
    183  
     149            END DO 
     150         END DO 
     151      END DO 
     152 
     153      !  Save the horizontal advective trends for diagnostic 
     154      ! ----------------------------------------------------- 
    184155      IF( l_trdtra ) THEN 
    185          ztrdt(:,:,:) = 0.e0   ;   ztrds(:,:,:) = 0.e0 
    186          !  
    187          ! T/S ZONAL advection trends 
    188          DO jk = 1, jpkm1 
    189             DO jj = 2, jpjm1 
    190                DO ji = fs_2, fs_jpim1   ! vector opt. 
    191                   zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
    192                   ztrdt(ji,jj,jk) = - ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) ) * zbtr 
    193                   ztrds(ji,jj,jk) = - ( zsu(ji,jj,jk) - zsu(ji-1,jj,jk) ) * zbtr 
    194                END DO 
    195             END DO 
    196          END DO 
    197          CALL trd_mod(ztrdt, ztrds, jptra_trd_xad, 'TRA', kt)    ! save the trends 
    198          ! 
    199          ! T/S MERIDIONAL advection trends 
    200          DO jk = 1, jpkm1 
    201             DO jj = 2, jpjm1 
    202                DO ji = fs_2, fs_jpim1   ! vector opt. 
    203                   zbtr      = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
    204                   ztrdt(ji,jj,jk) = - ( ztv(ji,jj,jk) - ztv(ji,jj-1,jk) ) * zbtr 
    205                   ztrds(ji,jj,jk) = - ( zsv(ji,jj,jk) - zsv(ji,jj-1,jk) ) * zbtr 
    206                END DO 
    207             END DO 
    208          END DO 
    209          CALL trd_mod(ztrdt, ztrds, jptra_trd_yad, 'TRA', kt)     ! save the trends 
    210          ! 
    211          ! T/S VERTICAL advection trends 
    212          DO jk = 1, jpkm1 
    213             DO jj = 2, jpjm1 
    214                DO ji = fs_2, fs_jpim1   ! vector opt.          
    215                   zbtr      = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
    216                   ztrdt(ji,jj,jk) = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * zbtr 
    217                   ztrds(ji,jj,jk) = - ( zsw(ji,jj,jk) - zsw(ji,jj,jk+1) ) * zbtr 
    218                END DO 
    219             END DO 
    220          END DO 
    221          CALL trd_mod(ztrdt, ztrds, jptra_trd_zad, 'TRA', kt)     ! save the trends 
    222          ! 
     156         CALL trd_tra_adv( kt, ktra, jpt_trd_xad, cdtype, ztu, pun, ptn ) 
     157         CALL trd_tra_adv( kt, ktra, jpt_trd_yad, cdtype, ztv, pvn, ptn ) 
     158         CALL trd_tra_adv( kt, ktra, jpt_trd_zad, cdtype, ztw, pwn, ptn ) 
    223159      ENDIF 
    224160 
     
    228164         DO jj = 2, jpjm1 
    229165            DO ji = fs_2, fs_jpim1   ! vector opt. 
    230                ta(ji,jj,jk) =  ta(ji,jj,jk) + zti(ji,jj,jk) 
    231                sa(ji,jj,jk) =  sa(ji,jj,jk) + zsi(ji,jj,jk) 
    232                zti (ji,jj,jk) = ( tb(ji,jj,jk) + z2dtt * zti(ji,jj,jk) ) * tmask(ji,jj,jk) 
    233                zsi (ji,jj,jk) = ( sb(ji,jj,jk) + z2dtt * zsi(ji,jj,jk) ) * tmask(ji,jj,jk) 
    234             END DO 
    235          END DO 
    236       END DO 
    237  
    238       ! Lateral boundary conditions on zti, zsi   (unchanged sign) 
     166               pta(ji,jj,jk) =  pta(ji,jj,jk) + zti(ji,jj,jk) 
     167               zti (ji,jj,jk) = ( ptb(ji,jj,jk) + z2dtt * zti(ji,jj,jk) ) * tmask(ji,jj,jk) 
     168            END DO 
     169         END DO 
     170      END DO 
     171 
     172      ! Lateral boundary conditions on zti   (unchanged sign) 
    239173      CALL lbc_lnk( zti, 'T', 1. ) 
    240       CALL lbc_lnk( zsi, 'T', 1. ) 
    241174 
    242175 
     
    249182               zeu = 0.5 * e2u(ji,jj) * fse3u(ji,jj,jk) * pun(ji,jj,jk) 
    250183               zev = 0.5 * e1v(ji,jj) * fse3v(ji,jj,jk) * pvn(ji,jj,jk) 
    251                ztu(ji,jj,jk) = zeu * ( tn(ji,jj,jk) + tn(ji+1,jj,jk) ) - ztu(ji,jj,jk) 
    252                zsu(ji,jj,jk) = zeu * ( sn(ji,jj,jk) + sn(ji+1,jj,jk) ) - zsu(ji,jj,jk) 
    253                ztv(ji,jj,jk) = zev * ( tn(ji,jj,jk) + tn(ji,jj+1,jk) ) - ztv(ji,jj,jk) 
    254                zsv(ji,jj,jk) = zev * ( sn(ji,jj,jk) + sn(ji,jj+1,jk) ) - zsv(ji,jj,jk) 
     184               ztu(ji,jj,jk) = zeu * ( ptn(ji,jj,jk) + ptn(ji+1,jj,jk) ) - ztu(ji,jj,jk) 
     185               ztv(ji,jj,jk) = zev * ( ptn(ji,jj,jk) + ptn(ji,jj+1,jk) ) - ztv(ji,jj,jk) 
    255186            END DO 
    256187         END DO 
     
    260191      ! Surface value 
    261192      ztw(:,:,1) = 0.e0 
    262       zsw(:,:,1) = 0.e0 
    263193 
    264194      ! Interior value 
     
    267197            DO ji = 1, jpi 
    268198               zew = 0.5 * e1t(ji,jj) * e2t(ji,jj) * pwn(ji,jj,jk) 
    269                ztw(ji,jj,jk) = zew * ( tn(ji,jj,jk) + tn(ji,jj,jk-1) ) - ztw(ji,jj,jk) 
    270                zsw(ji,jj,jk) = zew * ( sn(ji,jj,jk) + sn(ji,jj,jk-1) ) - zsw(ji,jj,jk) 
     199               ztw(ji,jj,jk) = zew * ( ptn(ji,jj,jk) + ptn(ji,jj,jk-1) ) - ztw(ji,jj,jk) 
    271200            END DO 
    272201         END DO 
     
    274203 
    275204      ! Lateral bondary conditions 
    276       CALL lbc_lnk( ztu, 'U', -1. )   ;   CALL lbc_lnk( zsu, 'U', -1. ) 
    277       CALL lbc_lnk( ztv, 'V', -1. )   ;   CALL lbc_lnk( zsv, 'V', -1. ) 
    278       CALL lbc_lnk( ztw, 'W',  1. )   ;   CALL lbc_lnk( zsw, 'W',  1. ) 
     205      CALL lbc_lnk( ztu, 'U', -1. ) 
     206      CALL lbc_lnk( ztv, 'V', -1. ) 
     207      CALL lbc_lnk( ztw, 'W',  1. ) 
    279208 
    280209      ! 4. monotonicity algorithm 
    281210      ! ------------------------- 
    282       CALL nonosc( tb, ztu, ztv, ztw, zti, z2 ) 
    283       CALL nonosc( sb, zsu, zsv, zsw, zsi, z2 ) 
     211      CALL nonosc( ptb, ztu, ztv, ztw, zti, z2 ) 
    284212 
    285213 
     
    294222               ztaj = - ( ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk  )) * zbtr 
    295223               ztak = - ( ztw(ji,jj,jk) - ztw(ji  ,jj  ,jk+1)) * zbtr 
    296                zsai = - ( zsu(ji,jj,jk) - zsu(ji-1,jj  ,jk  )) * zbtr 
    297                zsaj = - ( zsv(ji,jj,jk) - zsv(ji  ,jj-1,jk  )) * zbtr 
    298                zsak = - ( zsw(ji,jj,jk) - zsw(ji  ,jj  ,jk+1)) * zbtr 
    299224 
    300225               ! add them to the general tracer trends 
    301                ta(ji,jj,jk) = ta(ji,jj,jk) + ztai + ztaj + ztak 
    302                sa(ji,jj,jk) = sa(ji,jj,jk) + zsai + zsaj + zsak 
    303             END DO 
    304          END DO 
    305       END DO 
    306  
    307  
    308       ! Save the advective trends for diagnostics 
    309       ! -------------------------------------------- 
    310  
     226               pta(ji,jj,jk) = pta(ji,jj,jk) + ztai + ztaj + ztak 
     227            END DO 
     228         END DO 
     229      END DO 
     230 
     231!!gm  the transport computation is wrong, the upstream part is missing ! 
     232      ! "zonal" mean advective heat and salt transport 
     233      IF( cdtype == 'TRA' .AND. ln_diaptr .AND. ( MOD( kt, nf_ptr ) == 0 ) ) THEN 
     234         IF( ktra == jp_tem)   pht_adv(:) = ptr_vj( ztv(:,:,:) ) 
     235         IF( ktra == jp_sal)   pst_adv(:) = ptr_vj( ztv(:,:,:) ) 
     236      ENDIF 
     237 
     238      !  Save the horizontal advective trends for diagnostic 
     239      ! ----------------------------------------------------- 
    311240      IF( l_trdtra ) THEN 
    312          ztrdt(:,:,:) = 0.e0   ;   ztrds(:,:,:) = 0.e0 
    313          ! 
    314          ! T/S ZONAL advection trends 
    315          DO jk = 1, jpkm1 
    316             DO jj = 2, jpjm1 
    317                DO ji = fs_2, fs_jpim1   ! vector opt. 
    318                   !-- Compute zonal divergence by splitting hdivn (see divcur.F90) 
    319                   !   N.B. This computation is not valid along OBCs (if any) 
    320                   zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
    321                   z_hdivn_x = (  e2u(ji  ,jj) * fse3u(ji  ,jj,jk) * pun(ji  ,jj,jk)          & 
    322                      &         - e2u(ji-1,jj) * fse3u(ji-1,jj,jk) * pun(ji-1,jj,jk) ) * zbtr 
    323                   !-- Compute T/S zonal advection trends 
    324                   ztrdt(ji,jj,jk) = - ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) ) * zbtr + tn(ji,jj,jk) * z_hdivn_x 
    325                   ztrds(ji,jj,jk) = - ( zsu(ji,jj,jk) - zsu(ji-1,jj,jk) ) * zbtr + sn(ji,jj,jk) * z_hdivn_x 
    326                END DO 
    327             END DO 
    328          END DO 
    329          CALL trd_mod(ztrdt, ztrds, jptra_trd_xad, 'TRA', kt, cnbpas='bis')   ! <<< ADD TO PREVIOUSLY COMPUTED 
    330          ! 
    331          ! T/S MERIDIONAL advection trends 
    332          DO jk = 1, jpkm1 
    333             DO jj = 2, jpjm1 
    334                DO ji = fs_2, fs_jpim1   ! vector opt. 
    335                   !-- Compute merid. divergence by splitting hdivn (see divcur.F90) 
    336                   !   N.B. This computation is not valid along OBCs (if any) 
    337                   zbtr      = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
    338                   z_hdivn_y = (  e1v(ji,  jj) * fse3v(ji,jj  ,jk) * pvn(ji,jj  ,jk)          & 
    339                      &         - e1v(ji,jj-1) * fse3v(ji,jj-1,jk) * pvn(ji,jj-1,jk) ) * zbtr 
    340                   !-- Compute T/S meridional advection trends 
    341                   ztrdt(ji,jj,jk) = - ( ztv(ji,jj,jk) - ztv(ji,jj-1,jk) ) * zbtr + tn(ji,jj,jk) * z_hdivn_y           
    342                   ztrds(ji,jj,jk) = - ( zsv(ji,jj,jk) - zsv(ji,jj-1,jk) ) * zbtr + sn(ji,jj,jk) * z_hdivn_y           
    343                END DO 
    344             END DO 
    345          END DO 
    346          CALL trd_mod(ztrdt, ztrds, jptra_trd_yad, 'TRA', kt, cnbpas='bis')   ! <<< ADD TO PREVIOUSLY COMPUTED 
    347          ! 
    348          ! T/S VERTICAL advection trends 
    349          DO jk = 1, jpkm1 
    350             DO jj = 2, jpjm1 
    351                DO ji = fs_2, fs_jpim1   ! vector opt. 
    352                   zbtr1     = 1. / ( e1t(ji,jj) * e2t(ji,jj) ) 
    353 #if defined key_zco 
    354                   zbtr      = zbtr1 
    355                   z_hdivn_x = e2u(ji,jj)*pun(ji,jj,jk) - e2u(ji-1,jj)*pun(ji-1,jj,jk) 
    356                   z_hdivn_y = e1v(ji,jj)*pvn(ji,jj,jk) - e1v(ji,jj-1)*pvn(ji,jj-1,jk) 
    357 #else 
    358                   zbtr      = zbtr1 / fse3t(ji,jj,jk) 
    359                   z_hdivn_x = e2u(ji,jj)*fse3u(ji,jj,jk)*pun(ji,jj,jk) - e2u(ji-1,jj)*fse3u(ji-1,jj,jk)*pun(ji-1,jj,jk) 
    360                   z_hdivn_y = e1v(ji,jj)*fse3v(ji,jj,jk)*pvn(ji,jj,jk) - e1v(ji,jj-1)*fse3v(ji,jj-1,jk)*pvn(ji,jj-1,jk) 
    361 #endif 
    362                   z_hdivn   = (z_hdivn_x + z_hdivn_y) * zbtr 
    363                   zbtr      = zbtr1 / fse3t(ji,jj,jk) 
    364                   ztrdt(ji,jj,jk) = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * zbtr - tn(ji,jj,jk) * z_hdivn 
    365                   ztrds(ji,jj,jk) = - ( zsw(ji,jj,jk) - zsw(ji,jj,jk+1) ) * zbtr - sn(ji,jj,jk) * z_hdivn 
    366                END DO 
    367             END DO 
    368          END DO 
    369          CALL trd_mod(ztrdt, ztrds, jptra_trd_zad, 'TRA', kt, cnbpas='bis')   ! <<< ADD TO PREVIOUSLY COMPUTED 
    370          ! 
    371       ENDIF 
    372  
    373       IF(ln_ctl)   CALL prt_ctl( tab3d_1=ta, clinfo1=' tvd adv  - Ta: ', mask1=tmask,   & 
    374          &                       tab3d_2=sa, clinfo2=           ' Sa: ', mask2=tmask, clinfo3='tra' ) 
    375  
    376       ! "zonal" mean advective heat and salt transport 
    377       IF( ln_diaptr .AND. ( MOD( kt, nf_ptr ) == 0 ) ) THEN 
    378          pht_adv(:) = ptr_vj( ztv(:,:,:) ) 
    379          pst_adv(:) = ptr_vj( zsv(:,:,:) ) 
    380       ENDIF 
     241         CALL trd_tra_adv( kt, ktra, jpt_trd_xad, cdtype, ztu, pun, ptn, cnbpas='bis' )   ! <<< Add to iad trend 
     242         CALL trd_tra_adv( kt, ktra, jpt_trd_yad, cdtype, ztv, pvn, ptn, cnbpas='bis' )   ! <<< Add to jad trend 
     243         CALL trd_tra_adv( kt, ktra, jpt_trd_zad, cdtype, ztw, pwn, ptn, cnbpas='bis' )   ! <<< Add to zad trend 
     244      ENDIF 
     245 
     246      IF(ln_ctl)   CALL prt_ctl( tab3d_1=pta, clinfo1=' tvd - adv: ', mask1=tmask, clinfo3=cdtype ) 
    381247      ! 
    382248   END SUBROUTINE tra_adv_tvd 
     
    396262      !!       in-space based differencing for fluid 
    397263      !!---------------------------------------------------------------------- 
    398       REAL(wp), INTENT( in ) ::   prdt                               ! ??? 
    399       REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT( inout ) ::   & 
    400          pbef,                            & ! before field 
    401          paft,                            & ! after field 
    402          paa,                             & ! monotonic flux in the i direction 
    403          pbb,                             & ! monotonic flux in the j direction 
    404          pcc                                ! monotonic flux in the k direction 
     264      REAL(wp), INTENT(in   ) ::   prdt                               ! ??? 
     265      REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) ::   pbef, paft       ! before & after field 
     266      REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) ::   paa, pbb, pcc    ! monotonic flux in the 3 directions 
    405267      !! 
    406268      INTEGER ::   ji, jj, jk               ! dummy loop indices 
Note: See TracChangeset for help on using the changeset viewer.