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 792 – NEMO

Changeset 792


Ignore:
Timestamp:
2008-01-13T13:05:05+01:00 (16 years ago)
Author:
gm
Message:

dev_001_GM - switch from velocity to transport in traadv ==> update the trend computation

Location:
branches/dev_001_GM/NEMO/OPA_SRC/TRD
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • branches/dev_001_GM/NEMO/OPA_SRC/TRD/trdicp.F90

    r790 r792  
    44   !! Ocean diagnostics:  ocean tracers and dynamic trends 
    55   !!===================================================================== 
    6    !! History :       !  91-12 (G. Madec) 
    7    !!                 !  92-06 (M. Imbard) add time step frequency 
    8    !!                 !  96-01 (G. Madec)  terrain following coordinates 
    9    !!            8.5  !  02-06 (G. Madec)  F90: Free form and module 
    10    !!            9.0  !  04-08 (C. Talandier) New trends organization 
     6   !! History :  OPA  !  1991-12  (G. Madec) 
     7   !!            4.0  !  1992-06  (M. Imbard) add time step frequency 
     8   !!            8.0  !  1996-01  (G. Madec)  terrain following coordinates 
     9   !!   NEMO     1.0  !  2002-06  (G. Madec)  F90: Free form and module 
     10   !!            1.0  !  2004-08  (C. Talandier) New trends organization 
     11   !!            2.4  !  2008-01  (G. Madec)  merge TRC-TRA 
    1112   !!---------------------------------------------------------------------- 
    1213#if  defined key_trdtra   ||   defined key_trddyn   ||   defined key_esopa 
  • branches/dev_001_GM/NEMO/OPA_SRC/TRD/trdmod.F90

    r790 r792  
    149149      !!                  ***  ROUTINE trd_mod  *** 
    150150      !!  
    151       !! ** Purpose :   transformed the i-advective flux into the i-advective trends 
    152       !! ** Method  :   i-advective trends = un. di[T] = di[fi] - tn di[un] 
    153       !!---------------------------------------------------------------------- 
    154       INTEGER         , INTENT(in   )                         ::   kt      ! time step 
    155       INTEGER         , INTENT(in   )                         ::   ktra    ! tracer index 
    156       INTEGER         , INTENT(in   )                         ::   ktrd    ! tracer trend index 
    157       CHARACTER(len=3), INTENT(in   )                         ::   ctype   ! tracers type 'TRA' or 'TRC' 
    158       REAL(wp)        , INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   pf      ! advective flux in one direction 
    159       REAL(wp)        , INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   pun     ! now velocity  in one direction 
    160       REAL(wp)        , INTENT(in   ), DIMENSION(jpi,jpj,jpk) ::   ptn     ! now or before tracer  
    161       CHARACTER(len=3), INTENT(in   )                        , OPTIONAL ::   cnbpas  ! number of passage 
     151      !! ** Purpose :   transformed the i-, j- or k-advective flux into thes 
     152      !!              i-, j- or k-advective trends, resp. 
     153      !! ** Method  :   i-advective trends = -un. di-1[T] = -( di-1[fi] - tn di-1[un] ) 
     154      !!                k-advective trends = -un. di-1[T] = -( dj-1[fi] - tn dj-1[un] ) 
     155      !!                k-advective trends = -un. di+1[T] = -( dk+1[fi] - tn dk+1[un] ) 
     156      !!   NB: for TVD advection scheme, this routine is called twice (for upstream 
     157      !!       advective fluxes and anti-diffusive fluxes, resp.), but the 
     158      !!       divergence is only removed once, at the first call. 
     159      !!---------------------------------------------------------------------- 
     160      INTEGER         , INTENT(in)                                   ::   kt      ! time step 
     161      INTEGER         , INTENT(in)                                   ::   ktra    ! tracer index 
     162      INTEGER         , INTENT(in)                                   ::   ktrd    ! tracer trend index 
     163      CHARACTER(len=3), INTENT(in)                                   ::   ctype   ! tracers type 'TRA' or 'TRC' 
     164      REAL(wp)        , INTENT(in), DIMENSION(jpi,jpj,jpk)           ::   pf      ! advective flux in one direction 
     165      REAL(wp)        , INTENT(in), DIMENSION(jpi,jpj,jpk)           ::   pun     ! now velocity  in one direction 
     166      REAL(wp)        , INTENT(in), DIMENSION(jpi,jpj,jpk)           ::   ptn     ! now or before tracer  
     167      CHARACTER(len=3), INTENT(in)                        , OPTIONAL ::   cnbpas  ! number of passage 
    162168      !! 
    163       INTEGER                          ::   ji, jj, jk      ! dummy loop indices 
    164       CHARACTER(len=3) ::   ccpas                           ! number of passage 
    165       REAL(wp)                         ::   zbtr, z_hdivn   !  
    166       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   ztrdt           ! 3D workspace  
     169      INTEGER                          ::   ji, jj, jk   ! dummy loop indices 
     170      INTEGER                          ::   ii, ij, ik   ! index shift function of the direction 
     171      CHARACTER(len=3)                 ::   ccpas        ! number of passage 
     172      REAL(wp)                         ::   zbtr         ! temporary scalar 
     173      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   ztrdt        ! 3D workspace  
    167174      !!---------------------------------------------------------------------- 
    168175 
     
    170177      IF( PRESENT(cnbpas) )  ccpas = cnbpas 
    171178 
    172       ztrdt(:,:,:) = 0.e0 
     179 
     180      SELECT CASE( ktrd )            ! shift depending on the direction 
     181      CASE( jpt_trd_xad )   ;   ii = 1   ; ij = 0   ;   ik = 0      ! i-advective trend 
     182      CASE( jpt_trd_yad )   ;   ii = 0   ; ij = 1   ;   ik = 0      ! j-advective trend 
     183      CASE( jpt_trd_zad )   ;   ii = 0   ; ij = 0   ;   ik =-1      ! k-advective trend 
     184      END SELECT 
     185 
     186      !                              ! set to zero uncomputed values 
     187      ztrdt(jpi,:,:) = 0.e0   ;   ztrdt(1,:,:) = 0.e0 
     188      ztrdt(:,jpj,:) = 0.e0   ;   ztrdt(:,1,:) = 0.e0 
     189      ztrdt(:,:,jpk) = 0.e0 
     190      ! 
    173191      ! 
    174192      IF( ccpas == 'fst' ) THEN      ! first treatment : remove the divergence 
    175          SELECT CASE( ktrd )  
    176          CASE( jpt_trd_xad )      ! i-advective trend 
    177             DO jk = 1, jpkm1 
    178                DO jj = 2, jpjm1 
    179                   DO ji = fs_2, fs_jpim1   ! vector opt. 
    180 #if defined key_zco 
    181                      zbtr    = 1.e0 / ( e1t(ji,jj) * e2t(ji,jj) ) 
    182                      z_hdivn = (  e2u(ji  ,jj) * pun(ji  ,jj,jk)                              & 
    183                         &       - e2u(ji-1,jj) * pun(ji-1,jj,jk) ) 
    184 #else 
    185                      zbtr    = 1.e0/ ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
    186                      z_hdivn = (  e2u(ji  ,jj) * fse3u(ji  ,jj,jk) * pun(ji  ,jj,jk)          & 
    187                         &       - e2u(ji-1,jj) * fse3u(ji-1,jj,jk) * pun(ji-1,jj,jk) ) 
    188 #endif 
    189                      ztrdt(ji,jj,jk) = - zbtr * (  pf(ji,jj,jk) - pf(ji-1,jj,jk)  -  ptn(ji,jj,jk) * z_hdivn  ) 
    190                   END DO 
     193         DO jk = 1, jpkm1 
     194            DO jj = 2, jpjm1 
     195               DO ji = fs_2, fs_jpim1   ! vector opt. 
     196                  zbtr    = 1.e0/ ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
     197                  ztrdt(ji,jj,jk) = - zbtr * (      pf (ji,jj,jk) - pf (ji-ii,jj-ij,jk-ik)                    & 
     198                     &                          - ( pun(ji,jj,jk) - pun(ji-ii,jj-ij,jk-ik) ) * ptn(ji,jj,jk)  ) 
    191199               END DO 
    192200            END DO 
    193             ! 
    194          CASE( jpt_trd_yad )      ! j-advective trend 
    195             DO jk = 1, jpkm1 
    196                DO jj = 2, jpjm1 
    197                   DO ji = fs_2, fs_jpim1   ! vector opt. 
    198 #if defined key_zco 
    199                      zbtr    = 1.e0 / ( e1t(ji,jj) * e2t(ji,jj) ) 
    200                      z_hdivn = (  e1v(ji,jj  ) * pun(ji,jj  ,jk)                              & 
    201                         &       - e1v(ji,jj-1) * pun(ji,jj-1,jk) ) 
    202 #else 
    203                      zbtr    = 1.e0/ ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
    204                      z_hdivn = (  e1v(ji  ,jj) * fse3v(ji,jj  ,jk) * pun(ji,jj  ,jk)          & 
    205                         &       - e1v(ji-1,jj) * fse3v(ji,jj-1,jk) * pun(ji,jj-1,jk) ) 
    206 #endif 
    207                      ztrdt(ji,jj,jk) = - zbtr * (  pf(ji,jj,jk) - pf(ji,jj-1,jk)  -  ptn(ji,jj,jk) * z_hdivn  ) 
    208                   END DO 
     201         END DO 
     202         ! 
     203      ELSE                           ! second call : just compute the trend   (TVD scheme) 
     204         DO jk = 1, jpkm1 
     205            DO jj = 2, jpjm1 
     206               DO ji = fs_2, fs_jpim1   ! vector opt. 
     207                  zbtr    = 1.e0/ ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
     208                  ztrdt(ji,jj,jk) = - zbtr * (      pf (ji,jj,jk) - pf (ji-ii,jj-ij,jk-ik)     ) 
    209209               END DO 
    210210            END DO 
    211             ! 
    212          CASE( jpt_trd_zad )      ! z-advective trend 
    213             DO jk = 1, jpkm1 
    214                DO jj = 2, jpjm1 
    215                   DO ji = fs_2, fs_jpim1   ! vector opt. 
    216                      zbtr    = 1.e0 / fse3t(ji,jj,jk) 
    217                      z_hdivn = pun(ji,jj,jk) - pun(ji,jj,jk+1) 
    218                      ztrdt(ji,jj,jk) = - zbtr * (  pf(ji,jj,jk) - pf(ji,jj,jk+1)  -  ptn(ji,jj,jk) * z_hdivn  ) 
    219                   END DO 
    220                END DO 
    221             END DO 
    222             ! 
    223          END SELECT 
    224          ! 
    225       ELSE                ! second call : just compute the trend   (TVD scheme) 
    226          SELECT CASE( ktrd )  
    227          CASE( jpt_trd_xad )      ! i-advective trend 
    228             DO jk = 1, jpkm1       
    229                DO jj = 2, jpjm1 
    230                   DO ji = fs_2, fs_jpim1   ! vector opt. 
    231                      zbtr    = 1.e0 / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
    232                      ztrdt(ji,jj,jk) = - zbtr * (  pf(ji,jj,jk) - pf(ji-1,jj,jk)  ) 
    233                   END DO 
    234                END DO 
    235             END DO 
    236             ! 
    237          CASE( jpt_trd_yad )      ! j-advective trend 
    238             DO jk = 1, jpkm1 
    239                DO jj = 2, jpjm1 
    240                   DO ji = fs_2, fs_jpim1   ! vector opt. 
    241                      zbtr    = 1.e0 / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) ) 
    242                      ztrdt(ji,jj,jk) = - zbtr * (  pf(ji,jj,jk) - pf(ji,jj-1,jk)  ) 
    243                   END DO 
    244                END DO 
    245             END DO 
    246             ! 
    247          CASE( jpt_trd_zad )      ! z-advective trend 
    248             DO jk = 1, jpkm1 
    249                DO jj = 2, jpjm1 
    250                   DO ji = fs_2, fs_jpim1   ! vector opt. 
    251                      zbtr    = 1.e0 / fse3t(ji,jj,jk) 
    252                      ztrdt(ji,jj,jk) = - zbtr * (  pf(ji,jj,jk) - pf(ji,jj,jk+1)  ) 
    253                   END DO 
    254                END DO 
    255             END DO 
    256             ! 
    257          END SELECT 
    258          ! 
    259       ENDIF 
    260       ! 
    261       CALL trd_tra( kt, ktra, ktrd, ctype, ptrd3d=ztrdt)       ! trend diagnostics 
     211         END DO 
     212         ! 
     213      ENDIF 
     214      ! 
     215      !                              ! trend diagnostics 
     216      CALL trd_tra( kt, ktra, ktrd, ctype, ptrd3d=ztrdt, cnbpas=ccpas ) 
    262217      ! 
    263218   END SUBROUTINE trd_tra_adv 
Note: See TracChangeset for help on using the changeset viewer.