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 791 for branches/dev_001_GM/NEMO/OPA_SRC/TRA/traadv_eiv.F90 – NEMO

Ignore:
Timestamp:
2008-01-12T21:33:34+01:00 (16 years ago)
Author:
gm
Message:

dev_001_GM - TRA/traadv : switch from velocity to transport + optimised traadv_eiv2 introduced - compilation OK

File:
1 edited

Legend:

Unmodified
Added
Removed

  • branches/dev_001_GM/NEMO/OPA_SRC/TRA/traadv_eiv.F90

    r719 r791  
    55   !!====================================================================== 
    66   !! History :  9.0  !  05-11  (G. Madec)  Original code, from traldf and zdf _iso 
     7   !!            2.4  !  2008-01  (G. Madec)  merge TRC-TRA + switch from velocity to transport 
    78   !!---------------------------------------------------------------------- 
    89#if defined key_traldf_eiv   ||   defined key_esopa 
    … …  
    3233#  include "vectopt_loop_substitute.h90" 
    3334   !!---------------------------------------------------------------------- 
    34    !!  OPA 9.0 , LOCEAN-IPSL (2006)  
    35    !! $Header$  
     35   !! NEMO/OPA 2.4 , LOCEAN-IPSL (2008)  
     36   !! $Id:$  
    3637   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    3738   !!---------------------------------------------------------------------- 
    … …  
    4344      !!                  ***  ROUTINE tra_adv_eiv  *** 
    4445      !!  
    45       !! ** Purpose :   Compute the before horizontal tracer (t & s) diffusive  
    46       !!      trend and add it to the general trend of tracer equation. 
    47       !! 
    48       !! ** Method  :   The eddy induced advection is computed from the slope 
     46      !! ** Purpose :   Compute the eddy induced transport and add it to the 
     47      !!              effective transport 
     48      !! 
     49      !! ** Method  :   The eddy induced transport is computed from the slope 
    4950      !!      of iso-neutral surfaces computed in routine ldf_slp as follows: 
    50       !!         zu_eiv =  1/(e2u e3u)   dk[ aeiu e2u mi(wslpi) ] 
    51       !!         zv_eiv =  1/(e1v e3v)   dk[ aeiv e1v mj(wslpj) 
    52       !!         zw_eiv = -1/(e1t e2t) { di[ aeiu e2u mi(wslpi) ] 
    53       !!                               + dj[ aeiv e1v mj(wslpj) ] } 
     51      !!         zu_eiv =   dk[ aeiu e2u mi(wslpi) ] 
     52      !!         zv_eiv =   dk[ aeiv e1v mj(wslpj) ] 
     53      !!         zw_eiv = - { di[ aeiu e2u mi(wslpi) ] + dj[ aeiv e1v mj(wslpj) ] } 
    5454      !!      add the eiv component to the model velocity: 
    5555      !!         p.n = p.n + z._eiv 
     56      !! CAUTION : the horizontal transports not updated along jpi column and jpj row 
     57      !!           the vertical   transports not updated along 1 & jpi columns and 1 & jpj rows 
    5658      !! 
    5759      !! ** Action  : - add to p.n the eiv component 
    5860      !!---------------------------------------------------------------------- 
    5961      INTEGER , INTENT(in   )                         ::   kt    ! ocean time-step index 
    60       REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pun   ! in : 3 ocean velocity components  
    61       REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pvn   ! out: 3 ocean velocity components 
     62      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pun   ! in : 3 ocean transport components  
     63      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pvn   ! out: 3 ocean transport components 
    6264      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pwn   !      increased by the eiv 
    6365      !! 
    … …  
    8890               zvwk1= ( wslpj(ji,jj,jk+1) + wslpj(ji,jj+1,jk+1) ) * fsaeiv(ji,jj,jk+1) * vmask(ji,jj,jk+1) 
    8991 
    90                zu_eiv = 0.5 * umask(ji,jj,jk) * ( zuwk - zuwk1 ) / fse3u(ji,jj,jk) 
    91                zv_eiv = 0.5 * vmask(ji,jj,jk) * ( zvwk - zvwk1 ) / fse3v(ji,jj,jk) 
     92               zu_eiv = 0.5 * umask(ji,jj,jk) * ( zuwk - zuwk1 ) 
     93               zv_eiv = 0.5 * vmask(ji,jj,jk) * ( zvwk - zvwk1 ) 
    9294    
    93                pun(ji,jj,jk) = pun(ji,jj,jk) + zu_eiv 
    94                pvn(ji,jj,jk) = pvn(ji,jj,jk) + zv_eiv 
    95 # if defined key_diaeiv 
    96                u_eiv(ji,jj,jk) = zu_eiv 
    97                v_eiv(ji,jj,jk) = zv_eiv 
     95               pun(ji,jj,jk) = pun(ji,jj,jk) + e2u(ji,jj) * zu_eiv 
     96               pvn(ji,jj,jk) = pvn(ji,jj,jk) + e1v(ji,jj) * zv_eiv 
     97# if defined key_diaeiv 
     98               u_eiv(ji,jj,jk) = zu_eiv / fse3u(ji,jj,jk) 
     99               v_eiv(ji,jj,jk) = zv_eiv / fse3v(ji,jj,jk) 
    98100# endif 
    99101            END DO 
    … …  
    115117                  zvwj1 = ( wslpj(ji,jj,jk) + wslpj(ji,jj+1,jk) ) * e1v(ji  ,jj) * vmask(ji  ,jj,jk) 
    116118 
    117                   zw_eiv = - 0.5 * tmask(ji,jj,jk) * fsaeiw(ji,jj,jk) * ( zuwi1 - zuwi + zvwj1 - zvwj )  & 
    118                      &                                                / ( e1t(ji,jj)*e2t(ji,jj) ) 
     119                  zw_eiv = - 0.5 * tmask(ji,jj,jk) * fsaeiw(ji,jj,jk) * ( zuwi1 - zuwi + zvwj1 - zvwj ) 
    119120# endif 
    120121                  pwn(ji,jj,jk) = pwn(ji,jj,jk) + zw_eiv 
    121122 
    122123# if defined key_diaeiv 
    123                   w_eiv(ji,jj,jk) = zw_eiv 
     124                  w_eiv(ji,jj,jk) = zw_eiv / ( e1t(ji,jj)*e2t(ji,jj) ) 
    124125# endif 
    125126               END DO 
    … …  
    130131      !                                             ! ================= 
    131132   END SUBROUTINE tra_adv_eiv 
     133 
     134 
     135!!gm   test tra_adv_eiv better (faster) coded?  to be verified 
     136 
     137   SUBROUTINE tra_adv_eiv2( kt, pun, pvn, pwn ) 
     138      !!---------------------------------------------------------------------- 
     139      !!                  ***  ROUTINE tra_adv_eiv  *** 
     140      !!  
     141      !! ** Purpose :   Compute the eddy induced transport and add it to the 
     142      !!              effective transport 
     143      !! 
     144      !! ** Method  :   The eddy induced transport is computed from the slope 
     145      !!              of iso-neutral surfaces (see ldfslp.F90) as follows: 
     146      !!                   zu_eiv =   dk[ aeiu e2u mi(wslpi) ] 
     147      !!                   zv_eiv =   dk[ aeiv e1v mj(wslpj) ] 
     148      !!                   zw_eiv = - { di[ aeiu e2u mi(wslpi) ] + dj[ aeiv e1v mj(wslpj) ] } 
     149      !!              add the eiv component to the model velocity: 
     150      !!                   p.n = p.n + z._eiv 
     151      !! CAUTION : the horizontal transports not updated along jpi column and jpj row 
     152      !!           the vertical   transports not updated along 1 & jpi columns and 1 & jpj rows 
     153      !! 
     154      !! ** Action  : - add to p.n the eiv transport component 
     155      !!---------------------------------------------------------------------- 
     156      INTEGER , INTENT(in   )                         ::   kt    ! ocean time-step index 
     157      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pun   ! in : 3 ocean transport components  
     158      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pvn   ! out: 3 ocean transport components 
     159      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) ::   pwn   !      increased by the eiv 
     160      !! 
     161      INTEGER  ::   ji, jj, jk                 ! dummy loop indices 
     162      REAL(wp) ::   zuwslpi, zvwslpj           ! temporary scalar 
     163      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zsfu, zsfv   ! eiv stream-function in u and v directions 
     164      !!---------------------------------------------------------------------- 
     165 
     166      IF( kt == nit000 ) THEN 
     167         IF(lwp) WRITE(numout,*) 
     168         IF(lwp) WRITE(numout,*) 'tra_adv_eiv : eddy induced advection :' 
     169         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~   add to velocity fields the eiv component' 
     170      ENDIF 
     171 
     172      ! eiv stream-function in u- and v-directions 
     173      ! NB: UW-point mask at level k is umask(:,:,k) idem form VW-point mask 
     174      zsfu(:,:, 1 ) = 0.e0   ;   zsfv(:,:, 1 ) = 0.e0      ! surface value set to zero 
     175      zsfu(:,:,jpk) = 0.e0   ;   zsfv(:,:,jpk) = 0.e0      ! bottom  value set to zero 
     176      DO jk = 2, jpkm1 
     177         DO jj = 1, jpjm1 
     178            DO ji = 1, fs_jpim1   ! vector opt. 
     179               zuwslpi = 0.5 * ( wslpi(ji+1,jj,jk) + wslpi(ji,jj,jk) ) 
     180               zvwslpj = 0.5 * ( wslpj(ji,jj+1,jk) + wslpj(ji,jj,jk) )  
     181               zsfu(ji,jj,jk) = zuwslpi * e2u(ji,jj) * fsaeiu(ji,jj,jk) * umask(ji,jj,jk) 
     182               zsfv(ji,jj,jk) = zvwslpj * e1v(ji,jj) * fsaeiv(ji,jj,jk) * vmask(ji,jj,jk) 
     183            END DO 
     184         END DO 
     185      END DO 
     186# if defined key_diaeiv 
     187      ! save eiv stream function in the output 
     188!!gm  to be done, u_sfeiv and v_sfeiv not defined    ==> new IOM.... 
     189!!gm  and zsfu, zfv notdefined for jpi column and jpj row 
     190!     u_sfeiv(:,:,:) = zsfu(:,:,:) 
     191!     v_sfeiv(:,:,:) = zsfu(:,:,:) 
     192# endif 
     193 
     194      ! increase the transport with the eiv transport 
     195      DO jk = 1, jpkm1 
     196         DO jj = 1, jpjm1 
     197            DO ji = 1, fs_jpim1   ! vector opt. 
     198               pun(ji,jj,jk) = pun(ji,jj,jk) + ( zsfu(ji,jj,jk) - zsfu(ji,jj,jk+1) ) 
     199               pvn(ji,jj,jk) = pvn(ji,jj,jk) + ( zsfv(ji,jj,jk) - zsfv(ji,jj,jk+1) ) 
     200            END DO 
     201         END DO 
     202      END DO 
     203      DO jk = 2, jpkm1 
     204         DO jj = 2, jpjm1 
     205            DO ji = fs_2, fs_jpim1   ! vector opt. 
     206               pwn(ji,jj,jk) = pwn(ji,jj,jk) - (  zsfu(ji,jj,jk) - zsfu(ji-1,jj  ,jk)                        & 
     207                  &                             + zsfv(ji,jj,jk) - zsfv(ji  ,jj-1,jk)  ) * tmask(ji,jj,jk) 
     208            END DO 
     209         END DO 
     210      END DO 
     211      ! 
     212   END SUBROUTINE tra_adv_eiv2 
    132213 
    133214#else 
Note: See TracChangeset for help on using the changeset viewer.