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_qck.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_qck.F90

    r719 r786  
    11MODULE traadv_qck 
    2    !!============================================================================== 
     2   !!====================================================================== 
    33   !!                       ***  MODULE  traadv_qck  *** 
    4    !! Ocean active tracers:  horizontal & vertical advective trend 
    5    !!============================================================================== 
    6  
     4   !! Ocean tracers:  horizontal & vertical advective trend using QUICKEST scheme 
     5   !!====================================================================== 
     6   !! History :  1.0  !  06-09  (G. Reffray)  Original code 
     7   !!            2.4  !  08-01  (G. Madec)  Merge TRA-TRC 
    78   !!---------------------------------------------------------------------- 
    8    !!   tra_adv_qck : update the tracer trend with the horizontal 
    9    !!                      advection trends using a 3st order 
    10    !!                      finite difference scheme 
    11    !!                      The vertical advection scheme is the 2nd centered scheme 
     9 
    1210   !!---------------------------------------------------------------------- 
    13    !! * Modules used 
    14    USE oce             ! ocean dynamics and active tracers 
     11   !!   tra_adv_qck : update the tracer trend with the horizontal advection 
     12   !!                 trends using a 3rd order finite difference scheme        ????? 
     13   !!                 The vertical advection scheme is the 2nd centered scheme ????? 
     14   !!---------------------------------------------------------------------- 
    1515   USE dom_oce         ! ocean space and time domain 
    1616   USE trdmod          ! ocean active tracers trends  
    1717   USE trdmod_oce      ! ocean variables trends 
    1818   USE flxrnf          ! 
    19    USE trabbl          ! advective term in the BBL 
    2019   USE ocfzpt          ! 
    2120   USE lib_mpp 
     
    2928   PRIVATE 
    3029 
    31    !! * Accessibility 
    32    PUBLIC tra_adv_qck    ! routine called by step.F90 
    33  
    34    !! * Module variables 
    35    REAL(wp), DIMENSION(jpi,jpj),     SAVE ::   & 
    36          zbtr2 
    37    REAL(wp), DIMENSION(jpi,jpj,jpk), SAVE ::   & 
    38          sl 
    39    REAL(wp) ::                                 & 
    40          cst1, cst2, dt, coef1                  ! temporary scalars 
    41    INTEGER  ::                                 &  
    42          ji, jj, jk                             ! dummy loop indices 
    43    !!---------------------------------------------------------------------- 
     30   PUBLIC tra_adv_qck    ! routine called by traadv.F90 
     31 
     32   REAL(wp), DIMENSION(jpi,jpj)     ::   zbtr2 
     33   REAL(wp), DIMENSION(jpi,jpj,jpk) ::   sl 
     34   REAL(wp) ::   cst1, cst2, dt, coef1                  ! temporary scalars 
     35   INTEGER  ::   ji, jj, jk                             ! dummy loop indices 
     36 
    4437   !! * Substitutions 
    4538#  include "domzgr_substitute.h90" 
    4639#  include "vectopt_loop_substitute.h90" 
    4740   !!---------------------------------------------------------------------- 
    48    !!   OPA 9.0 , LOCEAN-IPSL (2005)  
    49    !! $Header$  
    50    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     41   !! NEMO/OPA 2.4 , LOCEAN-IPSL (2008)  
     42   !! $Id:$  
     43   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    5144   !!---------------------------------------------------------------------- 
    5245 
    5346CONTAINS 
    5447 
    55 #if ! defined key_mpp_omp 
    56    !!---------------------------------------------------------------------- 
    57    !!   Default option :             quickest advection scheme (k-j-i loop) 
    58    !!---------------------------------------------------------------------- 
    59  
    60    SUBROUTINE tra_adv_qck( kt, pun, pvn, pwn ) 
     48   SUBROUTINE tra_adv_qck( kt, cdtype, ktra, pun, pvn, pwn,   & 
     49      &                                      ptb, ptn, pta ) 
    6150      !!---------------------------------------------------------------------- 
    6251      !!                  ***  ROUTINE tra_adv_qck  *** 
     
    6756      !! ** Method :   The advection is evaluated by a third order scheme 
    6857      !!               For a positive velocity u : 
    69       !! 
    7058      !! 
    7159      !!                  i-1    i      i+1   i+2 
     
    8371      !!                FC is the central point (or the first upwind point) 
    8472      !! 
    85       !!      Flux(i) = u(i) * {0.5(FC+FD)  -0.5C(i)(FD-FC)  -((1-C(i)å?)/6)(FU+FD-2FC)} 
     73      !!      Flux(i) = u(i) * {0.5(FC+FD)  -0.5C(i)(FD-FC)  -((1-C(i)Â?)/6)(FU+FD-2FC)} 
    8674      !!                with C(i)=|u(i)|dx(i)/dt (Courant number) 
    8775      !! 
     
    9785      !!             - save the trends in (ttrdh,strdh) ('key_trdtra') 
    9886      !! 
    99       !! ** Reference : Leonard (1979, 1991) 
    100       !! History : 
    101       !!  9.0     !  06-09  (G. Reffray)  Original code 
    102       !!---------------------------------------------------------------------- 
    103       !! * Arguments 
    104       INTEGER, INTENT( in ) ::   kt             ! ocean time-step index 
    105       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::  pun   ! effective ocean velocity, u_component 
    106       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::  pvn   ! effective ocean velocity, v_component 
    107       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::  pwn   ! effective ocean velocity, w_component 
    108       !! 
    109       REAL(wp) :: z2                                          ! temporary scalar  
    110       REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztrdt, ztrds       ! temporary 3D workspace 
     87      !! References : Leonard (1979, 1991) 
     88      !!---------------------------------------------------------------------- 
     89      INTEGER         , INTENT(in   )                         ::   kt              ! ocean time-step index 
     90      CHARACTER(len=3), INTENT(in   )                         ::   cdtype          ! =TRA or TRC (tracer indicator) 
     91      INTEGER         , INTENT(in   )                         ::   ktra            ! tracer index 
     92      REAL(wp)        , INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   pun, pvn, pwn   ! 3 ocean velocity components 
     93      REAL(wp)        , INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   ptb, ptn        ! before and now tracer fields 
     94      REAL(wp)        , INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pta             ! tracer trend  
     95      !! 
     96      REAL(wp) :: z2                                   ! temporary scalar  
     97      !!---------------------------------------------------------------------- 
    11198 
    11299      IF( kt == nit000 ) THEN 
    113100         IF(lwp) WRITE(numout,*) 
    114101         IF(lwp) WRITE(numout,*) 'tra_adv_qck : 3st order quickest advection scheme' 
    115          IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~   Vector optimization case' 
    116          IF(lwp) WRITE(numout,*) 
     102         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 
    117103 
    118104         zbtr2(:,:) = 1. / ( e1t(:,:) * e2t(:,:) ) 
    119105         cst1 = 1./12. 
    120106         cst2 = 2./3. 
    121          IF (l_trdtra ) THEN 
    122          CALL ctl_warn( ' Trends not yet implemented for PPM advection scheme ' ) 
    123          ENDIF 
    124107      ENDIF 
    125108 
     
    128111      ENDIF 
    129112 
    130       ! Save ta and sa trends 
    131       IF( l_trdtra )   THEN      ! to be done 
    132          ztrdt(:,:,:) = ta(:,:,:)  
    133          ztrds(:,:,:) = sa(:,:,:)  
    134          l_adv = 'qst' 
    135       ENDIF 
    136113 
    137114      ! I. Slope estimation at the T-point for the limiter ULTIMATE 
     
    163140      CALL lbc_lnk(   sl(:,:,:), 'T', 1. ) 
    164141 
     142 
    165143      ! II. The horizontal fluxes are computed with the QUICKEST + ULTIMATE scheme 
    166144      !--------------------------------------------------------------------------- 
    167145 
    168       CALL tra_adv_qck_hor( kt , pun, pvn, tb , ta , pht_adv , z2) 
    169       CALL tra_adv_qck_hor( kt , pun, pvn, sb , sa , pst_adv , z2)   
    170  
    171       ! Save the horizontal advective trends for diagnostic 
    172       ! --------------------------------------------------- 
    173 !      IF( l_trdtra )   THEN    ! to be done 
    174 !         ! T/S ZONAL advection trends 
    175 !      ENDIF 
    176  
    177       IF(ln_ctl)   THEN 
    178           CALL prt_ctl(tab3d_1=ta, clinfo1=' centered2 had  - Ta: ', mask1=tmask, & 
    179              &         tab3d_2=sa, clinfo2=' Sa: ', mask2=tmask, clinfo3='tra') 
    180       ENDIF 
     146      CALL tra_adv_qck_hor( kt, cdtype, ktra, pun, pvn, ptb, pta, z2 ) 
     147 
     148      IF(ln_ctl)   CALL prt_ctl( tab3d_1=pta, clinfo1=' qck - had: ', mask1=tmask, clinfo3=cdtype ) 
     149 
    181150 
    182151      ! III. The vertical fluxes are computed with the 2nd order centered scheme 
    183152      !------------------------------------------------------------------------- 
    184153 
    185       CALL tra_adv_qck_ver( pwn, tn , ta, z2 ) 
    186       CALL tra_adv_qck_ver( pwn, sn , sa, z2 ) 
    187  
    188       ! Save the vertical advective trends for diagnostic 
    189       ! ------------------------------------------------- 
    190 !      IF( l_trdtra )   THEN    ! to be done 
    191          ! Recompute the vertical advection zta & zsa trends computed 
    192          ! at the step 2. above in making the difference between the new 
    193          ! trends and the previous one: ta()/sa - ztdta()/ztdsa() and substract 
    194          ! the term tn()/sn()*hdivn() to recover the W gradz(T/S) trends 
    195 !     ENDIF 
    196  
    197       IF(ln_ctl)   THEN 
    198           CALL prt_ctl(tab3d_1=ta, clinfo1=' centered2 zad  - Ta: ', mask1=tmask, & 
    199              &         tab3d_2=sa, clinfo2=' Sa: ', mask2=tmask, clinfo3='tra') 
    200       ENDIF 
    201  
     154      CALL tra_adv_qck_ver( pwn, ptn , pta, z2 ) 
     155 
     156      IF(ln_ctl)   CALL prt_ctl( tab3d_1=pta, clinfo1=' qck - zad: ', mask1=tmask, clinfo3=cdtype ) 
     157      ! 
    202158   END SUBROUTINE tra_adv_qck 
    203159 
    204    SUBROUTINE tra_adv_qck_hor ( kt , pun, pvn, tra , traa , phtra_adv ,z2 ) 
    205       !!---------------------------------------------------------------------- 
    206       !! 
    207       !!---------------------------------------------------------------------- 
    208       !! * Arguments 
    209       INTEGER, INTENT( in ) ::   kt             ! ocean time-step index 
    210       REAL, INTENT( in )    ::   z2 
    211       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pun, pvn            ! horizontal effective velocity 
    212  
    213       REAL(wp), INTENT ( out   ), DIMENSION(jpj)          ::     & 
    214          phtra_adv 
    215  
    216       REAL(wp), INTENT ( inout ), DIMENSION(jpi,jpj,jpk)  ::     & 
    217          tra, traa 
    218  
    219       REAL(wp) ::                                & 
    220          za, zbtr, e1, e2, c, dir, fu, fc, fd,   & ! temporary scalars 
    221          coef2, coef3, fho, mask, dx 
    222  
    223       REAL(wp), DIMENSION(jpi,jpj) ::            & 
    224          zee 
    225  
    226       REAL(wp), DIMENSION(jpi,jpj,jpk) ::  & 
    227          zmask, zlap, dwst, lim  
    228  
    229  
     160 
     161   SUBROUTINE tra_adv_qck_hor ( kt, cdtype, ktra, pun, pvn, ptb, pta, z2 ) 
     162      !!---------------------------------------------------------------------- 
     163      !! 
     164      !!---------------------------------------------------------------------- 
     165      INTEGER         , INTENT(in   )                         ::   kt          ! ocean time-step index 
     166      CHARACTER(len=3), INTENT(in   )                         ::   cdtype      ! =TRA or TRC (tracer indicator) 
     167      INTEGER         , INTENT(in   )                         ::   ktra        ! tracer index 
     168      REAL(wp)        , INTENT(in   )                         ::   z2          ! ??? 
     169      REAL(wp)        , INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   pun, pvn    ! horizontal effective velocity 
     170      REAL(wp)        , INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   ptb         ! before tracer field 
     171      REAL(wp)        , INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pta         ! tracer trend 
     172 
     173      REAL(wp) ::   za, zbtr, e1, e2, c, dir, fu, fc, fd   ! temporary scalars 
     174      REAL(wp) ::   coef2, coef3, fho, mask, dx 
     175      REAL(wp), DIMENSION(jpi,jpj) ::  zee 
     176      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zmask, zlap, dwst, lim  
     177      !!---------------------------------------------------------------------- 
    230178 
    231179      !---------------------------------------------------------------------- 
    232180      ! 0. Initialization (should ot be needed on the whole array ???) 
    233181      !---------------------------------------------------------------------- 
    234  
    235       zmask = 0.0                                                 
    236       zlap  = 0.0                                                   
    237       dwst  = 0.0                                                   
    238       lim   = 0.0      
     182      zmask(:,:,:)= 0.e0                                                 
     183      zlap (:,:,:)= 0.e0                                                   
     184      dwst (:,:,:)= 0.e0                                                   
     185      lim  (:,:,:)= 0.e0      
    239186                                             
    240187      !---------------------------------------------------------------------- 
     
    264211         DO jj = 1, jpjm1 
    265212            DO ji = 1, fs_jpim1   ! vector opt. 
    266                zmask(ji,jj,jk) = zee(ji,jj) * ( tra(ji+1,jj  ,jk) - tra(ji,jj,jk) ) 
     213               zmask(ji,jj,jk) = zee(ji,jj) * ( ptb(ji+1,jj  ,jk) - ptb(ji,jj,jk) ) 
    267214            END DO 
    268215         END DO 
     
    279226         !--- Function lim=FU+SL*(FC-FU) used by the limiter 
    280227         !--- Computation of the ustream and downstream lim at the T-points 
     228!!gm bug : fs_2  instead of 2 ... 
     229!!gm  a lot of optimisation to be done in this routine.... 
    281230         DO jj = 2, jpjm1 
    282231            DO ji = 2, fs_jpim1   ! vector opt. 
    283232               ! Upstream in the x-direction for the tracer 
    284                zmask(ji,jj,jk)=tra(ji-1,jj,jk)+sl(ji,jj,jk)*(tra(ji,jj,jk)-tra(ji-1,jj,jk)) 
     233               zmask(ji,jj,jk) =ptb(ji-1,jj,jk)+sl(ji,jj,jk)*(ptb(ji,jj,jk)-ptb(ji-1,jj,jk)) 
    285234               ! Downstream in the x-direction for the tracer 
    286                dwst (ji,jj,jk)=tra(ji+1,jj,jk)+sl(ji,jj,jk)*(tra(ji,jj,jk)-tra(ji+1,jj,jk)) 
     235               dwst (ji,jj,jk) =ptb(ji+1,jj,jk)+sl(ji,jj,jk)*(ptb(ji,jj,jk)-ptb(ji+1,jj,jk)) 
    287236            ENDDO 
    288237         ENDDO 
     
    329278 
    330279               fu  = lim(ji,jj,jk)                                ! FU + sl(FC-FU) in the x-direction for T 
    331                fc  = dir*tra(ji  ,jj,jk)+(1-dir)*tra(ji+1,jj,jk)  ! FC in the x-direction for T 
    332                fd  = dir*tra(ji+1,jj,jk)+(1-dir)*tra(ji  ,jj,jk)  ! FD in the x-direction for T 
     280               fc  = dir*ptb(ji  ,jj,jk)+(1-dir)*ptb(ji+1,jj,jk)  ! FC in the x-direction for T 
     281               fd  = dir*ptb(ji+1,jj,jk)+(1-dir)*ptb(ji  ,jj,jk)  ! FD in the x-direction for T 
    333282 
    334283               !--- QUICKEST scheme 
     
    358307               za = - zbtr * ( dwst(ji,jj,jk) - dwst(ji-1,jj  ,jk) ) 
    359308               !--- add it to the general tracer trends 
    360                traa(ji,jj,jk) = traa(ji,jj,jk) + za 
     309               pta(ji,jj,jk) = pta(ji,jj,jk) + za 
    361310            END DO 
    362311         END DO 
     
    364313      END DO                                           !   End of slab 
    365314      !                                                ! =============== 
     315 
     316      ! Save the horizontal advective trends for diagnostic 
     317      IF( l_trdtra )   CALL trd_tra_adv( kt, ktra, jpt_trd_xad, cdtype, dwst, pun, ptb ) 
     318 
     319 
    366320      !---------------------------------------------------------------------- 
    367321      ! I. Part 2 : y-direction 
     
    389343         DO jj = 1, jpjm1 
    390344            DO ji = 1, fs_jpim1   ! vector opt. 
    391                zmask(ji,jj,jk) = zee(ji,jj) * ( tra(ji  ,jj+1,jk) - tra(ji,jj,jk) ) 
     345               zmask(ji,jj,jk) = zee(ji,jj) * ( ptb(ji  ,jj+1,jk) - ptb(ji,jj,jk) ) 
    392346            END DO 
    393347         END DO 
     
    408362            DO ji = 2, fs_jpim1   ! vector opt. 
    409363               ! Upstream in the y-direction for the tracer 
    410                zmask(ji,jj,jk)=tra(ji,jj-1,jk)+sl(ji,jj,jk)*(tra(ji,jj,jk)-tra(ji,jj-1,jk)) 
     364               zmask(ji,jj,jk)=ptb(ji,jj-1,jk)+sl(ji,jj,jk)*(ptb(ji,jj,jk)-ptb(ji,jj-1,jk)) 
    411365               ! Downstream in the y-direction for the tracer 
    412                dwst (ji,jj,jk)=tra(ji,jj+1,jk)+sl(ji,jj,jk)*(tra(ji,jj,jk)-tra(ji,jj+1,jk)) 
     366               dwst (ji,jj,jk)=ptb(ji,jj+1,jk)+sl(ji,jj,jk)*(ptb(ji,jj,jk)-ptb(ji,jj+1,jk)) 
    413367            ENDDO 
    414368         ENDDO 
     
    455409 
    456410               fu  = lim(ji,jj,jk)                                ! FU + sl(FC-FU) in the y-direction for T 
    457                fc  = dir*tra(ji,jj  ,jk)+(1-dir)*tra(ji,jj+1,jk)  ! FC in the y-direction for T 
    458                fd  = dir*tra(ji,jj+1,jk)+(1-dir)*tra(ji,jj  ,jk)  ! FD in the y-direction for T 
     411               fc  = dir*ptb(ji,jj  ,jk)+(1-dir)*ptb(ji,jj+1,jk)  ! FC in the y-direction for T 
     412               fd  = dir*ptb(ji,jj+1,jk)+(1-dir)*ptb(ji,jj  ,jk)  ! FD in the y-direction for T 
    459413 
    460414               !--- QUICKEST scheme 
     
    484438               za = - zbtr * (  dwst(ji,jj,jk) - dwst(ji  ,jj-1,jk) ) 
    485439               !--- add it to the general tracer trends 
    486                traa(ji,jj,jk) = traa(ji,jj,jk) + za 
     440               pta(ji,jj,jk) = pta(ji,jj,jk) + za 
    487441            END DO 
    488442         END DO 
     
    491445      !                                                ! =============== 
    492446 
     447      ! Save the horizontal advective trends for diagnostic 
     448      IF( l_trdtra )   CALL trd_tra_adv( kt, ktra, jpt_trd_yad, cdtype, dwst, pvn, ptb ) 
     449 
    493450      ! "zonal" mean advective heat and salt transport  
    494       IF( ln_diaptr .AND. ( MOD( kt, nf_ptr ) == 0 ) ) THEN 
     451      IF( ln_diaptr .AND. cdtype == 'TRA' .AND. ( MOD( kt, nf_ptr ) == 0 ) ) THEN 
    495452#if defined key_zco 
    496453         DO jk = 1, jpkm1 
     
    501458            END DO 
    502459         END DO 
    503          phtra_adv(:) = ptr_vj( dwst(:,:,:) ) 
    504 #else 
    505          phtra_adv(:) = ptr_vj( dwst(:,:,:) ) 
    506460# endif 
     461         IF( ktra == jp_tem)   pht_adv(:) = ptr_vj( dwst(:,:,:) ) 
     462         IF( ktra == jp_sal)   pst_adv(:) = ptr_vj( dwst(:,:,:) ) 
    507463      ENDIF 
    508  
     464      ! 
    509465   END SUBROUTINE tra_adv_qck_hor 
    510466 
    511    SUBROUTINE tra_adv_qck_ver ( pwn, tra , traa, z2  )       
    512       !!---------------------------------------------------------------------- 
    513       !! 
    514       !!---------------------------------------------------------------------- 
    515       !! * Arguments 
    516  
    517       REAL(wp), INTENT ( in ) :: z2 
    518       REAL(wp), INTENT ( in ), DIMENSION(jpi,jpj,jpk)  ::     & 
    519          pwn 
    520       REAL(wp), INTENT ( inout ), DIMENSION(jpi,jpj,jpk)  ::     & 
    521          tra, traa 
    522  
    523       REAL(wp) ::                             & 
    524          za, ze3tr, dt, dir, fc, fd             ! temporary scalars 
     467 
     468   SUBROUTINE tra_adv_qck_ver( pwn, ptn , pta, z2  )       
     469      !!---------------------------------------------------------------------- 
     470      !! 
     471      !!---------------------------------------------------------------------- 
     472      REAL(wp), INTENT(in   )                         ::   z2 
     473      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   pwn 
     474      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   ptn 
     475      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pta 
     476      !! 
     477      REAL(wp) ::   za, ze3tr, dt, dir, fc, fd             ! temporary scalars 
    525478 
    526479      ! Vertical advection 
     
    530483      ! ---------------------------- 
    531484 
    532       !Bottom value : flux set to zero 
    533       sl(:,:,jpk) = 0.e0 
     485      sl(:,:,jpk) = 0.e0      !Bottom value : flux set to zero 
    534486 
    535487      ! Surface value 
    536       IF( lk_dynspg_rl .OR. lk_vvl ) THEN 
    537          ! rigid lid : flux set to zero 
     488      IF( lk_dynspg_rl .OR. lk_vvl ) THEN         ! rigid lid : flux set to zero 
    538489         sl(:,:, 1 ) = 0.e0 
    539       ELSE 
    540          ! free surface-constant volume 
    541          sl(:,:, 1 ) = pwn(:,:,1) * tra(:,:,1) 
     490      ELSE                                        ! free surface-constant volume 
     491         sl(:,:, 1 ) = pwn(:,:,1) * ptn(:,:,1) 
    542492      ENDIF 
    543493 
    544494      ! Second order centered tracer flux at w-point 
    545  
    546495      DO jk = 2, jpkm1 
    547496         dt  = z2 *  rdttra(jk) 
     
    549498            DO ji = fs_2, fs_jpim1   ! vector opt. 
    550499               dir = 0.5 + sign(0.5,pwn(ji,jj,jk))                                         ! if pwn>0 : dirw = 1 otherwise dirw = 0 
    551                fc = dir*tra(ji,jj,jk  )*fse3t(ji,jj,jk-1)+(1-dir)*tra(ji,jj,jk-1)*fse3t(ji,jj,jk  )   ! FC in the z-direction for T 
    552                fd = dir*tra(ji,jj,jk-1)*fse3t(ji,jj,jk  )+(1-dir)*tra(ji,jj,jk  )*fse3t(ji,jj,jk-1)   ! FD in the z-direction for T 
     500               fc = dir*ptn(ji,jj,jk  )*fse3t(ji,jj,jk-1)+(1-dir)*ptn(ji,jj,jk-1)*fse3t(ji,jj,jk  )   ! FC in the z-direction for T 
     501               fd = dir*ptn(ji,jj,jk-1)*fse3t(ji,jj,jk  )+(1-dir)*ptn(ji,jj,jk  )*fse3t(ji,jj,jk-1)   ! FD in the z-direction for T 
    553502               !--- Second order centered scheme 
    554503               sl(ji,jj,jk)=pwn(ji,jj,jk)*(fc+fd)/(fse3t(ji,jj,jk-1)+fse3t(ji,jj,jk)) 
     
    559508      ! 2. Tracer flux divergence at t-point added to the general trend 
    560509      ! --------------------------------------------------------------- 
    561  
    562510      DO jk = 1, jpkm1 
    563511         DO jj = 2, jpjm1 
     
    567515               za = - ze3tr * ( sl(ji,jj,jk) - sl(ji,jj,jk+1) ) 
    568516               ! add it to the general tracer trends 
    569                traa(ji,jj,jk) =  traa(ji,jj,jk) + za 
     517               pta(ji,jj,jk) =  pta(ji,jj,jk) + za 
    570518            END DO 
    571519         END DO 
    572520      END DO 
    573  
     521      ! 
    574522   END SUBROUTINE tra_adv_qck_ver 
    575523 
     524 
    576525   REAL FUNCTION bound(fu,fd,fc,fho) 
    577       real     ::  fu,fd,fc,fho,fref1,fref2 
     526      REAL(wp) ::  fu, fd, fc, fho, fref1, fref2 
    578527      fref1 = fu 
    579       fref2 = MAX(MIN(fc,fd),MIN(MAX(fc,fd),fref1)) 
    580       bound = MAX(MIN(fho,fc),MIN(MAX(fho,fc),fref2)) 
     528      fref2 = MAX(  MIN( fc , fd ), MIN( MAX( fc , fd ), fref1 )  ) 
     529      bound = MAX(  MIN( fho, fc ), MIN( MAX( fho, fc ), fref2 )  ) 
    581530   END FUNCTION 
    582  
    583 #else 
    584    !!---------------------------------------------------------------------- 
    585    !!   'key_mpp_omp' :      quickest advection (k- and j-slabs) 
    586    !!---------------------------------------------------------------------- 
    587    SUBROUTINE tra_adv_qck( kt, pun, pvn, pwn  ) 
    588       !!---------------------------------------------------------------------- 
    589       !! * Arguments 
    590       INTEGER, INTENT( in ) ::   kt             ! ocean time-step index 
    591       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::  pun   ! effective ocean velocity, u_component 
    592       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::  pvn   ! effective ocean velocity, v_component 
    593       REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) ::  pwn   ! effective ocean velocity, w_component 
    594       !!----------------------------------------------------------------------   
    595          IF(lwp) WRITE(numout,*) 
    596          IF(lwp) WRITE(numout,*) 'tra_adv_ qck:3st order quickest advection scheme'  
    597          IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~   Vector optimization case' 
    598          IF(lwp) WRITE(numout,*) 'WITH AUTOTASKING =>this routine doesn t exist for the moment' 
    599     IF(lwp) WRITE(numout,*) ' EMPTY ROUTINE!!!!!!'       
    600  
    601    END SUBROUTINE tra_adv_qck 
    602  
    603 #endif 
    604531 
    605532   !!====================================================================== 
Note: See TracChangeset for help on using the changeset viewer.