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 12377 for NEMO/trunk/src/OCE/TRA/trazdf.F90 – NEMO

Ignore:
Timestamp:
2020-02-12T15:39:06+01:00 (4 years ago)
Author:
acc
Message:

The big one. Merging all 2019 developments from the option 1 branch back onto the trunk.

This changeset reproduces 2019/dev_r11943_MERGE_2019 on the trunk using a 2-URL merge
onto a working copy of the trunk. I.e.:

svn merge --ignore-ancestry \

svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/trunk \
svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/branches/2019/dev_r11943_MERGE_2019 ./

The --ignore-ancestry flag avoids problems that may otherwise arise from the fact that
the merge history been trunk and branch may have been applied in a different order but
care has been taken before this step to ensure that all applicable fixes and updates
are present in the merge branch.

The trunk state just before this step has been branched to releases/release-4.0-HEAD
and that branch has been immediately tagged as releases/release-4.0.2. Any fixes
or additions in response to tickets on 4.0, 4.0.1 or 4.0.2 should be done on
releases/release-4.0-HEAD. From now on future 'point' releases (e.g. 4.0.2) will
remain unchanged with periodic releases as needs demand. Note release-4.0-HEAD is a
transitional naming convention. Future full releases, say 4.2, will have a release-4.2
branch which fulfills this role and the first point release (e.g. 4.2.0) will be made
immediately following the release branch creation.

2020 developments can be started from any trunk revision later than this one.

Location:
NEMO/trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk

    • Property svn:externals
      •  

        old new  
        33^/utils/build/mk@HEAD         mk 
        44^/utils/tools@HEAD            tools 
        5 ^/vendors/AGRIF/dev@HEAD      ext/AGRIF 
         5^/vendors/AGRIF/dev_r11615_ENHANCE-04_namelists_as_internalfiles_agrif@HEAD      ext/AGRIF 
        66^/vendors/FCM@HEAD            ext/FCM 
        77^/vendors/IOIPSL@HEAD         ext/IOIPSL 
  • NEMO/trunk/src/OCE/TRA/trazdf.F90

    r10425 r12377  
    3636 
    3737   !! * Substitutions 
    38 #  include "vectopt_loop_substitute.h90" 
     38#  include "do_loop_substitute.h90" 
    3939   !!---------------------------------------------------------------------- 
    4040   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    4444CONTAINS 
    4545 
    46    SUBROUTINE tra_zdf( kt ) 
     46   SUBROUTINE tra_zdf( kt, Kbb, Kmm, Krhs, pts, Kaa ) 
    4747      !!---------------------------------------------------------------------- 
    4848      !!                  ***  ROUTINE tra_zdf  *** 
     
    5050      !! ** Purpose :   compute the vertical ocean tracer physics. 
    5151      !!--------------------------------------------------------------------- 
    52       INTEGER, INTENT(in) ::   kt   ! ocean time-step index 
     52      INTEGER                                  , INTENT(in)    :: kt                  ! ocean time-step index 
     53      INTEGER                                  , INTENT(in)    :: Kbb, Kmm, Krhs, Kaa ! time level indices 
     54      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts                 ! active tracers and RHS of tracer equation 
    5355      ! 
    5456      INTEGER  ::   jk   ! Dummy loop indices 
     
    7072      IF( l_trdtra )   THEN                  !* Save ta and sa trends 
    7173         ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) ) 
    72          ztrdt(:,:,:) = tsa(:,:,:,jp_tem) 
    73          ztrds(:,:,:) = tsa(:,:,:,jp_sal) 
     74         ztrdt(:,:,:) = pts(:,:,:,jp_tem,Kaa) 
     75         ztrds(:,:,:) = pts(:,:,:,jp_sal,Kaa) 
    7476      ENDIF 
    7577      ! 
    7678      !                                      !* compute lateral mixing trend and add it to the general trend 
    77       CALL tra_zdf_imp( kt, nit000, 'TRA', r2dt, tsb, tsa, jpts )  
     79      CALL tra_zdf_imp( kt, nit000, 'TRA', r2dt, Kbb, Kmm, Krhs, pts, Kaa, jpts )  
    7880 
    7981!!gm WHY here !   and I don't like that ! 
     
    8183      ! JMM avoid negative salinities near river outlet ! Ugly fix 
    8284      ! JMM : restore negative salinities to small salinities: 
    83       WHERE( tsa(:,:,:,jp_sal) < 0._wp )   tsa(:,:,:,jp_sal) = 0.1_wp 
     85      WHERE( pts(:,:,:,jp_sal,Kaa) < 0._wp )   pts(:,:,:,jp_sal,Kaa) = 0.1_wp 
    8486!!gm 
    8587 
    8688      IF( l_trdtra )   THEN                      ! save the vertical diffusive trends for further diagnostics 
    8789         DO jk = 1, jpkm1 
    88             ztrdt(:,:,jk) = ( ( tsa(:,:,jk,jp_tem)*e3t_a(:,:,jk) - tsb(:,:,jk,jp_tem)*e3t_b(:,:,jk) ) & 
    89                &          / (e3t_n(:,:,jk)*r2dt) ) - ztrdt(:,:,jk) 
    90             ztrds(:,:,jk) = ( ( tsa(:,:,jk,jp_sal)*e3t_a(:,:,jk) - tsb(:,:,jk,jp_sal)*e3t_b(:,:,jk) ) & 
    91               &           / (e3t_n(:,:,jk)*r2dt) ) - ztrds(:,:,jk) 
     90            ztrdt(:,:,jk) = ( ( pts(:,:,jk,jp_tem,Kaa)*e3t(:,:,jk,Kaa) - pts(:,:,jk,jp_tem,Kbb)*e3t(:,:,jk,Kbb) ) & 
     91               &          / (e3t(:,:,jk,Kmm)*r2dt) ) - ztrdt(:,:,jk) 
     92            ztrds(:,:,jk) = ( ( pts(:,:,jk,jp_sal,Kaa)*e3t(:,:,jk,Kaa) - pts(:,:,jk,jp_sal,Kbb)*e3t(:,:,jk,Kbb) ) & 
     93              &           / (e3t(:,:,jk,Kmm)*r2dt) ) - ztrds(:,:,jk) 
    9294         END DO 
    9395!!gm this should be moved in trdtra.F90 and done on all trends 
    9496         CALL lbc_lnk_multi( 'trazdf', ztrdt, 'T', 1. , ztrds, 'T', 1. ) 
    9597!!gm 
    96          CALL trd_tra( kt, 'TRA', jp_tem, jptra_zdf, ztrdt ) 
    97          CALL trd_tra( kt, 'TRA', jp_sal, jptra_zdf, ztrds ) 
     98         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_zdf, ztrdt ) 
     99         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_zdf, ztrds ) 
    98100         DEALLOCATE( ztrdt , ztrds ) 
    99101      ENDIF 
    100102      !                                          ! print mean trends (used for debugging) 
    101       IF(ln_ctl)   CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' zdf  - Ta: ', mask1=tmask,               & 
    102          &                       tab3d_2=tsa(:,:,:,jp_sal), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' ) 
     103      IF(sn_cfctl%l_prtctl)   CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Kaa), clinfo1=' zdf  - Ta: ', mask1=tmask,               & 
     104         &                                  tab3d_2=pts(:,:,:,jp_sal,Kaa), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' ) 
    103105      ! 
    104106      IF( ln_timing )   CALL timing_stop('tra_zdf') 
     
    107109 
    108110  
    109    SUBROUTINE tra_zdf_imp( kt, kit000, cdtype, p2dt, ptb, pta, kjpt )  
     111   SUBROUTINE tra_zdf_imp( kt, kit000, cdtype, p2dt, Kbb, Kmm, Krhs, pt, Kaa, kjpt )  
    110112      !!---------------------------------------------------------------------- 
    111113      !!                  ***  ROUTINE tra_zdf_imp  *** 
     
    125127      !!      If iso-neutral mixing, add to avt the contribution due to lateral mixing. 
    126128      !! 
    127       !! ** Action  : - pta  becomes the after tracer 
    128       !!--------------------------------------------------------------------- 
    129       INTEGER                              , INTENT(in   ) ::   kt       ! ocean time-step index 
    130       INTEGER                              , INTENT(in   ) ::   kit000   ! first time step index 
    131       CHARACTER(len=3)                     , INTENT(in   ) ::   cdtype   ! =TRA or TRC (tracer indicator) 
    132       INTEGER                              , INTENT(in   ) ::   kjpt     ! number of tracers 
    133       REAL(wp)                             , INTENT(in   ) ::   p2dt     ! tracer time-step 
    134       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb      ! before and now tracer fields 
    135       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta      ! in: tracer trend ; out: after tracer field 
     129      !! ** Action  : - pt(:,:,:,:,Kaa)  becomes the after tracer 
     130      !!--------------------------------------------------------------------- 
     131      INTEGER                                  , INTENT(in   ) ::   kt       ! ocean time-step index 
     132      INTEGER                                  , INTENT(in   ) ::   Kbb, Kmm, Krhs, Kaa  ! ocean time level indices 
     133      INTEGER                                  , INTENT(in   ) ::   kit000   ! first time step index 
     134      CHARACTER(len=3)                         , INTENT(in   ) ::   cdtype   ! =TRA or TRC (tracer indicator) 
     135      INTEGER                                  , INTENT(in   ) ::   kjpt     ! number of tracers 
     136      REAL(wp)                                 , INTENT(in   ) ::   p2dt     ! tracer time-step 
     137      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) ::   pt       ! tracers and RHS of tracer equation 
    136138      ! 
    137139      INTEGER  ::  ji, jj, jk, jn   ! dummy loop indices 
     
    158160            IF( l_ldfslp ) THEN            ! isoneutral diffusion: add the contribution  
    159161               IF( ln_traldf_msc  ) THEN     ! MSC iso-neutral operator  
    160                   DO jk = 2, jpkm1 
    161                      DO jj = 2, jpjm1 
    162                         DO ji = fs_2, fs_jpim1   ! vector opt. 
    163                            zwt(ji,jj,jk) = zwt(ji,jj,jk) + akz(ji,jj,jk)   
    164                         END DO 
    165                      END DO 
    166                   END DO 
     162                  DO_3D_00_00( 2, jpkm1 ) 
     163                     zwt(ji,jj,jk) = zwt(ji,jj,jk) + akz(ji,jj,jk)   
     164                  END_3D 
    167165               ELSE                          ! standard or triad iso-neutral operator 
    168                   DO jk = 2, jpkm1 
    169                      DO jj = 2, jpjm1 
    170                         DO ji = fs_2, fs_jpim1   ! vector opt. 
    171                            zwt(ji,jj,jk) = zwt(ji,jj,jk) + ah_wslp2(ji,jj,jk) 
    172                         END DO 
    173                      END DO 
    174                   END DO 
     166                  DO_3D_00_00( 2, jpkm1 ) 
     167                     zwt(ji,jj,jk) = zwt(ji,jj,jk) + ah_wslp2(ji,jj,jk) 
     168                  END_3D 
    175169               ENDIF 
    176170            ENDIF 
     
    178172            ! Diagonal, lower (i), upper (s)  (including the bottom boundary condition since avt is masked) 
    179173            IF( ln_zad_Aimp ) THEN         ! Adaptive implicit vertical advection 
    180                DO jk = 1, jpkm1 
    181                   DO jj = 2, jpjm1 
    182                      DO ji = fs_2, fs_jpim1   ! vector opt. (ensure same order of calculation as below if wi=0.) 
    183                         zzwi = - p2dt * zwt(ji,jj,jk  ) / e3w_n(ji,jj,jk  ) 
    184                         zzws = - p2dt * zwt(ji,jj,jk+1) / e3w_n(ji,jj,jk+1) 
    185                         zwd(ji,jj,jk) = e3t_a(ji,jj,jk) - zzwi - zzws   & 
    186                            &                 + p2dt * ( MAX( wi(ji,jj,jk  ) , 0._wp ) - MIN( wi(ji,jj,jk+1) , 0._wp ) ) 
    187                         zwi(ji,jj,jk) = zzwi + p2dt *   MIN( wi(ji,jj,jk  ) , 0._wp ) 
    188                         zws(ji,jj,jk) = zzws - p2dt *   MAX( wi(ji,jj,jk+1) , 0._wp ) 
    189                     END DO 
    190                   END DO 
    191                END DO 
     174               DO_3D_00_00( 1, jpkm1 ) 
     175                  zzwi = - p2dt * zwt(ji,jj,jk  ) / e3w(ji,jj,jk  ,Kmm) 
     176                  zzws = - p2dt * zwt(ji,jj,jk+1) / e3w(ji,jj,jk+1,Kmm) 
     177                  zwd(ji,jj,jk) = e3t(ji,jj,jk,Kaa) - zzwi - zzws   & 
     178                     &                 + p2dt * ( MAX( wi(ji,jj,jk  ) , 0._wp ) - MIN( wi(ji,jj,jk+1) , 0._wp ) ) 
     179                  zwi(ji,jj,jk) = zzwi + p2dt *   MIN( wi(ji,jj,jk  ) , 0._wp ) 
     180                  zws(ji,jj,jk) = zzws - p2dt *   MAX( wi(ji,jj,jk+1) , 0._wp ) 
     181               END_3D 
    192182            ELSE 
    193                DO jk = 1, jpkm1 
    194                   DO jj = 2, jpjm1 
    195                      DO ji = fs_2, fs_jpim1   ! vector opt. 
    196                         zwi(ji,jj,jk) = - p2dt * zwt(ji,jj,jk  ) / e3w_n(ji,jj,jk) 
    197                         zws(ji,jj,jk) = - p2dt * zwt(ji,jj,jk+1) / e3w_n(ji,jj,jk+1) 
    198                         zwd(ji,jj,jk) = e3t_a(ji,jj,jk) - zwi(ji,jj,jk) - zws(ji,jj,jk) 
    199                     END DO 
    200                   END DO 
    201                END DO 
     183               DO_3D_00_00( 1, jpkm1 ) 
     184                  zwi(ji,jj,jk) = - p2dt * zwt(ji,jj,jk  ) / e3w(ji,jj,jk,Kmm) 
     185                  zws(ji,jj,jk) = - p2dt * zwt(ji,jj,jk+1) / e3w(ji,jj,jk+1,Kmm) 
     186                  zwd(ji,jj,jk) = e3t(ji,jj,jk,Kaa) - zwi(ji,jj,jk) - zws(ji,jj,jk) 
     187               END_3D 
    202188            ENDIF 
    203189            ! 
     
    218204            !   The solution will be in the 4d array pta. 
    219205            !   The 3d array zwt is used as a work space array. 
    220             !   En route to the solution pta is used a to evaluate the rhs and then  
     206            !   En route to the solution pt(:,:,:,:,Kaa) is used a to evaluate the rhs and then  
    221207            !   used as a work space array: its value is modified. 
    222208            ! 
    223             DO jj = 2, jpjm1        !* 1st recurrence:   Tk = Dk - Ik Sk-1 / Tk-1   (increasing k) 
    224                DO ji = fs_2, fs_jpim1            ! done one for all passive tracers (so included in the IF instruction) 
    225                   zwt(ji,jj,1) = zwd(ji,jj,1) 
    226                END DO 
    227             END DO 
    228             DO jk = 2, jpkm1 
    229                DO jj = 2, jpjm1 
    230                   DO ji = fs_2, fs_jpim1 
    231                      zwt(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) / zwt(ji,jj,jk-1) 
    232                   END DO 
    233                END DO 
    234             END DO 
     209            DO_2D_00_00 
     210               zwt(ji,jj,1) = zwd(ji,jj,1) 
     211            END_2D 
     212            DO_3D_00_00( 2, jpkm1 ) 
     213               zwt(ji,jj,jk) = zwd(ji,jj,jk) - zwi(ji,jj,jk) * zws(ji,jj,jk-1) / zwt(ji,jj,jk-1) 
     214            END_3D 
    235215            ! 
    236216         ENDIF  
    237217         !          
    238          DO jj = 2, jpjm1           !* 2nd recurrence:    Zk = Yk - Ik / Tk-1  Zk-1 
    239             DO ji = fs_2, fs_jpim1 
    240                pta(ji,jj,1,jn) = e3t_b(ji,jj,1) * ptb(ji,jj,1,jn) + p2dt * e3t_n(ji,jj,1) * pta(ji,jj,1,jn) 
    241             END DO 
    242          END DO 
    243          DO jk = 2, jpkm1 
    244             DO jj = 2, jpjm1 
    245                DO ji = fs_2, fs_jpim1 
    246                   zrhs = e3t_b(ji,jj,jk) * ptb(ji,jj,jk,jn) + p2dt * e3t_n(ji,jj,jk) * pta(ji,jj,jk,jn)   ! zrhs=right hand side 
    247                   pta(ji,jj,jk,jn) = zrhs - zwi(ji,jj,jk) / zwt(ji,jj,jk-1) * pta(ji,jj,jk-1,jn) 
    248                END DO 
    249             END DO 
    250          END DO 
     218         DO_2D_00_00 
     219            pt(ji,jj,1,jn,Kaa) = e3t(ji,jj,1,Kbb) * pt(ji,jj,1,jn,Kbb) + p2dt * e3t(ji,jj,1,Kmm) * pt(ji,jj,1,jn,Krhs) 
     220         END_2D 
     221         DO_3D_00_00( 2, jpkm1 ) 
     222            zrhs = e3t(ji,jj,jk,Kbb) * pt(ji,jj,jk,jn,Kbb) + p2dt * e3t(ji,jj,jk,Kmm) * pt(ji,jj,jk,jn,Krhs)   ! zrhs=right hand side 
     223            pt(ji,jj,jk,jn,Kaa) = zrhs - zwi(ji,jj,jk) / zwt(ji,jj,jk-1) * pt(ji,jj,jk-1,jn,Kaa) 
     224         END_3D 
    251225         ! 
    252          DO jj = 2, jpjm1           !* 3d recurrence:    Xk = (Zk - Sk Xk+1 ) / Tk   (result is the after tracer) 
    253             DO ji = fs_2, fs_jpim1 
    254                pta(ji,jj,jpkm1,jn) = pta(ji,jj,jpkm1,jn) / zwt(ji,jj,jpkm1) * tmask(ji,jj,jpkm1) 
    255             END DO 
    256          END DO 
    257          DO jk = jpk-2, 1, -1 
    258             DO jj = 2, jpjm1 
    259                DO ji = fs_2, fs_jpim1 
    260                   pta(ji,jj,jk,jn) = ( pta(ji,jj,jk,jn) - zws(ji,jj,jk) * pta(ji,jj,jk+1,jn) )   & 
    261                      &             / zwt(ji,jj,jk) * tmask(ji,jj,jk) 
    262                END DO 
    263             END DO 
    264          END DO 
     226         DO_2D_00_00 
     227            pt(ji,jj,jpkm1,jn,Kaa) = pt(ji,jj,jpkm1,jn,Kaa) / zwt(ji,jj,jpkm1) * tmask(ji,jj,jpkm1) 
     228         END_2D 
     229         DO_3DS_00_00( jpk-2, 1, -1 ) 
     230            pt(ji,jj,jk,jn,Kaa) = ( pt(ji,jj,jk,jn,Kaa) - zws(ji,jj,jk) * pt(ji,jj,jk+1,jn,Kaa) )   & 
     231               &             / zwt(ji,jj,jk) * tmask(ji,jj,jk) 
     232         END_3D 
    265233         !                                            ! ================= ! 
    266234      END DO                                          !  end tracer loop  ! 
Note: See TracChangeset for help on using the changeset viewer.