Changeset 10880


Ignore:
Timestamp:
2019-04-17T12:02:14+02:00 (19 months ago)
Author:
davestorkey
Message:

branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps:

  1. Move time indices from dom_oce.F90 to step.F90.
  2. Implement time dimension for passive tracers with independent set of time indices in trcstp.F90.
  3. Update all the traadv and trcadv modules.
Location:
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src
Files:
15 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DOM/dom_oce.F90

    r10876 r10880  
    3838   LOGICAL , PUBLIC ::   ln_iscpl       !: coupling with ice sheet 
    3939   LOGICAL , PUBLIC ::   ln_crs         !: Apply grid coarsening to dynamical model output or online passive tracers 
    40  
    41    !!---------------------------------------------------------------------- 
    42    !! time level indices 
    43    !!---------------------------------------------------------------------- 
    44    INTEGER, PUBLIC :: Nbb, Nnn, Naa, Nrhs 
    4540 
    4641   !! Free surface parameters 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv.F90

    r10874 r10880  
    7575CONTAINS 
    7676 
    77    SUBROUTINE tra_adv( kt ) 
     77   SUBROUTINE tra_adv( kt, Kbb, Kmm, pts, Krhs ) 
    7878      !!---------------------------------------------------------------------- 
    7979      !!                  ***  ROUTINE tra_adv  *** 
     
    8181      !! ** Purpose :   compute the ocean tracer advection trend. 
    8282      !! 
    83       !! ** Method  : - Update (ua,va) with the advection term following nadv 
    84       !!---------------------------------------------------------------------- 
    85       INTEGER, INTENT(in) ::   kt   ! ocean time-step index 
     83      !! ** Method  : - Update (uu(:,:,:,Krhs),vv(:,:,:,Krhs)) with the advection term following nadv 
     84      !!---------------------------------------------------------------------- 
     85      INTEGER                                  , INTENT(in)    :: kt             ! ocean time-step index 
     86      INTEGER                                  , INTENT(in)    :: Kbb, Kmm, Krhs ! time level indices 
     87      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts,jpt), INTENT(inout) :: pts            ! active tracers and RHS of tracer equation 
    8688      ! 
    8789      INTEGER ::   jk   ! dummy loop index 
    88       REAL(wp), DIMENSION(jpi,jpj,jpk)        :: zun, zvn, zwn   ! 3D workspace 
     90      REAL(wp), DIMENSION(jpi,jpj,jpk)        :: zuu, zvv, zww   ! 3D workspace 
    8991      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ztrdt, ztrds 
    9092      !!---------------------------------------------------------------------- 
     
    98100      ! 
    99101      !                                         !==  effective transport  ==! 
    100       zun(:,:,jpk) = 0._wp 
    101       zvn(:,:,jpk) = 0._wp 
    102       zwn(:,:,jpk) = 0._wp 
     102      zuu(:,:,jpk) = 0._wp 
     103      zvv(:,:,jpk) = 0._wp 
     104      zww(:,:,jpk) = 0._wp 
    103105      IF( ln_wave .AND. ln_sdw )  THEN 
    104106         DO jk = 1, jpkm1                                                       ! eulerian transport + Stokes Drift 
    105             zun(:,:,jk) = e2u  (:,:) * e3u_n(:,:,jk) * ( un(:,:,jk) + usd(:,:,jk) ) 
    106             zvn(:,:,jk) = e1v  (:,:) * e3v_n(:,:,jk) * ( vn(:,:,jk) + vsd(:,:,jk) ) 
    107             zwn(:,:,jk) = e1e2t(:,:)                 * ( wn(:,:,jk) + wsd(:,:,jk) ) 
     107            zuu(:,:,jk) = e2u  (:,:) * e3u(:,:,jk,Kmm) * ( uu(:,:,jk,Kmm) + usd(:,:,jk) ) 
     108            zvv(:,:,jk) = e1v  (:,:) * e3v(:,:,jk,Kmm) * ( vv(:,:,jk,Kmm) + vsd(:,:,jk) ) 
     109            zww(:,:,jk) = e1e2t(:,:)                 * ( ww(:,:,jk) + wsd(:,:,jk) ) 
    108110         END DO 
    109111      ELSE 
    110112         DO jk = 1, jpkm1 
    111             zun(:,:,jk) = e2u  (:,:) * e3u_n(:,:,jk) * un(:,:,jk)               ! eulerian transport only 
    112             zvn(:,:,jk) = e1v  (:,:) * e3v_n(:,:,jk) * vn(:,:,jk) 
    113             zwn(:,:,jk) = e1e2t(:,:)                 * wn(:,:,jk) 
     113            zuu(:,:,jk) = e2u  (:,:) * e3u(:,:,jk,Kmm) * uu(:,:,jk,Kmm)               ! eulerian transport only 
     114            zvv(:,:,jk) = e1v  (:,:) * e3v(:,:,jk,Kmm) * vv(:,:,jk,Kmm) 
     115            zww(:,:,jk) = e1e2t(:,:)                 * ww(:,:,jk) 
    114116         END DO 
    115117      ENDIF 
    116118      ! 
    117119      IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN                                ! add z-tilde and/or vvl corrections 
    118          zun(:,:,:) = zun(:,:,:) + un_td(:,:,:) 
    119          zvn(:,:,:) = zvn(:,:,:) + vn_td(:,:,:) 
    120       ENDIF 
    121       ! 
    122       zun(:,:,jpk) = 0._wp                                                      ! no transport trough the bottom 
    123       zvn(:,:,jpk) = 0._wp 
    124       zwn(:,:,jpk) = 0._wp 
     120         zuu(:,:,:) = zuu(:,:,:) + un_td(:,:,:) 
     121         zvv(:,:,:) = zvv(:,:,:) + vn_td(:,:,:) 
     122      ENDIF 
     123      ! 
     124      zuu(:,:,jpk) = 0._wp                                                      ! no transport trough the bottom 
     125      zvv(:,:,jpk) = 0._wp 
     126      zww(:,:,jpk) = 0._wp 
    125127      ! 
    126128      IF( ln_ldfeiv .AND. .NOT. ln_traldf_triad )   & 
    127          &              CALL ldf_eiv_trp( kt, nit000, zun, zvn, zwn, 'TRA' )   ! add the eiv transport (if necessary) 
    128       ! 
    129       IF( ln_mle    )   CALL tra_mle_trp( kt, nit000, zun, zvn, zwn, 'TRA' )   ! add the mle transport (if necessary) 
    130       ! 
    131       CALL iom_put( "uocetr_eff", zun )                                        ! output effective transport       
    132       CALL iom_put( "vocetr_eff", zvn ) 
    133       CALL iom_put( "wocetr_eff", zwn ) 
     129         &              CALL ldf_eiv_trp( kt, nit000, zuu, zvv, zww, 'TRA' )   ! add the eiv transport (if necessary) 
     130      ! 
     131      IF( ln_mle    )   CALL tra_mle_trp( kt, nit000, zuu, zvv, zww, 'TRA' )   ! add the mle transport (if necessary) 
     132      ! 
     133      CALL iom_put( "uocetr_eff", zuu )                                        ! output effective transport       
     134      CALL iom_put( "vocetr_eff", zvv ) 
     135      CALL iom_put( "wocetr_eff", zww ) 
    134136      ! 
    135137!!gm ??? 
    136       IF( ln_diaptr )   CALL dia_ptr( zvn )                                    ! diagnose the effective MSF  
     138      IF( ln_diaptr )   CALL dia_ptr( zvv )                                    ! diagnose the effective MSF  
    137139!!gm ??? 
    138140      ! 
    139141      IF( l_trdtra )   THEN                    !* Save ta and sa trends 
    140142         ALLOCATE( ztrdt(jpi,jpj,jpk), ztrds(jpi,jpj,jpk) ) 
    141          ztrdt(:,:,:) = tsa(:,:,:,jp_tem) 
    142          ztrds(:,:,:) = tsa(:,:,:,jp_sal) 
     143         ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) 
     144         ztrds(:,:,:) = pts(:,:,:,jp_sal,Krhs) 
    143145      ENDIF 
    144146      ! 
     
    146148      ! 
    147149      CASE ( np_CEN )                                 ! Centered scheme : 2nd / 4th order 
    148          CALL tra_adv_cen    ( kt, nit000, 'TRA',         zun, zvn, zwn     , tsn, tsa, jpts, nn_cen_h, nn_cen_v ) 
     150         CALL tra_adv_cen    ( kt, nit000, 'TRA',         zuu, zvv, zww, Kmm, pts, jpts, Krhs, nn_cen_h, nn_cen_v ) 
    149151      CASE ( np_FCT )                                 ! FCT scheme      : 2nd / 4th order 
    150          CALL tra_adv_fct    ( kt, nit000, 'TRA', r2dt, zun, zvn, zwn, tsb, tsn, tsa, jpts, nn_fct_h, nn_fct_v ) 
     152         CALL tra_adv_fct    ( kt, nit000, 'TRA', r2dt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, nn_fct_h, nn_fct_v ) 
    151153      CASE ( np_MUS )                                 ! MUSCL 
    152          CALL tra_adv_mus    ( kt, nit000, 'TRA', r2dt, zun, zvn, zwn, tsb,      tsa, jpts        , ln_mus_ups )  
     154         CALL tra_adv_mus    ( kt, nit000, 'TRA', r2dt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, ln_mus_ups )  
    153155      CASE ( np_UBS )                                 ! UBS 
    154          CALL tra_adv_ubs    ( kt, nit000, 'TRA', r2dt, zun, zvn, zwn, tsb, tsn, tsa, jpts        , nn_ubs_v   ) 
     156         CALL tra_adv_ubs    ( kt, nit000, 'TRA', r2dt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs, nn_ubs_v   ) 
    155157      CASE ( np_QCK )                                 ! QUICKEST 
    156          CALL tra_adv_qck    ( kt, nit000, 'TRA', r2dt, zun, zvn, zwn, tsb, tsn, tsa, jpts                    ) 
     158         CALL tra_adv_qck    ( kt, nit000, 'TRA', r2dt, zuu, zvv, zww, Kbb, Kmm, pts, jpts, Krhs ) 
    157159      ! 
    158160      END SELECT 
     
    160162      IF( l_trdtra )   THEN                      ! save the advective trends for further diagnostics 
    161163         DO jk = 1, jpkm1 
    162             ztrdt(:,:,jk) = tsa(:,:,jk,jp_tem) - ztrdt(:,:,jk) 
    163             ztrds(:,:,jk) = tsa(:,:,jk,jp_sal) - ztrds(:,:,jk) 
     164            ztrdt(:,:,jk) = pts(:,:,jk,jp_tem,Krhs) - ztrdt(:,:,jk) 
     165            ztrds(:,:,jk) = pts(:,:,jk,jp_sal,Krhs) - ztrds(:,:,jk) 
    164166         END DO 
    165167         CALL trd_tra( kt, 'TRA', jp_tem, jptra_totad, ztrdt ) 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_cen.F90

    r10874 r10880  
    4444CONTAINS 
    4545 
    46    SUBROUTINE tra_adv_cen( kt, kit000, cdtype, pun, pvn, pwn,     & 
    47       &                                        ptn, pta, kjpt, kn_cen_h, kn_cen_v )  
     46   SUBROUTINE tra_adv_cen( kt, kit000, cdtype, pu_mm, pv_mm, pww,     & 
     47      &                    Kmm, pt, kjpt, Krhs, kn_cen_h, kn_cen_v )  
    4848      !!---------------------------------------------------------------------- 
    4949      !!                  ***  ROUTINE tra_adv_cen  *** 
     
    5959      !!                = 4  ==>> 4th order COMPACT  scheme     -      - 
    6060      !! 
    61       !! ** Action : - update pta  with the now advective tracer trends 
     61      !! ** Action : - update pt(:,:,:,:,Krhs)  with the now advective tracer trends 
    6262      !!             - send trends to trdtra module for further diagnostcs (l_trdtra=T) 
    6363      !!             - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 
    6464      !!---------------------------------------------------------------------- 
    65       INTEGER                              , INTENT(in   ) ::   kt              ! ocean time-step index 
    66       INTEGER                              , INTENT(in   ) ::   kit000          ! first time step index 
    67       CHARACTER(len=3)                     , INTENT(in   ) ::   cdtype          ! =TRA or TRC (tracer indicator) 
    68       INTEGER                              , INTENT(in   ) ::   kjpt            ! number of tracers 
    69       INTEGER                              , INTENT(in   ) ::   kn_cen_h        ! =2/4 (2nd or 4th order scheme) 
    70       INTEGER                              , INTENT(in   ) ::   kn_cen_v        ! =2/4 (2nd or 4th order scheme) 
    71       REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ) ::   pun, pvn, pwn   ! 3 ocean velocity components 
    72       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptn             ! now tracer fields 
    73       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta             ! tracer trend  
     65      INTEGER                                  , INTENT(in   ) ::   kt              ! ocean time-step index 
     66      INTEGER                                  , INTENT(in   ) ::   Kmm, Krhs       ! ocean time level indices 
     67      INTEGER                                  , INTENT(in   ) ::   kit000          ! first time step index 
     68      CHARACTER(len=3)                         , INTENT(in   ) ::   cdtype          ! =TRA or TRC (tracer indicator) 
     69      INTEGER                                  , INTENT(in   ) ::   kjpt            ! number of tracers 
     70      INTEGER                                  , INTENT(in   ) ::   kn_cen_h        ! =2/4 (2nd or 4th order scheme) 
     71      INTEGER                                  , INTENT(in   ) ::   kn_cen_v        ! =2/4 (2nd or 4th order scheme) 
     72      REAL(wp), DIMENSION(jpi,jpj,jpk         ), INTENT(in   ) ::   pu_mm, pv_mm, pww   ! 3 ocean velocity components 
     73      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) ::   pt              ! active tracers and RHS of tracer equation 
    7474      ! 
    7575      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices 
     
    106106               DO jj = 1, jpjm1 
    107107                  DO ji = 1, fs_jpim1   ! vector opt. 
    108                      zwx(ji,jj,jk) = 0.5_wp * pun(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji+1,jj  ,jk,jn) ) 
    109                      zwy(ji,jj,jk) = 0.5_wp * pvn(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji  ,jj+1,jk,jn) ) 
     108                     zwx(ji,jj,jk) = 0.5_wp * pu_mm(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj  ,jk,jn,Kmm) ) 
     109                     zwy(ji,jj,jk) = 0.5_wp * pv_mm(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji  ,jj+1,jk,jn,Kmm) ) 
    110110                  END DO 
    111111               END DO 
     
    118118               DO jj = 2, jpjm1 
    119119                  DO ji = fs_2, fs_jpim1   ! vector opt. 
    120                      ztu(ji,jj,jk) = ( ptn(ji+1,jj  ,jk,jn) - ptn(ji,jj,jk,jn) ) * umask(ji,jj,jk) 
    121                      ztv(ji,jj,jk) = ( ptn(ji  ,jj+1,jk,jn) - ptn(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 
     120                     ztu(ji,jj,jk) = ( pt(ji+1,jj  ,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * umask(ji,jj,jk) 
     121                     ztv(ji,jj,jk) = ( pt(ji  ,jj+1,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * vmask(ji,jj,jk) 
    122122                  END DO 
    123123               END DO 
     
    128128               DO jj = 2, jpjm1 
    129129                  DO ji = 1, fs_jpim1   ! vector opt. 
    130                      zC2t_u = ptn(ji,jj,jk,jn) + ptn(ji+1,jj  ,jk,jn)   ! C2 interpolation of T at u- & v-points (x2) 
    131                      zC2t_v = ptn(ji,jj,jk,jn) + ptn(ji  ,jj+1,jk,jn) 
     130                     zC2t_u = pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj  ,jk,jn,Kmm)   ! C2 interpolation of T at u- & v-points (x2) 
     131                     zC2t_v = pt(ji,jj,jk,jn,Kmm) + pt(ji  ,jj+1,jk,jn,Kmm) 
    132132                     !                                                  ! C4 interpolation of T at u- & v-points (x2) 
    133133                     zC4t_u =  zC2t_u + r1_6 * ( ztu(ji-1,jj,jk) - ztu(ji+1,jj,jk) ) 
    134134                     zC4t_v =  zC2t_v + r1_6 * ( ztv(ji,jj-1,jk) - ztv(ji,jj+1,jk) ) 
    135135                     !                                                  ! C4 fluxes 
    136                      zwx(ji,jj,jk) =  0.5_wp * pun(ji,jj,jk) * zC4t_u 
    137                      zwy(ji,jj,jk) =  0.5_wp * pvn(ji,jj,jk) * zC4t_v 
     136                     zwx(ji,jj,jk) =  0.5_wp * pu_mm(ji,jj,jk) * zC4t_u 
     137                     zwy(ji,jj,jk) =  0.5_wp * pv_mm(ji,jj,jk) * zC4t_v 
    138138                  END DO 
    139139               END DO 
     
    150150               DO jj = 2, jpjm1 
    151151                  DO ji = fs_2, fs_jpim1   ! vector opt. 
    152                      zwz(ji,jj,jk) = 0.5 * pwn(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji,jj,jk-1,jn) ) * wmask(ji,jj,jk) 
     152                     zwz(ji,jj,jk) = 0.5 * pww(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji,jj,jk-1,jn,Kmm) ) * wmask(ji,jj,jk) 
    153153                  END DO 
    154154               END DO 
     
    156156            ! 
    157157         CASE(  4  )                         !* 4th order compact 
    158             CALL interp_4th_cpt( ptn(:,:,:,jn) , ztw )      ! ztw = interpolated value of T at w-point 
     158            CALL interp_4th_cpt( pt(:,:,:,jn,Kmm) , ztw )      ! ztw = interpolated value of T at w-point 
    159159            DO jk = 2, jpkm1 
    160160               DO jj = 2, jpjm1 
    161161                  DO ji = fs_2, fs_jpim1 
    162                      zwz(ji,jj,jk) = pwn(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk) 
     162                     zwz(ji,jj,jk) = pww(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk) 
    163163                  END DO 
    164164               END DO 
     
    171171               DO jj = 1, jpj 
    172172                  DO ji = 1, jpi 
    173                      zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * ptn(ji,jj,mikt(ji,jj),jn)  
     173                     zwz(ji,jj, mikt(ji,jj) ) = pww(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kmm)  
    174174                  END DO 
    175175               END DO    
    176176            ELSE                                   ! no ice-shelf cavities (only ocean surface) 
    177                zwz(:,:,1) = pwn(:,:,1) * ptn(:,:,1,jn) 
     177               zwz(:,:,1) = pww(:,:,1) * pt(:,:,1,jn,Kmm) 
    178178            ENDIF 
    179179         ENDIF 
     
    182182            DO jj = 2, jpjm1 
    183183               DO ji = fs_2, fs_jpim1   ! vector opt. 
    184                   pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn)    & 
     184                  pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs)    & 
    185185                     &             - (  zwx(ji,jj,jk) - zwx(ji-1,jj  ,jk  )    & 
    186186                     &                + zwy(ji,jj,jk) - zwy(ji  ,jj-1,jk  )    & 
    187                      &                + zwz(ji,jj,jk) - zwz(ji  ,jj  ,jk+1)  ) * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 
     187                     &                + zwz(ji,jj,jk) - zwz(ji  ,jj  ,jk+1)  ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    188188               END DO 
    189189            END DO 
     
    191191         !                             ! trend diagnostics 
    192192         IF( l_trd ) THEN 
    193             CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pun, ptn(:,:,:,jn) ) 
    194             CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pvn, ptn(:,:,:,jn) ) 
    195             CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pwn, ptn(:,:,:,jn) ) 
     193            CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pu_mm, pt(:,:,:,jn,Kmm) ) 
     194            CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pv_mm, pt(:,:,:,jn,Kmm) ) 
     195            CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pww, pt(:,:,:,jn,Kmm) ) 
    196196         END IF 
    197197         !                                 ! "Poleward" heat and salt transports  
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_fct.F90

    r10874 r10880  
    5252CONTAINS 
    5353 
    54    SUBROUTINE tra_adv_fct( kt, kit000, cdtype, p2dt, pun, pvn, pwn,       & 
    55       &                                              ptb, ptn, pta, kjpt, kn_fct_h, kn_fct_v ) 
     54   SUBROUTINE tra_adv_fct( kt, kit000, cdtype, p2dt, pu_mm, pv_mm, pww,       & 
     55      &                    Kbb, Kmm, pt, kjpt, Krhs, kn_fct_h, kn_fct_v ) 
    5656      !!---------------------------------------------------------------------- 
    5757      !!                  ***  ROUTINE tra_adv_fct  *** 
     
    6565      !!               - corrected flux (monotonic correction)  
    6666      !! 
    67       !! ** Action : - update pta  with the now advective tracer trends 
     67      !! ** Action : - update pt(:,:,:,:,Krhs)  with the now advective tracer trends 
    6868      !!             - send trends to trdtra module for further diagnostics (l_trdtra=T) 
    6969      !!             - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 
    7070      !!---------------------------------------------------------------------- 
    71       INTEGER                              , INTENT(in   ) ::   kt              ! ocean time-step index 
    72       INTEGER                              , INTENT(in   ) ::   kit000          ! first time step index 
    73       CHARACTER(len=3)                     , INTENT(in   ) ::   cdtype          ! =TRA or TRC (tracer indicator) 
    74       INTEGER                              , INTENT(in   ) ::   kjpt            ! number of tracers 
    75       INTEGER                              , INTENT(in   ) ::   kn_fct_h        ! order of the FCT scheme (=2 or 4) 
    76       INTEGER                              , INTENT(in   ) ::   kn_fct_v        ! order of the FCT scheme (=2 or 4) 
    77       REAL(wp)                             , INTENT(in   ) ::   p2dt            ! tracer time-step 
    78       REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ) ::   pun, pvn, pwn   ! 3 ocean velocity components 
    79       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb, ptn        ! before and now tracer fields 
    80       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta             ! tracer trend  
     71      INTEGER                                  , INTENT(in   ) ::   kt              ! ocean time-step index 
     72      INTEGER                                  , INTENT(in   ) ::   Kbb, Kmm, Krhs  ! ocean time level indices 
     73      INTEGER                                  , INTENT(in   ) ::   kit000          ! first time step index 
     74      CHARACTER(len=3)                         , INTENT(in   ) ::   cdtype          ! =TRA or TRC (tracer indicator) 
     75      INTEGER                                  , INTENT(in   ) ::   kjpt            ! number of tracers 
     76      INTEGER                                  , INTENT(in   ) ::   kn_fct_h        ! order of the FCT scheme (=2 or 4) 
     77      INTEGER                                  , INTENT(in   ) ::   kn_fct_v        ! order of the FCT scheme (=2 or 4) 
     78      REAL(wp)                                 , INTENT(in   ) ::   p2dt            ! tracer time-step 
     79      REAL(wp), DIMENSION(jpi,jpj,jpk         ), INTENT(in   ) ::   pu_mm, pv_mm, pww   ! 3 ocean velocity components 
     80      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) ::   pt              ! active tracers and RHS of tracer equation 
    8181      ! 
    8282      INTEGER  ::   ji, jj, jk, jn                           ! dummy loop indices   
     
    125125               DO ji = 1, fs_jpim1   ! vector opt. 
    126126                  ! upstream scheme 
    127                   zfp_ui = pun(ji,jj,jk) + ABS( pun(ji,jj,jk) ) 
    128                   zfm_ui = pun(ji,jj,jk) - ABS( pun(ji,jj,jk) ) 
    129                   zfp_vj = pvn(ji,jj,jk) + ABS( pvn(ji,jj,jk) ) 
    130                   zfm_vj = pvn(ji,jj,jk) - ABS( pvn(ji,jj,jk) ) 
    131                   zwx(ji,jj,jk) = 0.5 * ( zfp_ui * ptb(ji,jj,jk,jn) + zfm_ui * ptb(ji+1,jj  ,jk,jn) ) 
    132                   zwy(ji,jj,jk) = 0.5 * ( zfp_vj * ptb(ji,jj,jk,jn) + zfm_vj * ptb(ji  ,jj+1,jk,jn) ) 
     127                  zfp_ui = pu_mm(ji,jj,jk) + ABS( pu_mm(ji,jj,jk) ) 
     128                  zfm_ui = pu_mm(ji,jj,jk) - ABS( pu_mm(ji,jj,jk) ) 
     129                  zfp_vj = pv_mm(ji,jj,jk) + ABS( pv_mm(ji,jj,jk) ) 
     130                  zfm_vj = pv_mm(ji,jj,jk) - ABS( pv_mm(ji,jj,jk) ) 
     131                  zwx(ji,jj,jk) = 0.5 * ( zfp_ui * pt(ji,jj,jk,jn,Kbb) + zfm_ui * pt(ji+1,jj  ,jk,jn,Kbb) ) 
     132                  zwy(ji,jj,jk) = 0.5 * ( zfp_vj * pt(ji,jj,jk,jn,Kbb) + zfm_vj * pt(ji  ,jj+1,jk,jn,Kbb) ) 
    133133               END DO 
    134134            END DO 
     
    138138            DO jj = 1, jpj 
    139139               DO ji = 1, jpi 
    140                   zfp_wk = pwn(ji,jj,jk) + ABS( pwn(ji,jj,jk) ) 
    141                   zfm_wk = pwn(ji,jj,jk) - ABS( pwn(ji,jj,jk) ) 
    142                   zwz(ji,jj,jk) = 0.5 * ( zfp_wk * ptb(ji,jj,jk,jn) + zfm_wk * ptb(ji,jj,jk-1,jn) ) * wmask(ji,jj,jk) 
     140                  zfp_wk = pww(ji,jj,jk) + ABS( pww(ji,jj,jk) ) 
     141                  zfm_wk = pww(ji,jj,jk) - ABS( pww(ji,jj,jk) ) 
     142                  zwz(ji,jj,jk) = 0.5 * ( zfp_wk * pt(ji,jj,jk,jn,Kbb) + zfm_wk * pt(ji,jj,jk-1,jn,Kbb) ) * wmask(ji,jj,jk) 
    143143               END DO 
    144144            END DO 
     
    148148               DO jj = 1, jpj 
    149149                  DO ji = 1, jpi 
    150                      zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * ptb(ji,jj,mikt(ji,jj),jn)   ! linear free surface  
     150                     zwz(ji,jj, mikt(ji,jj) ) = pww(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kbb)   ! linear free surface  
    151151                  END DO 
    152152               END DO    
    153153            ELSE                             ! no cavities: only at the ocean surface 
    154                zwz(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn) 
     154               zwz(:,:,1) = pww(:,:,1) * pt(:,:,1,jn,Kbb) 
    155155            ENDIF 
    156156         ENDIF 
     
    164164                     &      + zwz(ji,jj,jk) - zwz(ji  ,jj  ,jk+1) ) * r1_e1e2t(ji,jj) 
    165165                  !                             ! update and guess with monotonic sheme 
    166                   pta(ji,jj,jk,jn) =                     pta(ji,jj,jk,jn) +        ztra   / e3t_n(ji,jj,jk) * tmask(ji,jj,jk) 
    167                   zwi(ji,jj,jk)    = ( e3t_b(ji,jj,jk) * ptb(ji,jj,jk,jn) + p2dt * ztra ) / e3t_a(ji,jj,jk) * tmask(ji,jj,jk) 
     166                  pt(ji,jj,jk,jn,Krhs) =                     pt(ji,jj,jk,jn,Krhs) +        ztra   / e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk) 
     167                  zwi(ji,jj,jk)    = ( e3t(ji,jj,jk,Kbb) * pt(ji,jj,jk,jn,Kbb) + p2dt * ztra ) / e3t(ji,jj,jk,Krhs) * tmask(ji,jj,jk) 
    168168               END DO 
    169169            END DO 
     
    184184               DO jj = 1, jpjm1 
    185185                  DO ji = 1, fs_jpim1   ! vector opt. 
    186                      zwx(ji,jj,jk) = 0.5_wp * pun(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji+1,jj,jk,jn) ) - zwx(ji,jj,jk) 
    187                      zwy(ji,jj,jk) = 0.5_wp * pvn(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji,jj+1,jk,jn) ) - zwy(ji,jj,jk) 
     186                     zwx(ji,jj,jk) = 0.5_wp * pu_mm(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj,jk,jn,Kmm) ) - zwx(ji,jj,jk) 
     187                     zwy(ji,jj,jk) = 0.5_wp * pv_mm(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji,jj+1,jk,jn,Kmm) ) - zwy(ji,jj,jk) 
    188188                  END DO 
    189189               END DO 
     
    196196               DO jj = 1, jpjm1                    ! 1st derivative (gradient) 
    197197                  DO ji = 1, fs_jpim1   ! vector opt. 
    198                      ztu(ji,jj,jk) = ( ptn(ji+1,jj  ,jk,jn) - ptn(ji,jj,jk,jn) ) * umask(ji,jj,jk) 
    199                      ztv(ji,jj,jk) = ( ptn(ji  ,jj+1,jk,jn) - ptn(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 
     198                     ztu(ji,jj,jk) = ( pt(ji+1,jj  ,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * umask(ji,jj,jk) 
     199                     ztv(ji,jj,jk) = ( pt(ji  ,jj+1,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * vmask(ji,jj,jk) 
    200200                  END DO 
    201201               END DO 
     
    212212               DO jj = 1, jpjm1 
    213213                  DO ji = 1, fs_jpim1   ! vector opt. 
    214                      zC2t_u = ptn(ji,jj,jk,jn) + ptn(ji+1,jj  ,jk,jn)   ! 2 x C2 interpolation of T at u- & v-points 
    215                      zC2t_v = ptn(ji,jj,jk,jn) + ptn(ji  ,jj+1,jk,jn) 
     214                     zC2t_u = pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj  ,jk,jn,Kmm)   ! 2 x C2 interpolation of T at u- & v-points 
     215                     zC2t_v = pt(ji,jj,jk,jn,Kmm) + pt(ji  ,jj+1,jk,jn,Kmm) 
    216216                     !                                                  ! C4 minus upstream advective fluxes  
    217                      zwx(ji,jj,jk) =  0.5_wp * pun(ji,jj,jk) * ( zC2t_u + zltu(ji,jj,jk) - zltu(ji+1,jj,jk) ) - zwx(ji,jj,jk) 
    218                      zwy(ji,jj,jk) =  0.5_wp * pvn(ji,jj,jk) * ( zC2t_v + zltv(ji,jj,jk) - zltv(ji,jj+1,jk) ) - zwy(ji,jj,jk) 
     217                     zwx(ji,jj,jk) =  0.5_wp * pu_mm(ji,jj,jk) * ( zC2t_u + zltu(ji,jj,jk) - zltu(ji+1,jj,jk) ) - zwx(ji,jj,jk) 
     218                     zwy(ji,jj,jk) =  0.5_wp * pv_mm(ji,jj,jk) * ( zC2t_v + zltv(ji,jj,jk) - zltv(ji,jj+1,jk) ) - zwy(ji,jj,jk) 
    219219                  END DO 
    220220               END DO 
     
    227227               DO jj = 1, jpjm1 
    228228                  DO ji = 1, fs_jpim1   ! vector opt. 
    229                      ztu(ji,jj,jk) = ( ptn(ji+1,jj  ,jk,jn) - ptn(ji,jj,jk,jn) ) * umask(ji,jj,jk) 
    230                      ztv(ji,jj,jk) = ( ptn(ji  ,jj+1,jk,jn) - ptn(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 
     229                     ztu(ji,jj,jk) = ( pt(ji+1,jj  ,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * umask(ji,jj,jk) 
     230                     ztv(ji,jj,jk) = ( pt(ji  ,jj+1,jk,jn,Kmm) - pt(ji,jj,jk,jn,Kmm) ) * vmask(ji,jj,jk) 
    231231                  END DO 
    232232               END DO 
     
    237237               DO jj = 2, jpjm1 
    238238                  DO ji = 2, fs_jpim1   ! vector opt. 
    239                      zC2t_u = ptn(ji,jj,jk,jn) + ptn(ji+1,jj  ,jk,jn)   ! 2 x C2 interpolation of T at u- & v-points (x2) 
    240                      zC2t_v = ptn(ji,jj,jk,jn) + ptn(ji  ,jj+1,jk,jn) 
     239                     zC2t_u = pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj  ,jk,jn,Kmm)   ! 2 x C2 interpolation of T at u- & v-points (x2) 
     240                     zC2t_v = pt(ji,jj,jk,jn,Kmm) + pt(ji  ,jj+1,jk,jn,Kmm) 
    241241                     !                                                  ! C4 interpolation of T at u- & v-points (x2) 
    242242                     zC4t_u =  zC2t_u + r1_6 * ( ztu(ji-1,jj  ,jk) - ztu(ji+1,jj  ,jk) ) 
    243243                     zC4t_v =  zC2t_v + r1_6 * ( ztv(ji  ,jj-1,jk) - ztv(ji  ,jj+1,jk) ) 
    244244                     !                                                  ! C4 minus upstream advective fluxes  
    245                      zwx(ji,jj,jk) =  0.5_wp * pun(ji,jj,jk) * zC4t_u - zwx(ji,jj,jk) 
    246                      zwy(ji,jj,jk) =  0.5_wp * pvn(ji,jj,jk) * zC4t_v - zwy(ji,jj,jk) 
     245                     zwx(ji,jj,jk) =  0.5_wp * pu_mm(ji,jj,jk) * zC4t_u - zwx(ji,jj,jk) 
     246                     zwy(ji,jj,jk) =  0.5_wp * pv_mm(ji,jj,jk) * zC4t_v - zwy(ji,jj,jk) 
    247247                  END DO 
    248248               END DO 
     
    257257               DO jj = 2, jpjm1 
    258258                  DO ji = fs_2, fs_jpim1 
    259                      zwz(ji,jj,jk) =  (  pwn(ji,jj,jk) * 0.5_wp * ( ptn(ji,jj,jk,jn) + ptn(ji,jj,jk-1,jn) )   & 
     259                     zwz(ji,jj,jk) =  (  pww(ji,jj,jk) * 0.5_wp * ( pt(ji,jj,jk,jn,Kmm) + pt(ji,jj,jk-1,jn,Kmm) )   & 
    260260                        &              - zwz(ji,jj,jk)  ) * wmask(ji,jj,jk) 
    261261                  END DO 
     
    264264            ! 
    265265         CASE(  4  )                   !- 4th order COMPACT 
    266             CALL interp_4th_cpt( ptn(:,:,:,jn) , ztw )   ! zwt = COMPACT interpolation of T at w-point 
     266            CALL interp_4th_cpt( pt(:,:,:,jn,Kmm) , ztw )   ! zwt = COMPACT interpolation of T at w-point 
    267267            DO jk = 2, jpkm1 
    268268               DO jj = 2, jpjm1 
    269269                  DO ji = fs_2, fs_jpim1 
    270                      zwz(ji,jj,jk) = ( pwn(ji,jj,jk) * ztw(ji,jj,jk) - zwz(ji,jj,jk) ) * wmask(ji,jj,jk) 
     270                     zwz(ji,jj,jk) = ( pww(ji,jj,jk) * ztw(ji,jj,jk) - zwz(ji,jj,jk) ) * wmask(ji,jj,jk) 
    271271                  END DO 
    272272               END DO 
     
    282282         !        !==  monotonicity algorithm  ==! 
    283283         ! 
    284          CALL nonosc( ptb(:,:,:,jn), zwx, zwy, zwz, zwi, p2dt ) 
     284         CALL nonosc( Kmm, pt(:,:,:,jn,Kbb), zwx, zwy, zwz, zwi, p2dt ) 
    285285         ! 
    286286         !        !==  final trend with corrected fluxes  ==! 
     
    289289            DO jj = 2, jpjm1 
    290290               DO ji = fs_2, fs_jpim1   ! vector opt.   
    291                   pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - (  zwx(ji,jj,jk) - zwx(ji-1,jj  ,jk  )   & 
     291                  pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) - (  zwx(ji,jj,jk) - zwx(ji-1,jj  ,jk  )   & 
    292292                     &                                   + zwy(ji,jj,jk) - zwy(ji  ,jj-1,jk  )   & 
    293293                     &                                   + zwz(ji,jj,jk) - zwz(ji  ,jj  ,jk+1) ) & 
    294                      &                                * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 
     294                     &                                * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    295295               END DO 
    296296            END DO 
     
    303303            ! 
    304304            IF( l_trd ) THEN              ! trend diagnostics 
    305                CALL trd_tra( kt, cdtype, jn, jptra_xad, ztrdx, pun, ptn(:,:,:,jn) ) 
    306                CALL trd_tra( kt, cdtype, jn, jptra_yad, ztrdy, pvn, ptn(:,:,:,jn) ) 
    307                CALL trd_tra( kt, cdtype, jn, jptra_zad, ztrdz, pwn, ptn(:,:,:,jn) ) 
     305               CALL trd_tra( kt, cdtype, jn, jptra_xad, ztrdx, pu_mm, pt(:,:,:,jn,Kmm) ) 
     306               CALL trd_tra( kt, cdtype, jn, jptra_yad, ztrdy, pv_mm, pt(:,:,:,jn,Kmm) ) 
     307               CALL trd_tra( kt, cdtype, jn, jptra_zad, ztrdz, pww, pt(:,:,:,jn,Kmm) ) 
    308308            ENDIF 
    309309            !                             ! heat/salt transport 
     
    328328 
    329329 
    330    SUBROUTINE nonosc( pbef, paa, pbb, pcc, paft, p2dt ) 
     330   SUBROUTINE nonosc( Kmm, pbef, paa, pbb, pcc, paft, p2dt ) 
    331331      !!--------------------------------------------------------------------- 
    332332      !!                    ***  ROUTINE nonosc  *** 
     
    341341      !!       in-space based differencing for fluid 
    342342      !!---------------------------------------------------------------------- 
     343      INTEGER                          , INTENT(in   ) ::   Kmm             ! time level index  
    343344      REAL(wp)                         , INTENT(in   ) ::   p2dt            ! tracer time-step 
    344345      REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in   ) ::   pbef, paft      ! before & after field 
     
    392393 
    393394               ! up & down beta terms 
    394                zbt = e1e2t(ji,jj) * e3t_n(ji,jj,jk) / p2dt 
     395               zbt = e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) / p2dt 
    395396               zbetup(ji,jj,jk) = ( zup            - paft(ji,jj,jk) ) / ( zpos + zrtrn ) * zbt 
    396397               zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zdo            ) / ( zneg + zrtrn ) * zbt 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_mus.F90

    r10874 r10880  
    5454CONTAINS 
    5555 
    56    SUBROUTINE tra_adv_mus( kt, kit000, cdtype, p2dt, pun, pvn, pwn,             & 
    57       &                                              ptb, pta, kjpt, ld_msc_ups ) 
     56   SUBROUTINE tra_adv_mus( kt, kit000, cdtype, p2dt, pu_mm, pv_mm, pww,             & 
     57      &                    Kbb, Kmm, pt, kjpt, Krhs, ld_msc_ups ) 
    5858      !!---------------------------------------------------------------------- 
    5959      !!                    ***  ROUTINE tra_adv_mus  *** 
     
    6666      !!              ld_msc_ups=T :  
    6767      !! 
    68       !! ** Action : - update pta  with the now advective tracer trends 
     68      !! ** Action : - update pt(:,:,:,:,Krhs)  with the now advective tracer trends 
    6969      !!             - send trends to trdtra module for further diagnostcs (l_trdtra=T) 
    7070      !!             - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 
     
    7373      !!              IPSL, Sept. 2000 (http://www.lodyc.jussieu.fr/opa) 
    7474      !!---------------------------------------------------------------------- 
    75       INTEGER                              , INTENT(in   ) ::   kt              ! ocean time-step index 
    76       INTEGER                              , INTENT(in   ) ::   kit000          ! first time step index 
    77       CHARACTER(len=3)                     , INTENT(in   ) ::   cdtype          ! =TRA or TRC (tracer indicator) 
    78       INTEGER                              , INTENT(in   ) ::   kjpt            ! number of tracers 
    79       LOGICAL                              , INTENT(in   ) ::   ld_msc_ups      ! use upstream scheme within muscl 
    80       REAL(wp)                             , INTENT(in   ) ::   p2dt            ! tracer time-step 
    81       REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ) ::   pun, pvn, pwn   ! 3 ocean velocity components 
    82       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb             ! before tracer field 
    83       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta             ! tracer trend  
     75      INTEGER                                  , INTENT(in   ) ::   kt              ! ocean time-step index 
     76      INTEGER                                  , INTENT(in   ) ::   Kbb, Kmm, Krhs  ! ocean time level indices 
     77      INTEGER                                  , INTENT(in   ) ::   kit000          ! first time step index 
     78      CHARACTER(len=3)                         , INTENT(in   ) ::   cdtype          ! =TRA or TRC (tracer indicator) 
     79      INTEGER                                  , INTENT(in   ) ::   kjpt            ! number of tracers 
     80      LOGICAL                                  , INTENT(in   ) ::   ld_msc_ups      ! use upstream scheme within muscl 
     81      REAL(wp)                                 , INTENT(in   ) ::   p2dt            ! tracer time-step 
     82      REAL(wp), DIMENSION(jpi,jpj,jpk         ), INTENT(in   ) ::   pu_mm, pv_mm, pww   ! 3 ocean velocity components 
     83      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) ::   pt              ! active tracers and RHS of tracer equation 
    8484      ! 
    8585      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices 
     
    134134            DO jj = 1, jpjm1       
    135135               DO ji = 1, fs_jpim1   ! vector opt. 
    136                   zwx(ji,jj,jk) = umask(ji,jj,jk) * ( ptb(ji+1,jj,jk,jn) - ptb(ji,jj,jk,jn) ) 
    137                   zwy(ji,jj,jk) = vmask(ji,jj,jk) * ( ptb(ji,jj+1,jk,jn) - ptb(ji,jj,jk,jn) ) 
     136                  zwx(ji,jj,jk) = umask(ji,jj,jk) * ( pt(ji+1,jj,jk,jn,Kbb) - pt(ji,jj,jk,jn,Kbb) ) 
     137                  zwy(ji,jj,jk) = vmask(ji,jj,jk) * ( pt(ji,jj+1,jk,jn,Kbb) - pt(ji,jj,jk,jn,Kbb) ) 
    138138               END DO 
    139139           END DO 
     
    172172               DO ji = fs_2, fs_jpim1   ! vector opt. 
    173173                  ! MUSCL fluxes 
    174                   z0u = SIGN( 0.5, pun(ji,jj,jk) ) 
     174                  z0u = SIGN( 0.5, pu_mm(ji,jj,jk) ) 
    175175                  zalpha = 0.5 - z0u 
    176                   zu  = z0u - 0.5 * pun(ji,jj,jk) * p2dt * r1_e1e2u(ji,jj) / e3u_n(ji,jj,jk) 
    177                   zzwx = ptb(ji+1,jj,jk,jn) + xind(ji,jj,jk) * zu * zslpx(ji+1,jj,jk) 
    178                   zzwy = ptb(ji  ,jj,jk,jn) + xind(ji,jj,jk) * zu * zslpx(ji  ,jj,jk) 
    179                   zwx(ji,jj,jk) = pun(ji,jj,jk) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) 
     176                  zu  = z0u - 0.5 * pu_mm(ji,jj,jk) * p2dt * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,Kmm) 
     177                  zzwx = pt(ji+1,jj,jk,jn,Kbb) + xind(ji,jj,jk) * zu * zslpx(ji+1,jj,jk) 
     178                  zzwy = pt(ji  ,jj,jk,jn,Kbb) + xind(ji,jj,jk) * zu * zslpx(ji  ,jj,jk) 
     179                  zwx(ji,jj,jk) = pu_mm(ji,jj,jk) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) 
    180180                  ! 
    181                   z0v = SIGN( 0.5, pvn(ji,jj,jk) ) 
     181                  z0v = SIGN( 0.5, pv_mm(ji,jj,jk) ) 
    182182                  zalpha = 0.5 - z0v 
    183                   zv  = z0v - 0.5 * pvn(ji,jj,jk) * p2dt * r1_e1e2v(ji,jj) / e3v_n(ji,jj,jk) 
    184                   zzwx = ptb(ji,jj+1,jk,jn) + xind(ji,jj,jk) * zv * zslpy(ji,jj+1,jk) 
    185                   zzwy = ptb(ji,jj  ,jk,jn) + xind(ji,jj,jk) * zv * zslpy(ji,jj  ,jk) 
    186                   zwy(ji,jj,jk) = pvn(ji,jj,jk) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) 
     183                  zv  = z0v - 0.5 * pv_mm(ji,jj,jk) * p2dt * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,Kmm) 
     184                  zzwx = pt(ji,jj+1,jk,jn,Kbb) + xind(ji,jj,jk) * zv * zslpy(ji,jj+1,jk) 
     185                  zzwy = pt(ji,jj  ,jk,jn,Kbb) + xind(ji,jj,jk) * zv * zslpy(ji,jj  ,jk) 
     186                  zwy(ji,jj,jk) = pv_mm(ji,jj,jk) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) 
    187187               END DO 
    188188            END DO 
     
    193193            DO jj = 2, jpjm1       
    194194               DO ji = fs_2, fs_jpim1   ! vector opt. 
    195                   pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( zwx(ji,jj,jk) - zwx(ji-1,jj  ,jk  )       & 
     195                  pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) - ( zwx(ji,jj,jk) - zwx(ji-1,jj  ,jk  )       & 
    196196                  &                                     + zwy(ji,jj,jk) - zwy(ji  ,jj-1,jk  ) )     & 
    197                   &                                   * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 
     197                  &                                   * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    198198               END DO 
    199199           END DO 
     
    201201         !                                ! trend diagnostics 
    202202         IF( l_trd )  THEN 
    203             CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pun, ptb(:,:,:,jn) ) 
    204             CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pvn, ptb(:,:,:,jn) ) 
     203            CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pu_mm, pt(:,:,:,jn,Kbb) ) 
     204            CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pv_mm, pt(:,:,:,jn,Kbb) ) 
    205205         END IF 
    206206         !                                 ! "Poleward" heat and salt transports  
     
    215215         zwx(:,:,jpk) = 0._wp 
    216216         DO jk = 2, jpkm1                       ! interior values 
    217             zwx(:,:,jk) = tmask(:,:,jk) * ( ptb(:,:,jk-1,jn) - ptb(:,:,jk,jn) ) 
     217            zwx(:,:,jk) = tmask(:,:,jk) * ( pt(:,:,jk-1,jn,Kbb) - pt(:,:,jk,jn,Kbb) ) 
    218218         END DO 
    219219         !                                !-- Slopes of tracer 
     
    239239            DO jj = 2, jpjm1       
    240240               DO ji = fs_2, fs_jpim1   ! vector opt. 
    241                   z0w = SIGN( 0.5, pwn(ji,jj,jk+1) ) 
     241                  z0w = SIGN( 0.5, pww(ji,jj,jk+1) ) 
    242242                  zalpha = 0.5 + z0w 
    243                   zw  = z0w - 0.5 * pwn(ji,jj,jk+1) * p2dt * r1_e1e2t(ji,jj) / e3w_n(ji,jj,jk+1) 
    244                   zzwx = ptb(ji,jj,jk+1,jn) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk+1) 
    245                   zzwy = ptb(ji,jj,jk  ,jn) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk  ) 
    246                   zwx(ji,jj,jk+1) = pwn(ji,jj,jk+1) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) * wmask(ji,jj,jk) 
     243                  zw  = z0w - 0.5 * pww(ji,jj,jk+1) * p2dt * r1_e1e2t(ji,jj) / e3w(ji,jj,jk+1,Kmm) 
     244                  zzwx = pt(ji,jj,jk+1,jn,Kbb) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk+1) 
     245                  zzwy = pt(ji,jj,jk  ,jn,Kbb) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk  ) 
     246                  zwx(ji,jj,jk+1) = pww(ji,jj,jk+1) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) * wmask(ji,jj,jk) 
    247247               END DO  
    248248            END DO 
     
    252252               DO jj = 1, jpj 
    253253                  DO ji = 1, jpi 
    254                      zwx(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * ptb(ji,jj,mikt(ji,jj),jn) 
     254                     zwx(ji,jj, mikt(ji,jj) ) = pww(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kbb) 
    255255                  END DO 
    256256               END DO    
    257257            ELSE                                      ! no cavities: only at the ocean surface 
    258                zwx(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn) 
     258               zwx(:,:,1) = pww(:,:,1) * pt(:,:,1,jn,Kbb) 
    259259            ENDIF 
    260260         ENDIF 
     
    263263            DO jj = 2, jpjm1       
    264264               DO ji = fs_2, fs_jpim1   ! vector opt. 
    265                   pta(ji,jj,jk,jn) =  pta(ji,jj,jk,jn) - ( zwx(ji,jj,jk) - zwx(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 
     265                  pt(ji,jj,jk,jn,Krhs) =  pt(ji,jj,jk,jn,Krhs) - ( zwx(ji,jj,jk) - zwx(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    266266               END DO 
    267267            END DO 
    268268         END DO 
    269269         !                                ! send trends for diagnostic 
    270          IF( l_trd )  CALL trd_tra( kt, cdtype, jn, jptra_zad, zwx, pwn, ptb(:,:,:,jn) ) 
     270         IF( l_trd )  CALL trd_tra( kt, cdtype, jn, jptra_zad, zwx, pww, pt(:,:,:,jn,Kbb) ) 
    271271         ! 
    272272      END DO                     ! end of tracer loop 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_qck.F90

    r10874 r10880  
    4747CONTAINS 
    4848 
    49    SUBROUTINE tra_adv_qck ( kt, kit000, cdtype, p2dt, pun, pvn, pwn,      & 
    50       &                                       ptb, ptn, pta, kjpt ) 
     49   SUBROUTINE tra_adv_qck ( kt, kit000, cdtype, p2dt, pu_mm, pv_mm, pww, Kbb, Kmm, pt, kjpt, Krhs ) 
    5150      !!---------------------------------------------------------------------- 
    5251      !!                  ***  ROUTINE tra_adv_qck  *** 
     
    7271      !!         dt = 2*rdtra and the scalar values are tb and sb 
    7372      !! 
    74       !!       On the vertical, the simple centered scheme used ptn 
     73      !!       On the vertical, the simple centered scheme used pt(:,:,:,:,Kmm) 
    7574      !! 
    7675      !!               The fluxes are bounded by the ULTIMATE limiter to 
     
    7877      !!            prevent the appearance of spurious numerical oscillations 
    7978      !! 
    80       !! ** Action : - update pta  with the now advective tracer trends 
     79      !! ** Action : - update pt(:,:,:,:,Krhs)  with the now advective tracer trends 
    8180      !!             - send trends to trdtra module for further diagnostcs (l_trdtra=T) 
    8281      !!             - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 
     
    8483      !! ** Reference : Leonard (1979, 1991) 
    8584      !!---------------------------------------------------------------------- 
    86       INTEGER                              , INTENT(in   ) ::   kt              ! ocean time-step index 
    87       INTEGER                              , INTENT(in   ) ::   kit000          ! first time step index 
    88       CHARACTER(len=3)                     , INTENT(in   ) ::   cdtype          ! =TRA or TRC (tracer indicator) 
    89       INTEGER                              , INTENT(in   ) ::   kjpt            ! number of tracers 
    90       REAL(wp)                             , INTENT(in   ) ::   p2dt            ! tracer time-step 
    91       REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ) ::   pun, pvn, pwn   ! 3 ocean velocity components 
    92       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb, ptn        ! before and now tracer fields 
    93       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta             ! tracer trend  
     85      INTEGER                                  , INTENT(in   ) ::   kt              ! ocean time-step index 
     86      INTEGER                                  , INTENT(in   ) ::   Kbb, Kmm, Krhs  ! ocean time level indices 
     87      INTEGER                                  , INTENT(in   ) ::   kit000          ! first time step index 
     88      CHARACTER(len=3)                         , INTENT(in   ) ::   cdtype          ! =TRA or TRC (tracer indicator) 
     89      INTEGER                                  , INTENT(in   ) ::   kjpt            ! number of tracers 
     90      REAL(wp)                                 , INTENT(in   ) ::   p2dt            ! tracer time-step 
     91      REAL(wp), DIMENSION(jpi,jpj,jpk         ), INTENT(in   ) ::   pu_mm, pv_mm, pww   ! 3 ocean velocity components 
     92      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) ::   pt              ! active tracers and RHS of tracer equation 
    9493      !!---------------------------------------------------------------------- 
    9594      ! 
     
    108107      ! 
    109108      !        ! horizontal fluxes are computed with the QUICKEST + ULTIMATE scheme 
    110       CALL tra_adv_qck_i( kt, cdtype, p2dt, pun, ptb, ptn, pta, kjpt )  
    111       CALL tra_adv_qck_j( kt, cdtype, p2dt, pvn, ptb, ptn, pta, kjpt )  
     109      CALL tra_adv_qck_i( kt, cdtype, p2dt, pu_mm, Kbb, Kmm, pt, kjpt, Krhs )  
     110      CALL tra_adv_qck_j( kt, cdtype, p2dt, pv_mm, Kbb, Kmm, pt, kjpt, Krhs )  
    112111 
    113112      !        ! vertical fluxes are computed with the 2nd order centered scheme 
    114       CALL tra_adv_cen2_k( kt, cdtype, pwn,         ptn, pta, kjpt ) 
     113      CALL tra_adv_cen2_k( kt, cdtype, pww, Kmm, pt, kjpt, Krhs ) 
    115114      ! 
    116115   END SUBROUTINE tra_adv_qck 
    117116 
    118117 
    119    SUBROUTINE tra_adv_qck_i( kt, cdtype, p2dt, pun,                  & 
    120       &                                        ptb, ptn, pta, kjpt   ) 
    121       !!---------------------------------------------------------------------- 
    122       !! 
    123       !!---------------------------------------------------------------------- 
    124       INTEGER                              , INTENT(in   ) ::   kt         ! ocean time-step index 
    125       CHARACTER(len=3)                     , INTENT(in   ) ::   cdtype     ! =TRA or TRC (tracer indicator) 
    126       INTEGER                              , INTENT(in   ) ::   kjpt       ! number of tracers 
    127       REAL(wp)                             , INTENT(in   ) ::   p2dt       ! tracer time-step 
    128       REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ) ::   pun        ! i-velocity components 
    129       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb, ptn   ! before and now tracer fields 
    130       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta        ! tracer trend  
     118   SUBROUTINE tra_adv_qck_i( kt, cdtype, p2dt, pu_mm, Kbb, Kmm, pt, kjpt, Krhs ) 
     119      !!---------------------------------------------------------------------- 
     120      !! 
     121      !!---------------------------------------------------------------------- 
     122      INTEGER                                  , INTENT(in   ) ::   kt         ! ocean time-step index 
     123      INTEGER                                  , INTENT(in   ) ::   Kbb, Kmm, Krhs  ! ocean time level indices 
     124      CHARACTER(len=3)                         , INTENT(in   ) ::   cdtype     ! =TRA or TRC (tracer indicator) 
     125      INTEGER                                  , INTENT(in   ) ::   kjpt       ! number of tracers 
     126      REAL(wp)                                 , INTENT(in   ) ::   p2dt       ! tracer time-step 
     127      REAL(wp), DIMENSION(jpi,jpj,jpk         ), INTENT(in   ) ::   pu_mm        ! i-velocity components 
     128      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) ::   pt              ! active tracers and RHS of tracer equation 
    131129      !! 
    132130      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices 
     
    145143            DO jj = 2, jpjm1 
    146144               DO ji = fs_2, fs_jpim1   ! vector opt. 
    147                   zfc(ji,jj,jk) = ptb(ji-1,jj,jk,jn)        ! Upstream   in the x-direction for the tracer 
    148                   zfd(ji,jj,jk) = ptb(ji+1,jj,jk,jn)        ! Downstream in the x-direction for the tracer 
     145                  zfc(ji,jj,jk) = pt(ji-1,jj,jk,jn,Kbb)        ! Upstream   in the x-direction for the tracer 
     146                  zfd(ji,jj,jk) = pt(ji+1,jj,jk,jn,Kbb)        ! Downstream in the x-direction for the tracer 
    149147               END DO 
    150148            END DO 
     
    158156            DO jj = 2, jpjm1 
    159157               DO ji = fs_2, fs_jpim1   ! vector opt.          
    160                   zdir = 0.5 + SIGN( 0.5, pun(ji,jj,jk) )   ! if pun > 0 : zdir = 1 otherwise zdir = 0  
     158                  zdir = 0.5 + SIGN( 0.5, pu_mm(ji,jj,jk) )   ! if pu_mm > 0 : zdir = 1 otherwise zdir = 0  
    161159                  zfu(ji,jj,jk) = zdir * zfc(ji,jj,jk ) + ( 1. - zdir ) * zfd(ji+1,jj,jk)  ! FU in the x-direction for T  
    162160               END DO 
     
    167165            DO jj = 2, jpjm1 
    168166               DO ji = fs_2, fs_jpim1   ! vector opt.    
    169                   zdir = 0.5 + SIGN( 0.5, pun(ji,jj,jk) )   ! if pun > 0 : zdir = 1 otherwise zdir = 0  
    170                   zdx = ( zdir * e1t(ji,jj) + ( 1. - zdir ) * e1t(ji+1,jj) ) * e2u(ji,jj) * e3u_n(ji,jj,jk) 
    171                   zwx(ji,jj,jk)  = ABS( pun(ji,jj,jk) ) * p2dt / zdx    ! (0<zc_cfl<1 : Courant number on x-direction) 
    172                   zfc(ji,jj,jk)  = zdir * ptb(ji  ,jj,jk,jn) + ( 1. - zdir ) * ptb(ji+1,jj,jk,jn)  ! FC in the x-direction for T 
    173                   zfd(ji,jj,jk)  = zdir * ptb(ji+1,jj,jk,jn) + ( 1. - zdir ) * ptb(ji  ,jj,jk,jn)  ! FD in the x-direction for T 
     167                  zdir = 0.5 + SIGN( 0.5, pu_mm(ji,jj,jk) )   ! if pu_mm > 0 : zdir = 1 otherwise zdir = 0  
     168                  zdx = ( zdir * e1t(ji,jj) + ( 1. - zdir ) * e1t(ji+1,jj) ) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm) 
     169                  zwx(ji,jj,jk)  = ABS( pu_mm(ji,jj,jk) ) * p2dt / zdx    ! (0<zc_cfl<1 : Courant number on x-direction) 
     170                  zfc(ji,jj,jk)  = zdir * pt(ji  ,jj,jk,jn,Kbb) + ( 1. - zdir ) * pt(ji+1,jj,jk,jn,Kbb)  ! FC in the x-direction for T 
     171                  zfd(ji,jj,jk)  = zdir * pt(ji+1,jj,jk,jn,Kbb) + ( 1. - zdir ) * pt(ji  ,jj,jk,jn,Kbb)  ! FD in the x-direction for T 
    174172               END DO 
    175173            END DO 
     
    197195            DO jj = 2, jpjm1 
    198196               DO ji = fs_2, fs_jpim1   ! vector opt.                
    199                   zdir = 0.5 + SIGN( 0.5, pun(ji,jj,jk) )   ! if pun > 0 : zdir = 1 otherwise zdir = 0  
     197                  zdir = 0.5 + SIGN( 0.5, pu_mm(ji,jj,jk) )   ! if pu_mm > 0 : zdir = 1 otherwise zdir = 0  
    200198                  !--- If the second ustream point is a land point 
    201199                  !--- the flux is computed by the 1st order UPWIND scheme 
    202200                  zmsk = zdir * zfu(ji,jj,jk) + ( 1. - zdir ) * zfu(ji+1,jj,jk) 
    203201                  zwx(ji,jj,jk) = zmsk * zwx(ji,jj,jk) + ( 1. - zmsk ) * zfc(ji,jj,jk) 
    204                   zwx(ji,jj,jk) = zwx(ji,jj,jk) * pun(ji,jj,jk) 
     202                  zwx(ji,jj,jk) = zwx(ji,jj,jk) * pu_mm(ji,jj,jk) 
    205203               END DO 
    206204            END DO 
     
    213211            DO jj = 2, jpjm1 
    214212               DO ji = fs_2, fs_jpim1   ! vector opt.   
    215                   zbtr = r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 
     213                  zbtr = r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    216214                  ! horizontal advective trends 
    217215                  ztra = - zbtr * ( zwx(ji,jj,jk) - zwx(ji-1,jj,jk) ) 
    218216                  !--- add it to the general tracer trends 
    219                   pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra 
     217                  pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) + ztra 
    220218               END DO 
    221219            END DO 
    222220         END DO 
    223221         !                                 ! trend diagnostics 
    224          IF( l_trd )   CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pun, ptn(:,:,:,jn) ) 
     222         IF( l_trd )   CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pu_mm, pt(:,:,:,jn,Kmm) ) 
    225223         ! 
    226224      END DO 
     
    229227 
    230228 
    231    SUBROUTINE tra_adv_qck_j( kt, cdtype, p2dt, pvn,                & 
    232       &                                        ptb, ptn, pta, kjpt ) 
    233       !!---------------------------------------------------------------------- 
    234       !! 
    235       !!---------------------------------------------------------------------- 
    236       INTEGER                              , INTENT(in   ) ::   kt         ! ocean time-step index 
    237       CHARACTER(len=3)                     , INTENT(in   ) ::   cdtype     ! =TRA or TRC (tracer indicator) 
    238       INTEGER                              , INTENT(in   ) ::   kjpt       ! number of tracers 
    239       REAL(wp)                             , INTENT(in   ) ::   p2dt       ! tracer time-step 
    240       REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ) ::   pvn        ! j-velocity components 
    241       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb, ptn   ! before and now tracer fields 
    242       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta        ! tracer trend  
     229   SUBROUTINE tra_adv_qck_j( kt, cdtype, p2dt, pv_mm, Kbb, Kmm, pt, kjpt, Krhs ) 
     230      !!---------------------------------------------------------------------- 
     231      !! 
     232      !!---------------------------------------------------------------------- 
     233      INTEGER                                  , INTENT(in   ) ::   kt         ! ocean time-step index 
     234      INTEGER                                  , INTENT(in   ) ::   Kbb, Kmm, Krhs  ! ocean time level indices 
     235      CHARACTER(len=3)                         , INTENT(in   ) ::   cdtype     ! =TRA or TRC (tracer indicator) 
     236      INTEGER                                  , INTENT(in   ) ::   kjpt       ! number of tracers 
     237      REAL(wp)                                 , INTENT(in   ) ::   p2dt       ! tracer time-step 
     238      REAL(wp), DIMENSION(jpi,jpj,jpk         ), INTENT(in   ) ::   pv_mm        ! j-velocity components 
     239      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) ::   pt              ! active tracers and RHS of tracer equation 
    243240      !! 
    244241      INTEGER  :: ji, jj, jk, jn                ! dummy loop indices 
     
    259256               DO ji = fs_2, fs_jpim1   ! vector opt. 
    260257                  ! Upstream in the x-direction for the tracer 
    261                   zfc(ji,jj,jk) = ptb(ji,jj-1,jk,jn) 
     258                  zfc(ji,jj,jk) = pt(ji,jj-1,jk,jn,Kbb) 
    262259                  ! Downstream in the x-direction for the tracer 
    263                   zfd(ji,jj,jk) = ptb(ji,jj+1,jk,jn) 
     260                  zfd(ji,jj,jk) = pt(ji,jj+1,jk,jn,Kbb) 
    264261               END DO 
    265262            END DO 
     
    275272            DO jj = 2, jpjm1 
    276273               DO ji = fs_2, fs_jpim1   ! vector opt.          
    277                   zdir = 0.5 + SIGN( 0.5, pvn(ji,jj,jk) )   ! if pun > 0 : zdir = 1 otherwise zdir = 0  
     274                  zdir = 0.5 + SIGN( 0.5, pv_mm(ji,jj,jk) )   ! if pu_mm > 0 : zdir = 1 otherwise zdir = 0  
    278275                  zfu(ji,jj,jk) = zdir * zfc(ji,jj,jk ) + ( 1. - zdir ) * zfd(ji,jj+1,jk)  ! FU in the x-direction for T  
    279276               END DO 
     
    284281            DO jj = 2, jpjm1 
    285282               DO ji = fs_2, fs_jpim1   ! vector opt.    
    286                   zdir = 0.5 + SIGN( 0.5, pvn(ji,jj,jk) )   ! if pun > 0 : zdir = 1 otherwise zdir = 0  
    287                   zdx = ( zdir * e2t(ji,jj) + ( 1. - zdir ) * e2t(ji,jj+1) ) * e1v(ji,jj) * e3v_n(ji,jj,jk) 
    288                   zwy(ji,jj,jk)  = ABS( pvn(ji,jj,jk) ) * p2dt / zdx    ! (0<zc_cfl<1 : Courant number on x-direction) 
    289                   zfc(ji,jj,jk)  = zdir * ptb(ji,jj  ,jk,jn) + ( 1. - zdir ) * ptb(ji,jj+1,jk,jn)  ! FC in the x-direction for T 
    290                   zfd(ji,jj,jk)  = zdir * ptb(ji,jj+1,jk,jn) + ( 1. - zdir ) * ptb(ji,jj  ,jk,jn)  ! FD in the x-direction for T 
     283                  zdir = 0.5 + SIGN( 0.5, pv_mm(ji,jj,jk) )   ! if pu_mm > 0 : zdir = 1 otherwise zdir = 0  
     284                  zdx = ( zdir * e2t(ji,jj) + ( 1. - zdir ) * e2t(ji,jj+1) ) * e1v(ji,jj) * e3v(ji,jj,jk,Kmm) 
     285                  zwy(ji,jj,jk)  = ABS( pv_mm(ji,jj,jk) ) * p2dt / zdx    ! (0<zc_cfl<1 : Courant number on x-direction) 
     286                  zfc(ji,jj,jk)  = zdir * pt(ji,jj  ,jk,jn,Kbb) + ( 1. - zdir ) * pt(ji,jj+1,jk,jn,Kbb)  ! FC in the x-direction for T 
     287                  zfd(ji,jj,jk)  = zdir * pt(ji,jj+1,jk,jn,Kbb) + ( 1. - zdir ) * pt(ji,jj  ,jk,jn,Kbb)  ! FD in the x-direction for T 
    291288               END DO 
    292289            END DO 
     
    314311            DO jj = 2, jpjm1 
    315312               DO ji = fs_2, fs_jpim1   ! vector opt.                
    316                   zdir = 0.5 + SIGN( 0.5, pvn(ji,jj,jk) )   ! if pun > 0 : zdir = 1 otherwise zdir = 0  
     313                  zdir = 0.5 + SIGN( 0.5, pv_mm(ji,jj,jk) )   ! if pu_mm > 0 : zdir = 1 otherwise zdir = 0  
    317314                  !--- If the second ustream point is a land point 
    318315                  !--- the flux is computed by the 1st order UPWIND scheme 
    319316                  zmsk = zdir * zfu(ji,jj,jk) + ( 1. - zdir ) * zfu(ji,jj+1,jk) 
    320317                  zwy(ji,jj,jk) = zmsk * zwy(ji,jj,jk) + ( 1. - zmsk ) * zfc(ji,jj,jk) 
    321                   zwy(ji,jj,jk) = zwy(ji,jj,jk) * pvn(ji,jj,jk) 
     318                  zwy(ji,jj,jk) = zwy(ji,jj,jk) * pv_mm(ji,jj,jk) 
    322319               END DO 
    323320            END DO 
     
    330327            DO jj = 2, jpjm1 
    331328               DO ji = fs_2, fs_jpim1   ! vector opt.   
    332                   zbtr = r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 
     329                  zbtr = r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    333330                  ! horizontal advective trends 
    334331                  ztra = - zbtr * ( zwy(ji,jj,jk) - zwy(ji,jj-1,jk) ) 
    335332                  !--- add it to the general tracer trends 
    336                   pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra 
     333                  pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) + ztra 
    337334               END DO 
    338335            END DO 
    339336         END DO 
    340337         !                                 ! trend diagnostics 
    341          IF( l_trd )   CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pvn, ptn(:,:,:,jn) ) 
     338         IF( l_trd )   CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pv_mm, pt(:,:,:,jn,Kmm) ) 
    342339         !                                 ! "Poleward" heat and salt transports (contribution of upstream fluxes) 
    343340         IF( l_ptr )   CALL dia_ptr_hst( jn, 'adv', zwy(:,:,:) ) 
     
    348345 
    349346 
    350    SUBROUTINE tra_adv_cen2_k( kt, cdtype, pwn,           & 
    351      &                                    ptn, pta, kjpt ) 
    352       !!---------------------------------------------------------------------- 
    353       !! 
    354       !!---------------------------------------------------------------------- 
    355       INTEGER                              , INTENT(in   ) ::   kt       ! ocean time-step index 
    356       CHARACTER(len=3)                     , INTENT(in   ) ::   cdtype   ! =TRA or TRC (tracer indicator) 
    357       INTEGER                              , INTENT(in   ) ::   kjpt     ! number of tracers 
    358       REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ) ::   pwn      ! vertical velocity  
    359       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptn      ! before and now tracer fields 
    360       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta      ! tracer trend  
     347   SUBROUTINE tra_adv_cen2_k( kt, cdtype, pww, Kmm, pt, kjpt, Krhs ) 
     348      !!---------------------------------------------------------------------- 
     349      !! 
     350      !!---------------------------------------------------------------------- 
     351      INTEGER                                  , INTENT(in   ) ::   kt       ! ocean time-step index 
     352      INTEGER                                  , INTENT(in   ) ::   Kmm, Krhs  ! ocean time level indices 
     353      CHARACTER(len=3)                         , INTENT(in   ) ::   cdtype   ! =TRA or TRC (tracer indicator) 
     354      INTEGER                                  , INTENT(in   ) ::   kjpt     ! number of tracers 
     355      REAL(wp), DIMENSION(jpi,jpj,jpk         ), INTENT(in   ) ::   pww      ! vertical velocity  
     356      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) ::   pt              ! active tracers and RHS of tracer equation 
    361357      ! 
    362358      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices 
     
    374370            DO jj = 2, jpjm1 
    375371               DO ji = fs_2, fs_jpim1   ! vector opt. 
    376                   zwz(ji,jj,jk) = 0.5 * pwn(ji,jj,jk) * ( ptn(ji,jj,jk-1,jn) + ptn(ji,jj,jk,jn) ) * wmask(ji,jj,jk) 
     372                  zwz(ji,jj,jk) = 0.5 * pww(ji,jj,jk) * ( pt(ji,jj,jk-1,jn,Kmm) + pt(ji,jj,jk,jn,Kmm) ) * wmask(ji,jj,jk) 
    377373               END DO 
    378374            END DO 
     
    382378               DO jj = 1, jpj 
    383379                  DO ji = 1, jpi 
    384                      zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * ptn(ji,jj,mikt(ji,jj),jn)   ! linear free surface  
     380                     zwz(ji,jj, mikt(ji,jj) ) = pww(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kmm)   ! linear free surface  
    385381                  END DO 
    386382               END DO    
    387383            ELSE                                   ! no ocean cavities (only ocean surface) 
    388                zwz(:,:,1) = pwn(:,:,1) * ptn(:,:,1,jn) 
     384               zwz(:,:,1) = pww(:,:,1) * pt(:,:,1,jn,Kmm) 
    389385            ENDIF 
    390386         ENDIF 
     
    393389            DO jj = 2, jpjm1 
    394390               DO ji = fs_2, fs_jpim1   ! vector opt. 
    395                   pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( zwz(ji,jj,jk) - zwz(ji,jj,jk+1) )   & 
    396                      &                                * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 
     391                  pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) - ( zwz(ji,jj,jk) - zwz(ji,jj,jk+1) )   & 
     392                     &                                * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    397393               END DO 
    398394            END DO 
    399395         END DO 
    400396         !                                 ! Send trends for diagnostic 
    401          IF( l_trd )  CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pwn, ptn(:,:,:,jn) ) 
     397         IF( l_trd )  CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pww, pt(:,:,:,jn,Kmm) ) 
    402398         ! 
    403399      END DO 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_ubs.F90

    r10874 r10880  
    4646CONTAINS 
    4747 
    48    SUBROUTINE tra_adv_ubs( kt, kit000, cdtype, p2dt, pun, pvn, pwn,          & 
    49       &                                                ptb, ptn, pta, kjpt, kn_ubs_v ) 
     48   SUBROUTINE tra_adv_ubs( kt, kit000, cdtype, p2dt, pu_mm, pv_mm, pww,          & 
     49      &                    Kbb, Kmm, pt, kjpt, Krhs, kn_ubs_v ) 
    5050      !!---------------------------------------------------------------------- 
    5151      !!                  ***  ROUTINE tra_adv_ubs  *** 
     
    5858      !!      It is only used in the horizontal direction. 
    5959      !!      For example the i-component of the advective fluxes are given by : 
    60       !!                !  e2u e3u un ( mi(Tn) - zltu(i  ) )   if un(i) >= 0 
     60      !!                !  e2u e3u uu ( mi(Tn) - zltu(i  ) ,Kmm)   if uu(i,Kmm) >= 0 
    6161      !!          ztu = !  or  
    62       !!                !  e2u e3u un ( mi(Tn) - zltu(i+1) )   if un(i) < 0 
     62      !!                !  e2u e3u uu ( mi(Tn) - zltu(i+1) ,Kmm)   if uu(i,Kmm) < 0 
    6363      !!      where zltu is the second derivative of the before temperature field: 
    6464      !!          zltu = 1/e3t di[ e2u e3u / e1u di[Tb] ] 
     
    7777      !!      scheme (kn_ubs_v=4). 
    7878      !! 
    79       !! ** Action : - update pta  with the now advective tracer trends 
     79      !! ** Action : - update pt(:,:,:,:,Krhs)  with the now advective tracer trends 
    8080      !!             - send trends to trdtra module for further diagnostcs (l_trdtra=T) 
    8181      !!             - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 
     
    8484      !!             Farrow, D.E., Stevens, D.P., 1995, J. Phys. Ocean. 25, 1731Ð1741.  
    8585      !!---------------------------------------------------------------------- 
    86       INTEGER                              , INTENT(in   ) ::   kt              ! ocean time-step index 
    87       INTEGER                              , INTENT(in   ) ::   kit000          ! first time step index 
    88       CHARACTER(len=3)                     , INTENT(in   ) ::   cdtype          ! =TRA or TRC (tracer indicator) 
    89       INTEGER                              , INTENT(in   ) ::   kjpt            ! number of tracers 
    90       INTEGER                              , INTENT(in   ) ::   kn_ubs_v        ! number of tracers 
    91       REAL(wp)                             , INTENT(in   ) ::   p2dt            ! tracer time-step 
    92       REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ) ::   pun, pvn, pwn   ! 3 ocean transport components 
    93       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb, ptn        ! before and now tracer fields 
    94       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta             ! tracer trend  
     86      INTEGER                                  , INTENT(in   ) ::   kt              ! ocean time-step index 
     87      INTEGER                                  , INTENT(in   ) ::   Kbb, Kmm, Krhs  ! ocean time level indices 
     88      INTEGER                                  , INTENT(in   ) ::   kit000          ! first time step index 
     89      CHARACTER(len=3)                         , INTENT(in   ) ::   cdtype          ! =TRA or TRC (tracer indicator) 
     90      INTEGER                                  , INTENT(in   ) ::   kjpt            ! number of tracers 
     91      INTEGER                                  , INTENT(in   ) ::   kn_ubs_v        ! number of tracers 
     92      REAL(wp)                                 , INTENT(in   ) ::   p2dt            ! tracer time-step 
     93      REAL(wp), DIMENSION(jpi,jpj,jpk         ), INTENT(in   ) ::   pu_mm, pv_mm, pww   ! 3 ocean transport components 
     94      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt,jpt), INTENT(inout) ::   pt              ! active tracers and RHS of tracer equation 
    9595      ! 
    9696      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices 
     
    126126            DO jj = 1, jpjm1              ! First derivative (masked gradient) 
    127127               DO ji = 1, fs_jpim1   ! vector opt. 
    128                   zeeu = e2_e1u(ji,jj) * e3u_n(ji,jj,jk) * umask(ji,jj,jk) 
    129                   zeev = e1_e2v(ji,jj) * e3v_n(ji,jj,jk) * vmask(ji,jj,jk) 
    130                   ztu(ji,jj,jk) = zeeu * ( ptb(ji+1,jj  ,jk,jn) - ptb(ji,jj,jk,jn) ) 
    131                   ztv(ji,jj,jk) = zeev * ( ptb(ji  ,jj+1,jk,jn) - ptb(ji,jj,jk,jn) ) 
     128                  zeeu = e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm) * umask(ji,jj,jk) 
     129                  zeev = e1_e2v(ji,jj) * e3v(ji,jj,jk,Kmm) * vmask(ji,jj,jk) 
     130                  ztu(ji,jj,jk) = zeeu * ( pt(ji+1,jj  ,jk,jn,Kbb) - pt(ji,jj,jk,jn,Kbb) ) 
     131                  ztv(ji,jj,jk) = zeev * ( pt(ji  ,jj+1,jk,jn,Kbb) - pt(ji,jj,jk,jn,Kbb) ) 
    132132               END DO 
    133133            END DO 
    134134            DO jj = 2, jpjm1              ! Second derivative (divergence) 
    135135               DO ji = fs_2, fs_jpim1   ! vector opt. 
    136                   zcoef = 1._wp / ( 6._wp * e3t_n(ji,jj,jk) ) 
     136                  zcoef = 1._wp / ( 6._wp * e3t(ji,jj,jk,Kmm) ) 
    137137                  zltu(ji,jj,jk) = (  ztu(ji,jj,jk) - ztu(ji-1,jj,jk)  ) * zcoef 
    138138                  zltv(ji,jj,jk) = (  ztv(ji,jj,jk) - ztv(ji,jj-1,jk)  ) * zcoef 
     
    146146            DO jj = 1, jpjm1 
    147147               DO ji = 1, fs_jpim1   ! vector opt. 
    148                   zfp_ui = pun(ji,jj,jk) + ABS( pun(ji,jj,jk) )      ! upstream transport (x2) 
    149                   zfm_ui = pun(ji,jj,jk) - ABS( pun(ji,jj,jk) ) 
    150                   zfp_vj = pvn(ji,jj,jk) + ABS( pvn(ji,jj,jk) ) 
    151                   zfm_vj = pvn(ji,jj,jk) - ABS( pvn(ji,jj,jk) ) 
     148                  zfp_ui = pu_mm(ji,jj,jk) + ABS( pu_mm(ji,jj,jk) )      ! upstream transport (x2) 
     149                  zfm_ui = pu_mm(ji,jj,jk) - ABS( pu_mm(ji,jj,jk) ) 
     150                  zfp_vj = pv_mm(ji,jj,jk) + ABS( pv_mm(ji,jj,jk) ) 
     151                  zfm_vj = pv_mm(ji,jj,jk) - ABS( pv_mm(ji,jj,jk) ) 
    152152                  !                                                  ! 2nd order centered advective fluxes (x2) 
    153                   zcenut = pun(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji+1,jj  ,jk,jn) ) 
    154                   zcenvt = pvn(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji  ,jj+1,jk,jn) ) 
     153                  zcenut = pu_mm(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji+1,jj  ,jk,jn,Kmm) ) 
     154                  zcenvt = pv_mm(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji  ,jj+1,jk,jn,Kmm) ) 
    155155                  !                                                  ! UBS advective fluxes 
    156156                  ztu(ji,jj,jk) = 0.5 * ( zcenut - zfp_ui * zltu(ji,jj,jk) - zfm_ui * zltu(ji+1,jj,jk) ) 
     
    160160         END DO          
    161161         ! 
    162          zltu(:,:,:) = pta(:,:,:,jn)      ! store the initial trends before its update 
     162         zltu(:,:,:) = pt(:,:,:,jn,Krhs)      ! store the initial trends before its update 
    163163         ! 
    164164         DO jk = 1, jpkm1        !==  add the horizontal advective trend  ==! 
    165165            DO jj = 2, jpjm1 
    166166               DO ji = fs_2, fs_jpim1   ! vector opt. 
    167                   pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn)                        & 
     167                  pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs)                        & 
    168168                     &             - (  ztu(ji,jj,jk) - ztu(ji-1,jj  ,jk)    & 
    169                      &                + ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk)  ) * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 
     169                     &                + ztv(ji,jj,jk) - ztv(ji  ,jj-1,jk)  ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    170170               END DO 
    171171            END DO 
     
    173173         END DO 
    174174         ! 
    175          zltu(:,:,:) = pta(:,:,:,jn) - zltu(:,:,:)    ! Horizontal advective trend used in vertical 2nd order FCT case 
     175         zltu(:,:,:) = pt(:,:,:,jn,Krhs) - zltu(:,:,:)    ! Horizontal advective trend used in vertical 2nd order FCT case 
    176176         !                                            ! and/or in trend diagnostic (l_trd=T)  
    177177         !                 
    178178         IF( l_trd ) THEN                  ! trend diagnostics 
    179              CALL trd_tra( kt, cdtype, jn, jptra_xad, ztu, pun, ptn(:,:,:,jn) ) 
    180              CALL trd_tra( kt, cdtype, jn, jptra_yad, ztv, pvn, ptn(:,:,:,jn) ) 
     179             CALL trd_tra( kt, cdtype, jn, jptra_xad, ztu, pu_mm, pt(:,:,:,jn,Kmm) ) 
     180             CALL trd_tra( kt, cdtype, jn, jptra_yad, ztv, pv_mm, pt(:,:,:,jn,Kmm) ) 
    181181         END IF 
    182182         !      
     
    193193         CASE(  2  )                   ! 2nd order FCT  
    194194            !          
    195             IF( l_trd )   zltv(:,:,:) = pta(:,:,:,jn)          ! store pta if trend diag. 
     195            IF( l_trd )   zltv(:,:,:) = pt(:,:,:,jn,Krhs)          ! store pt(:,:,:,:,Krhs) if trend diag. 
    196196            ! 
    197197            !                          !*  upstream advection with initial mass fluxes & intermediate update  ==! 
     
    199199               DO jj = 1, jpj 
    200200                  DO ji = 1, jpi 
    201                      zfp_wk = pwn(ji,jj,jk) + ABS( pwn(ji,jj,jk) ) 
    202                      zfm_wk = pwn(ji,jj,jk) - ABS( pwn(ji,jj,jk) ) 
    203                      ztw(ji,jj,jk) = 0.5_wp * (  zfp_wk * ptb(ji,jj,jk,jn) + zfm_wk * ptb(ji,jj,jk-1,jn)  ) * wmask(ji,jj,jk) 
     201                     zfp_wk = pww(ji,jj,jk) + ABS( pww(ji,jj,jk) ) 
     202                     zfm_wk = pww(ji,jj,jk) - ABS( pww(ji,jj,jk) ) 
     203                     ztw(ji,jj,jk) = 0.5_wp * (  zfp_wk * pt(ji,jj,jk,jn,Kbb) + zfm_wk * pt(ji,jj,jk-1,jn,Kbb)  ) * wmask(ji,jj,jk) 
    204204                  END DO 
    205205               END DO 
     
    209209                  DO jj = 1, jpj 
    210210                     DO ji = 1, jpi 
    211                         ztw(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * ptb(ji,jj,mikt(ji,jj),jn)   ! linear free surface  
     211                        ztw(ji,jj, mikt(ji,jj) ) = pww(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn,Kbb)   ! linear free surface  
    212212                     END DO 
    213213                  END DO    
    214214               ELSE                                ! no cavities: only at the ocean surface 
    215                   ztw(:,:,1) = pwn(:,:,1) * ptb(:,:,1,jn) 
     215                  ztw(:,:,1) = pww(:,:,1) * pt(:,:,1,jn,Kbb) 
    216216               ENDIF 
    217217            ENDIF 
     
    220220               DO jj = 2, jpjm1 
    221221                  DO ji = fs_2, fs_jpim1   ! vector opt. 
    222                      ztak = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 
    223                      pta(ji,jj,jk,jn) =   pta(ji,jj,jk,jn) +  ztak  
    224                      zti(ji,jj,jk)    = ( ptb(ji,jj,jk,jn) + p2dt * ( ztak + zltu(ji,jj,jk) ) ) * tmask(ji,jj,jk) 
     222                     ztak = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     223                     pt(ji,jj,jk,jn,Krhs) =   pt(ji,jj,jk,jn,Krhs) +  ztak  
     224                     zti(ji,jj,jk)    = ( pt(ji,jj,jk,jn,Kbb) + p2dt * ( ztak + zltu(ji,jj,jk) ) ) * tmask(ji,jj,jk) 
    225225                  END DO 
    226226               END DO 
     
    232232               DO jj = 1, jpj 
    233233                  DO ji = 1, jpi 
    234                      ztw(ji,jj,jk) = (   0.5_wp * pwn(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji,jj,jk-1,jn) )   & 
     234                     ztw(ji,jj,jk) = (   0.5_wp * pww(ji,jj,jk) * ( pt(ji,jj,jk,jn,Kmm) + pt(ji,jj,jk-1,jn,Kmm) )   & 
    235235                        &              - ztw(ji,jj,jk)   ) * wmask(ji,jj,jk) 
    236236                  END DO 
     
    240240            IF( ln_linssh )   ztw(:,:, 1 ) = 0._wp       ! only ocean surface as interior zwz values have been w-masked 
    241241            ! 
    242             CALL nonosc_z( ptb(:,:,:,jn), ztw, zti, p2dt )      !  monotonicity algorithm 
     242            CALL nonosc_z( Kmm, pt(:,:,:,jn,Kbb), ztw, zti, p2dt )      !  monotonicity algorithm 
    243243            ! 
    244244         CASE(  4  )                               ! 4th order COMPACT 
    245             CALL interp_4th_cpt( ptn(:,:,:,jn) , ztw )         ! 4th order compact interpolation of T at w-point 
     245            CALL interp_4th_cpt( pt(:,:,:,jn,Kmm) , ztw )         ! 4th order compact interpolation of T at w-point 
    246246            DO jk = 2, jpkm1 
    247247               DO jj = 2, jpjm1 
    248248                  DO ji = fs_2, fs_jpim1 
    249                      ztw(ji,jj,jk) = pwn(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk) 
    250                   END DO 
    251                END DO 
    252             END DO 
    253             IF( ln_linssh )   ztw(:,:, 1 ) = pwn(:,:,1) * ptn(:,:,1,jn)     !!gm ISF & 4th COMPACT doesn't work 
     249                     ztw(ji,jj,jk) = pww(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk) 
     250                  END DO 
     251               END DO 
     252            END DO 
     253            IF( ln_linssh )   ztw(:,:, 1 ) = pww(:,:,1) * pt(:,:,1,jn,Kmm)     !!gm ISF & 4th COMPACT doesn't work 
    254254            ! 
    255255         END SELECT 
     
    258258            DO jj = 2, jpjm1  
    259259               DO ji = fs_2, fs_jpim1   ! vector opt.    
    260                   pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 
     260                  pt(ji,jj,jk,jn,Krhs) = pt(ji,jj,jk,jn,Krhs) - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    261261               END DO 
    262262            END DO 
     
    267267               DO jj = 2, jpjm1 
    268268                  DO ji = fs_2, fs_jpim1   ! vector opt. 
    269                      zltv(ji,jj,jk) = pta(ji,jj,jk,jn) - zltv(ji,jj,jk)                          & 
    270                         &           + ptn(ji,jj,jk,jn) * (  pwn(ji,jj,jk) - pwn(ji,jj,jk+1)  )   & 
    271                         &                              * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk) 
     269                     zltv(ji,jj,jk) = pt(ji,jj,jk,jn,Krhs) - zltv(ji,jj,jk)                          & 
     270                        &           + pt(ji,jj,jk,jn,Kmm) * (  pww(ji,jj,jk) - pww(ji,jj,jk+1)  )   & 
     271                        &                              * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    272272                  END DO 
    273273               END DO 
     
    281281 
    282282 
    283    SUBROUTINE nonosc_z( pbef, pcc, paft, p2dt ) 
     283   SUBROUTINE nonosc_z( Kmm, pbef, pcc, paft, p2dt ) 
    284284      !!--------------------------------------------------------------------- 
    285285      !!                    ***  ROUTINE nonosc_z  *** 
     
    294294      !!       in-space based differencing for fluid 
    295295      !!---------------------------------------------------------------------- 
     296      INTEGER , INTENT(in   )                          ::   Kmm    ! time level index 
    296297      REAL(wp), INTENT(in   )                          ::   p2dt   ! tracer time-step 
    297298      REAL(wp),                DIMENSION (jpi,jpj,jpk) ::   pbef   ! before field 
     
    352353               zneg = MAX( 0., pcc(ji  ,jj  ,jk  ) ) - MIN( 0., pcc(ji  ,jj  ,jk+1) ) 
    353354               ! up & down beta terms 
    354                zbt = e1e2t(ji,jj) * e3t_n(ji,jj,jk) / p2dt 
     355               zbt = e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) / p2dt 
    355356               zbetup(ji,jj,jk) = ( zbetup(ji,jj,jk) - paft(ji,jj,jk) ) / (zpos+zrtrn) * zbt 
    356357               zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zbetdo(ji,jj,jk) ) / (zneg+zrtrn) * zbt 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/step.F90

    r10877 r10880  
    4747 
    4848   !!---------------------------------------------------------------------- 
     49   !! time level indices 
     50   !!---------------------------------------------------------------------- 
     51   INTEGER, PUBLIC :: Nbb, Nnn, Naa, Nrhs          !! used by nemo_init 
     52 
     53   !!---------------------------------------------------------------------- 
    4954   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
    5055   !! $Id$ 
     
    247252               &         CALL Agrif_Sponge_tra        ! tracers sponge 
    248253#endif 
    249                          CALL tra_adv       ( kstp )  ! horizontal & vertical advection 
     254                         CALL tra_adv( kstp, Nbb, Nnn      , ts, Nrhs )  ! hor. + vert. advection  ==> RHS 
    250255      IF( ln_zdfosm  )   CALL tra_osm       ( kstp )  ! OSMOSIS non-local tracer fluxes 
    251256      IF( lrst_oce .AND. ln_zdfosm ) & 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/step_oce.F90

    r10068 r10880  
    107107#endif 
    108108#if defined key_top 
    109    USE trcstp           ! passive tracer time-stepping      (trc_stp routine) 
     109   USE trcstp, ONLY : trc_stp    ! passive tracer time-stepping      (trc_stp routine) 
    110110#endif 
    111111   !!---------------------------------------------------------------------- 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcadv.F90

    r10874 r10880  
    6868CONTAINS 
    6969 
    70    SUBROUTINE trc_adv( kt ) 
     70   SUBROUTINE trc_adv( kt, Kbb, Kmm, ptr, Krhs ) 
    7171      !!---------------------------------------------------------------------- 
    7272      !!                  ***  ROUTINE trc_adv  *** 
     
    7676      !! ** Method  : - Update after tracers (tra) with the advection term following nadv 
    7777      !!---------------------------------------------------------------------- 
    78       INTEGER, INTENT(in) ::   kt   ! ocean time-step index 
     78      INTEGER                                   , INTENT(in)    :: kt   ! ocean time-step index 
     79      INTEGER                                   , INTENT(in)    :: Kbb, Kmm, Krhs ! time level indices 
     80      REAL(wp), DIMENSION(jpi,jpj,jpk,jptra,jpt), INTENT(inout) :: ptr            ! passive tracers and RHS of tracer equation 
    7981      ! 
    8082      INTEGER ::   jk   ! dummy loop index 
     
    123125      ! 
    124126      CASE ( np_CEN )                                 ! Centered : 2nd / 4th order 
    125          CALL tra_adv_cen( kt, nittrc000,'TRC',          zun, zvn, zwn     , trn, tra, jptra, nn_cen_h, nn_cen_v ) 
     127         CALL tra_adv_cen( kt, nittrc000,'TRC',          zun, zvn, zwn,      Kmm, ptr, jptra, Krhs, nn_cen_h, nn_cen_v ) 
    126128      CASE ( np_FCT )                                 ! FCT      : 2nd / 4th order 
    127          CALL tra_adv_fct( kt, nittrc000,'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra, nn_fct_h, nn_fct_v ) 
     129         CALL tra_adv_fct( kt, nittrc000,'TRC', r2dttrc, zun, zvn, zwn, Kbb, Kmm, ptr, jptra, Krhs, nn_fct_h, nn_fct_v ) 
    128130      CASE ( np_MUS )                                 ! MUSCL 
    129          CALL tra_adv_mus( kt, nittrc000,'TRC', r2dttrc, zun, zvn, zwn, trb,      tra, jptra        , ln_mus_ups )  
     131         CALL tra_adv_mus( kt, nittrc000,'TRC', r2dttrc, zun, zvn, zwn, Kbb, Kmm, ptr, jptra, Krhs, ln_mus_ups        )  
    130132      CASE ( np_UBS )                                 ! UBS 
    131          CALL tra_adv_ubs( kt, nittrc000,'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra          , nn_ubs_v ) 
     133         CALL tra_adv_ubs( kt, nittrc000,'TRC', r2dttrc, zun, zvn, zwn, Kbb, Kmm, ptr, jptra, Krhs, nn_ubs_v          ) 
    132134      CASE ( np_QCK )                                 ! QUICKEST 
    133          CALL tra_adv_qck( kt, nittrc000,'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra                     ) 
     135         CALL tra_adv_qck( kt, nittrc000,'TRC', r2dttrc, zun, zvn, zwn, Kbb, Kmm, ptr, jptra, Krhs                     ) 
    134136      ! 
    135137      END SELECT 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trctrp.F90

    r10874 r10880  
    4444CONTAINS 
    4545 
    46    SUBROUTINE trc_trp( kt ) 
     46   SUBROUTINE trc_trp( kt, Kbb, Kmm, Krhs, Kaa ) 
    4747      !!---------------------------------------------------------------------- 
    4848      !!                     ***  ROUTINE trc_trp  *** 
     
    5353      !!              - Update the passive tracers 
    5454      !!---------------------------------------------------------------------- 
    55       INTEGER, INTENT( in ) ::  kt  ! ocean time-step index 
     55      INTEGER, INTENT( in ) :: kt                  ! ocean time-step index 
     56      INTEGER, INTENT( in ) :: Kbb, Kmm, Krhs, Kaa ! time level indices 
    5657      !! --------------------------------------------------------------------- 
    5758      ! 
     
    6465         IF( ln_trcdmp )        CALL trc_dmp    ( kt )      ! internal damping trends 
    6566         IF( ln_bdy )           CALL trc_bdy_dmp( kt )      ! BDY damping trends 
    66                                 CALL trc_adv    ( kt )      ! horizontal & vertical advection  
     67                                CALL trc_adv    ( kt, Kbb, Kmm, tr, Krhs )      ! horizontal & vertical advection  
    6768         !                                                         ! Partial top/bottom cell: GRADh( trb )   
    6869         IF( ln_zps ) THEN 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/oce_trc.F90

    r10351 r10880  
    88   !!---------------------------------------------------------------------- 
    99   !                                            !* Domain size * 
     10   USE par_oce , ONLY :   jpt      =>   jpt        !: time dimension 
    1011   USE par_oce , ONLY :   jpi      =>   jpi        !: first  dimension of grid --> i  
    1112   USE par_oce , ONLY :   jpj      =>   jpj        !: second dimension of grid --> j   
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/trc.F90

    r10425 r10880  
    3333   REAL(wp), PUBLIC                                        ::  areatot        !: total volume  
    3434   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:  ) ::  cvol           !: volume correction -degrad option-  
    35    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::  trn            !: tracer concentration for now time step 
    36    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::  tra            !: tracer concentration for next time step 
    37    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::  trb            !: tracer concentration for before time step 
     35   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:,:), TARGET ::  tr             !: tracer concentration  
    3836   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:  ) ::  sbc_trc_b      !: Before sbc fluxes for tracers 
    3937   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:  ) ::  sbc_trc        !: Now sbc fluxes for tracers 
     
    4240   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:  ) ::  trc_o          !: prescribed tracer concentration in ocean for SBC 
    4341   INTEGER             , PUBLIC                            ::  nn_ice_tr      !: handling of sea ice tracers 
     42 
     43   !! TEMPORARY POINTERS - TO BE DELETED AFTER IMMERSE DEVELOPMENT COMPLETE 
     44   REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:,:,:)   ::  trn            !: tracer concentration for now time step 
     45   REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:,:,:)   ::  tra            !: tracer concentration for next time step 
     46   REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:,:,:)   ::  trb            !: tracer concentration for before time step 
     47   !! TEMPORARY POINTERS - TO BE DELETED AFTER IMMERSE DEVELOPMENT COMPLETE 
    4448 
    4549   !! interpolated gradient 
     
    147151      ierr(:) = 0 
    148152      ! 
    149       ALLOCATE( trn(jpi,jpj,jpk,jptra), trb(jpi,jpj,jpk,jptra), tra(jpi,jpj,jpk,jptra),       &   
     153      ALLOCATE( tr(jpi,jpj,jpk,jptra,jpt)                                             ,       &   
    150154         &      trc_i(jpi,jpj,jptra)  , trc_o(jpi,jpj,jptra)                          ,       & 
    151155         &      gtru (jpi,jpj,jptra)  , gtrv (jpi,jpj,jptra)                          ,       & 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/trcini.F90

    r10570 r10880  
    2626   USE trcice          ! tracers in sea ice 
    2727   USE trcbc,   only : trc_bc_ini ! generalized Boundary Conditions 
     28   USE trcstp          ! for time level indices (to be initialised) 
    2829  
    2930   IMPLICIT NONE 
     
    6162      CALL trc_nam       ! read passive tracers namelists 
    6263      CALL top_alloc()   ! allocate TOP arrays 
     64 
     65      ! Initialise time level indices 
     66      Nbb = 1; Nnn = 2; Naa = 3; Nrhs = Naa 
     67 
     68      ! Initialisation of temporary pointers (to be deleted after development finished) 
     69      CALL update_pointers_trc() 
    6370      ! 
    6471      IF(.NOT.ln_trcdta )   ln_trc_ini(:) = .FALSE. 
  • NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/trcstp.F90

    r10570 r10880  
    3030 
    3131   PUBLIC   trc_stp    ! called by step 
     32   PUBLIC   update_pointers_trc ! called in initialisation 
     33 
     34   !!---------------------------------------------------------------------- 
     35   !! time level indices 
     36   !!---------------------------------------------------------------------- 
     37   INTEGER, PUBLIC :: Nbb, Nnn, Naa, Nrhs               !! used by trc_init 
    3238 
    3339   LOGICAL  ::   llnew                   ! ??? 
     
    100106                                   CALL trc_wri      ( kt )       ! output of passive tracers with iom I/O manager 
    101107                                   CALL trc_sms      ( kt )       ! tracers: sinks and sources 
    102                                    CALL trc_trp      ( kt )       ! transport of passive tracers 
     108                                   CALL trc_trp      ( kt, Nbb, Nnn, Nrhs, Naa )       ! transport of passive tracers 
    103109         IF( kt == nittrc000 ) THEN 
    104110            CALL iom_close( numrtr )       ! close input tracer restart file 
     
    125131   END SUBROUTINE trc_stp 
    126132 
     133   SUBROUTINE update_pointers_trc 
     134      !!---------------------------------------------------------------------- 
     135      !!                     ***  ROUTINE update_pointers_trc  *** 
     136      !! 
     137      !! ** Purpose :   Associate temporary pointer arrays. 
     138      !!                For IMMERSE development phase only - to be deleted 
     139      !! 
     140      !! ** Method  : 
     141      !!---------------------------------------------------------------------- 
     142 
     143      trb => tr(:,:,:,:,Nbb); trn => tr(:,:,:,:,Nnn); tra => tr(:,:,:,:,Naa) 
     144 
     145   END SUBROUTINE update_pointers_trc 
    127146 
    128147   SUBROUTINE trc_mean_qsr( kt ) 
Note: See TracChangeset for help on using the changeset viewer.