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 12928 for NEMO/branches/2019/dev_r11078_OSMOSIS_IMMERSE_Nurser/src/OCE/TRA/traldf_lap_blp.F90 – NEMO

Ignore:
Timestamp:
2020-05-14T21:46:00+02:00 (4 years ago)
Author:
smueller
Message:

Synchronizing with /NEMO/trunk@12925 (ticket #2170)

Location:
NEMO/branches/2019/dev_r11078_OSMOSIS_IMMERSE_Nurser
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r11078_OSMOSIS_IMMERSE_Nurser

    • Property svn:externals
      •  

        old new  
        66^/vendors/FCM@HEAD            ext/FCM 
        77^/vendors/IOIPSL@HEAD         ext/IOIPSL 
         8 
         9# SETTE 
         10^/utils/CI/sette@HEAD         sette 
  • NEMO/branches/2019/dev_r11078_OSMOSIS_IMMERSE_Nurser/src/OCE/TRA/traldf_lap_blp.F90

    r10425 r12928  
    3737 
    3838   !! * Substitutions 
    39 #  include "vectopt_loop_substitute.h90" 
     39#  include "do_loop_substitute.h90" 
    4040   !!---------------------------------------------------------------------- 
    4141   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    4545CONTAINS 
    4646 
    47    SUBROUTINE tra_ldf_lap( kt, kit000, cdtype, pahu, pahv, pgu , pgv ,   & 
    48       &                                                   pgui, pgvi,   & 
    49       &                                        ptb , pta , kjpt, kpass )  
     47   SUBROUTINE tra_ldf_lap( kt, Kmm, kit000, cdtype, pahu, pahv  ,               & 
     48      &                                             pgu , pgv   , pgui, pgvi,   & 
     49      &                                             pt  , pt_rhs, kjpt, kpass )  
    5050      !!---------------------------------------------------------------------- 
    5151      !!                  ***  ROUTINE tra_ldf_lap  *** 
     
    5959      !!          difft = 1/(e1e2t*e3t) {  di-1[ pahu e2u*e3u/e1u di(tb) ] 
    6060      !!                                 + dj-1[ pahv e1v*e3v/e2v dj(tb) ] } 
    61       !!      Add this trend to the general tracer trend pta : 
    62       !!          pta = pta + difft 
    63       !! 
    64       !! ** Action  : - Update pta arrays with the before iso-level  
     61      !!      Add this trend to the general tracer trend pt_rhs : 
     62      !!          pt_rhs = pt_rhs + difft 
     63      !! 
     64      !! ** Action  : - Update pt_rhs arrays with the before iso-level  
    6565      !!                harmonic mixing trend. 
    6666      !!---------------------------------------------------------------------- 
     
    7070      INTEGER                              , INTENT(in   ) ::   kjpt       ! number of tracers 
    7171      INTEGER                              , INTENT(in   ) ::   kpass      ! =1/2 first or second passage 
     72      INTEGER                              , INTENT(in   ) ::   Kmm        ! ocean time level index 
    7273      REAL(wp), DIMENSION(jpi,jpj,jpk)     , INTENT(in   ) ::   pahu, pahv ! eddy diffusivity at u- and v-points  [m2/s] 
    7374      REAL(wp), DIMENSION(jpi,jpj    ,kjpt), INTENT(in   ) ::   pgu, pgv   ! tracer gradient at pstep levels 
    7475      REAL(wp), DIMENSION(jpi,jpj,    kjpt), INTENT(in   ) ::   pgui, pgvi ! tracer gradient at top   levels 
    75       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb        ! before and now tracer fields 
    76       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta        ! tracer trend  
     76      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   pt         ! before tracer fields 
     77      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pt_rhs     ! tracer trend  
    7778      ! 
    7879      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices 
     
    8990      l_hst = .FALSE. 
    9091      l_ptr = .FALSE. 
    91       IF( cdtype == 'TRA' .AND. ln_diaptr )                                                l_ptr = .TRUE.  
     92      IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) )     l_ptr = .TRUE. 
    9293      IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. & 
    9394         &                        iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) )  l_hst = .TRUE. 
     
    9798      ELSE                    ;   zsign = -1._wp 
    9899      ENDIF 
    99       DO jk = 1, jpkm1 
    100          DO jj = 1, jpjm1 
    101             DO ji = 1, fs_jpim1   ! vector opt. 
    102                zaheeu(ji,jj,jk) = zsign * pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u_n(ji,jj,jk)   !!gm   * umask(ji,jj,jk) pah masked! 
    103                zaheev(ji,jj,jk) = zsign * pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v_n(ji,jj,jk)   !!gm   * vmask(ji,jj,jk) 
    104             END DO 
    105          END DO 
    106       END DO 
     100      DO_3D_10_10( 1, jpkm1 ) 
     101         zaheeu(ji,jj,jk) = zsign * pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm)   !!gm   * umask(ji,jj,jk) pah masked! 
     102         zaheev(ji,jj,jk) = zsign * pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v(ji,jj,jk,Kmm)   !!gm   * vmask(ji,jj,jk) 
     103      END_3D 
    107104      ! 
    108105      !                             ! =========== ! 
     
    110107         !                          ! =========== !     
    111108         !                                
    112          DO jk = 1, jpkm1              !== First derivative (gradient)  ==! 
    113             DO jj = 1, jpjm1 
    114                DO ji = 1, fs_jpim1 
    115                   ztu(ji,jj,jk) = zaheeu(ji,jj,jk) * ( ptb(ji+1,jj  ,jk,jn) - ptb(ji,jj,jk,jn) ) 
    116                   ztv(ji,jj,jk) = zaheev(ji,jj,jk) * ( ptb(ji  ,jj+1,jk,jn) - ptb(ji,jj,jk,jn) ) 
    117                END DO 
    118             END DO 
    119          END DO   
     109         DO_3D_10_10( 1, jpkm1 ) 
     110            ztu(ji,jj,jk) = zaheeu(ji,jj,jk) * ( pt(ji+1,jj  ,jk,jn) - pt(ji,jj,jk,jn) ) 
     111            ztv(ji,jj,jk) = zaheev(ji,jj,jk) * ( pt(ji  ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) 
     112         END_3D 
    120113         IF( ln_zps ) THEN                ! set gradient at bottom/top ocean level 
    121             DO jj = 1, jpjm1                    ! bottom 
    122                DO ji = 1, fs_jpim1 
    123                   ztu(ji,jj,mbku(ji,jj)) = zaheeu(ji,jj,mbku(ji,jj)) * pgu(ji,jj,jn) 
    124                   ztv(ji,jj,mbkv(ji,jj)) = zaheev(ji,jj,mbkv(ji,jj)) * pgv(ji,jj,jn) 
    125                END DO 
    126             END DO   
     114            DO_2D_10_10 
     115               ztu(ji,jj,mbku(ji,jj)) = zaheeu(ji,jj,mbku(ji,jj)) * pgu(ji,jj,jn) 
     116               ztv(ji,jj,mbkv(ji,jj)) = zaheev(ji,jj,mbkv(ji,jj)) * pgv(ji,jj,jn) 
     117            END_2D 
    127118            IF( ln_isfcav ) THEN                ! top in ocean cavities only 
    128                DO jj = 1, jpjm1 
    129                   DO ji = 1, fs_jpim1   ! vector opt. 
    130                      IF( miku(ji,jj) > 1 )   ztu(ji,jj,miku(ji,jj)) = zaheeu(ji,jj,miku(ji,jj)) * pgui(ji,jj,jn)  
    131                      IF( mikv(ji,jj) > 1 )   ztv(ji,jj,mikv(ji,jj)) = zaheev(ji,jj,mikv(ji,jj)) * pgvi(ji,jj,jn)  
    132                   END DO 
    133                END DO 
     119               DO_2D_10_10 
     120                  IF( miku(ji,jj) > 1 )   ztu(ji,jj,miku(ji,jj)) = zaheeu(ji,jj,miku(ji,jj)) * pgui(ji,jj,jn)  
     121                  IF( mikv(ji,jj) > 1 )   ztv(ji,jj,mikv(ji,jj)) = zaheev(ji,jj,mikv(ji,jj)) * pgvi(ji,jj,jn)  
     122               END_2D 
    134123            ENDIF 
    135124         ENDIF 
    136125         ! 
    137          DO jk = 1, jpkm1              !== Second derivative (divergence) added to the general tracer trends  ==! 
    138             DO jj = 2, jpjm1 
    139                DO ji = fs_2, fs_jpim1 
    140                   pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + (  ztu(ji,jj,jk) - ztu(ji-1,jj,jk)     & 
    141                      &                                   + ztv(ji,jj,jk) - ztv(ji,jj-1,jk) )   & 
    142                      &                                / ( e1e2t(ji,jj) * e3t_n(ji,jj,jk) ) 
    143                END DO 
    144             END DO 
    145          END DO   
     126         DO_3D_00_00( 1, jpkm1 ) 
     127            pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + (  ztu(ji,jj,jk) - ztu(ji-1,jj,jk)     & 
     128               &                                      +    ztv(ji,jj,jk) - ztv(ji,jj-1,jk) )   & 
     129               &                                      / ( e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) ) 
     130         END_3D 
    146131         ! 
    147132         !                             !== "Poleward" diffusive heat or salt transports  ==! 
     
    159144    
    160145 
    161    SUBROUTINE tra_ldf_blp( kt, kit000, cdtype, pahu, pahv, pgu , pgv ,   & 
    162       &                                                    pgui, pgvi,  & 
    163       &                                                    ptb , pta , kjpt, kldf ) 
     146   SUBROUTINE tra_ldf_blp( kt, Kmm, kit000, cdtype, pahu, pahv  ,             & 
     147      &                                             pgu , pgv   , pgui, pgvi, & 
     148      &                                             pt  , pt_rhs, kjpt, kldf ) 
    164149      !!---------------------------------------------------------------------- 
    165150      !!                 ***  ROUTINE tra_ldf_blp  *** 
     
    179164      INTEGER                              , INTENT(in   ) ::   kjpt       ! number of tracers 
    180165      INTEGER                              , INTENT(in   ) ::   kldf       ! type of operator used 
     166      INTEGER                              , INTENT(in   ) ::   Kmm        ! ocean time level indices 
    181167      REAL(wp), DIMENSION(jpi,jpj,jpk)     , INTENT(in   ) ::   pahu, pahv ! eddy diffusivity at u- and v-points  [m2/s] 
    182168      REAL(wp), DIMENSION(jpi,jpj    ,kjpt), INTENT(in   ) ::   pgu, pgv   ! tracer gradient at pstep levels 
    183169      REAL(wp), DIMENSION(jpi,jpj,    kjpt), INTENT(in   ) ::   pgui, pgvi ! tracer gradient at top levels 
    184       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb        ! before and now tracer fields 
    185       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta        ! tracer trend 
     170      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   pt         ! before and now tracer fields 
     171      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pt_rhs     ! tracer trend 
    186172      ! 
    187173      INTEGER ::   ji, jj, jk, jn   ! dummy loop indices 
     
    203189      zlap(:,:,:,:) = 0._wp 
    204190      ! 
    205       SELECT CASE ( kldf )       !==  1st laplacian applied to ptb (output in zlap)  ==! 
     191      SELECT CASE ( kldf )       !==  1st laplacian applied to pt (output in zlap)  ==! 
    206192      ! 
    207193      CASE ( np_blp    )               ! iso-level bilaplacian 
    208          CALL tra_ldf_lap  ( kt, kit000, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, ptb,      zlap, kjpt, 1 ) 
     194         CALL tra_ldf_lap  ( kt, Kmm, kit000, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, pt,     zlap, kjpt, 1 ) 
    209195      CASE ( np_blp_i  )               ! rotated   bilaplacian : standard operator (Madec) 
    210          CALL tra_ldf_iso  ( kt, kit000, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, ptb, ptb, zlap, kjpt, 1 ) 
     196         CALL tra_ldf_iso  ( kt, Kmm, kit000, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, pt, pt, zlap, kjpt, 1 ) 
    211197      CASE ( np_blp_it )               ! rotated  bilaplacian : triad operator (griffies) 
    212          CALL tra_ldf_triad( kt, kit000, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, ptb, ptb, zlap, kjpt, 1 ) 
     198         CALL tra_ldf_triad( kt, Kmm, kit000, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, pt, pt, zlap, kjpt, 1 ) 
    213199      END SELECT 
    214200      ! 
    215201      CALL lbc_lnk( 'traldf_lap_blp', zlap(:,:,:,:) , 'T', 1. )     ! Lateral boundary conditions (unchanged sign) 
    216202      !                                               ! Partial top/bottom cell: GRADh( zlap )   
    217       IF( ln_isfcav .AND. ln_zps ) THEN   ;   CALL zps_hde_isf( kt, kjpt, zlap, zglu, zglv, zgui, zgvi )  ! both top & bottom 
    218       ELSEIF(             ln_zps ) THEN   ;   CALL zps_hde    ( kt, kjpt, zlap, zglu, zglv )              ! only bottom  
    219       ENDIF 
    220       ! 
    221       SELECT CASE ( kldf )       !==  2nd laplacian applied to zlap (output in pta)  ==! 
     203      IF( ln_isfcav .AND. ln_zps ) THEN   ;   CALL zps_hde_isf( kt, Kmm, kjpt, zlap, zglu, zglv, zgui, zgvi )  ! both top & bottom 
     204      ELSEIF(             ln_zps ) THEN   ;   CALL zps_hde    ( kt, Kmm, kjpt, zlap, zglu, zglv )              ! only bottom  
     205      ENDIF 
     206      ! 
     207      SELECT CASE ( kldf )       !==  2nd laplacian applied to zlap (output in pt_rhs)  ==! 
    222208      ! 
    223209      CASE ( np_blp    )               ! iso-level bilaplacian 
    224          CALL tra_ldf_lap  ( kt, kit000, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, pta,      kjpt, 2 ) 
     210         CALL tra_ldf_lap  ( kt, Kmm, kit000, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, pt_rhs,         kjpt, 2 ) 
    225211      CASE ( np_blp_i  )               ! rotated   bilaplacian : standard operator (Madec) 
    226          CALL tra_ldf_iso  ( kt, kit000, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, ptb, pta, kjpt, 2 ) 
     212         CALL tra_ldf_iso  ( kt, Kmm, kit000, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, pt    , pt_rhs, kjpt, 2 ) 
    227213      CASE ( np_blp_it )               ! rotated  bilaplacian : triad operator (griffies) 
    228          CALL tra_ldf_triad( kt, kit000, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, ptb, pta, kjpt, 2 ) 
     214         CALL tra_ldf_triad( kt, Kmm, kit000, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, pt    , pt_rhs, kjpt, 2 ) 
    229215      END SELECT 
    230216      ! 
Note: See TracChangeset for help on using the changeset viewer.