New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 15188 for NEMO/branches/2021 – NEMO

Changeset 15188 for NEMO/branches/2021


Ignore:
Timestamp:
2021-08-13T14:47:00+02:00 (3 years ago)
Author:
techene
Message:

#2715 add advective velocities

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2021/dev_r14318_RK3_stage1/src/TOP/TRP/trcadv.F90

    r14086 r15188  
    88   !!            3.7  !  2014-05  (G. Madec, C. Ethe)  Add 2nd/4th order cases for CEN and FCT schemes  
    99   !!            4.0  !  2017-09  (G. Madec)  remove vertical time-splitting option 
     10   !!            4.5  !  2021-08  (G. Madec, S. Techene) add advective velocities as optional arguments 
    1011   !!---------------------------------------------------------------------- 
    1112#if defined key_top 
     
    3839   PRIVATE 
    3940 
    40    PUBLIC   trc_adv       ! called by trctrp.F90 
     41   PUBLIC   trc_adv       ! called by trctrp.F90 and stprk3_stg.F90 
    4142   PUBLIC   trc_adv_ini   ! called by trcini.F90 
    4243 
     
    6162   INTEGER, PARAMETER ::   np_UBS     = 4   ! 3rd order Upstream Biased Scheme 
    6263   INTEGER, PARAMETER ::   np_QCK     = 5   ! QUICK scheme 
    63     
     64 
     65   !! * Substitutions 
     66#  include "do_loop_substitute.h90" 
    6467#  include "domzgr_substitute.h90" 
    6568   !!---------------------------------------------------------------------- 
     
    7073CONTAINS 
    7174 
    72    SUBROUTINE trc_adv( kt, Kbb, Kmm, ptr, Krhs  ) 
     75   SUBROUTINE trc_adv( kt, Kbb, Kmm, ptr, Krhs, pau, pav, paw ) 
    7376      !!---------------------------------------------------------------------- 
    7477      !!                  ***  ROUTINE trc_adv  *** 
     
    7679      !! ** Purpose :   compute the ocean tracer advection trend. 
    7780      !! 
    78       !! ** Method  : - Update after tracers (tr(Krhs)) with the advection term following nadv 
    79       !!---------------------------------------------------------------------- 
    80       INTEGER                                   , INTENT(in)    :: kt   ! ocean time-step index 
    81       INTEGER                                   , INTENT(in)    :: Kbb, Kmm, Krhs ! time level indices 
    82       REAL(wp), DIMENSION(jpi,jpj,jpk,jptra,jpt), INTENT(inout) :: ptr            ! passive tracers and RHS of tracer equation 
    83       ! 
    84       INTEGER ::   jk   ! dummy loop index 
     81      !! ** Method  : - Update tr(Krhs) with the advective trend following nadv 
     82      !!---------------------------------------------------------------------- 
     83      INTEGER                                     , INTENT(in   ) ::   kt             ! ocean time-step index 
     84      INTEGER                                     , INTENT(in   ) ::   Kbb, Kmm, Krhs ! time level indices 
     85      REAL(wp), DIMENSION(:,:,:), OPTIONAL, TARGET, INTENT(in   ) ::   pau, pav, paw  ! advective velocity 
     86      REAL(wp), DIMENSION(jpi,jpj,jpk,jptra,jpt)  , INTENT(inout) ::   ptr            ! passive tracers and RHS of tracer equation 
     87      ! 
     88      INTEGER ::  ji, jj, jk   ! dummy loop index 
    8589      CHARACTER (len=22) ::   charout 
     90      REAL(wp), DIMENSION(:,:,:), POINTER ::   zptu, zptv, zptw 
    8691      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zuu, zvv, zww  ! effective velocity 
    8792      !!---------------------------------------------------------------------- 
     
    9196      !                                         !==  effective transport  ==! 
    9297      IF( l_offline ) THEN 
    93          zuu(:,:,:) = uu(:,:,:,Kmm)                ! already in (uu(Kmm),vv(Kmm),ww) 
     98         zuu(:,:,:) = uu(:,:,:,Kmm)                != already in (uu(Kmm),vv(Kmm),ww) 
    9499         zvv(:,:,:) = vv(:,:,:,Kmm) 
    95100         zww(:,:,:) = ww(:,:,:) 
    96       ELSE                                         ! build the effective transport 
    97          zuu(:,:,jpk) = 0._wp 
     101      ELSE                                         != build the effective transport 
     102         zuu(:,:,jpk) = 0._wp                            ! no transport trough the bottom 
    98103         zvv(:,:,jpk) = 0._wp 
    99104         zww(:,:,jpk) = 0._wp 
    100          IF( ln_wave .AND. ln_sdw )  THEN 
    101             DO jk = 1, jpkm1                                                       ! eulerian transport + Stokes Drift 
    102                zuu(:,:,jk) = e2u  (:,:) * e3u(:,:,jk,Kmm) * ( uu(:,:,jk,Kmm) + usd(:,:,jk) ) 
    103                zvv(:,:,jk) = e1v  (:,:) * e3v(:,:,jk,Kmm) * ( vv(:,:,jk,Kmm) + vsd(:,:,jk) ) 
    104                zww(:,:,jk) = e1e2t(:,:)                   * ( ww(:,:,jk) + wsd(:,:,jk) ) 
    105             END DO 
    106          ELSE 
    107             DO jk = 1, jpkm1 
    108                zuu(:,:,jk) = e2u  (:,:) * e3u(:,:,jk,Kmm) * uu(:,:,jk,Kmm)                   ! eulerian transport 
    109                zvv(:,:,jk) = e1v  (:,:) * e3v(:,:,jk,Kmm) * vv(:,:,jk,Kmm) 
    110                zww(:,:,jk) = e1e2t(:,:)                   * ww(:,:,jk) 
    111             END DO 
     105         ! 
     106         IF( PRESENT( pau ) ) THEN                       ! RK3: advective velocity (pau,pav,paw) /= advected velocity (uu,vv,ww) 
     107            zptu => pau(:,:,:) 
     108            zptv => pav(:,:,:) 
     109            zptw => paw(:,:,:) 
     110         ELSE                                            ! MLF: advective velocity = (uu,vv,ww) 
     111            zptu => uu(:,:,:,Kmm) 
     112            zptv => vv(:,:,:,Kmm) 
     113            zptw => ww(:,:,:    ) 
    112114         ENDIF 
    113115         ! 
    114          IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN                                 ! add z-tilde and/or vvl corrections 
    115             zuu(:,:,:) = zuu(:,:,:) + un_td(:,:,:) 
    116             zvv(:,:,:) = zvv(:,:,:) + vn_td(:,:,:) 
     116         IF( ln_wave .AND. ln_sdw )  THEN                ! eulerian transport + Stokes Drift 
     117            DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) 
     118               zuu(ji,jj,jk) = e2u  (ji,jj) * e3u(ji,jj,jk,Kmm) * ( zptu(ji,jj,jk) + usd(ji,jj,jk) ) 
     119               zvv(ji,jj,jk) = e1v  (ji,jj) * e3v(ji,jj,jk,Kmm) * ( zptv(ji,jj,jk) + vsd(ji,jj,jk) ) 
     120               zww(ji,jj,jk) = e1e2t(ji,jj)                     * ( zptw(ji,jj,jk) + wsd(ji,jj,jk) ) 
     121            END_3D 
     122         ELSE                                            ! eulerian transport only 
     123            DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) 
     124               zuu(ji,jj,jk) = e2u  (ji,jj) * e3u(ji,jj,jk,Kmm) * zptu(ji,jj,jk) 
     125               zvv(ji,jj,jk) = e1v  (ji,jj) * e3v(ji,jj,jk,Kmm) * zptv(ji,jj,jk) 
     126               zww(ji,jj,jk) = e1e2t(ji,jj)                     * zptw(ji,jj,jk) 
     127            END_3D 
     128         ENDIF 
     129         ! 
     130         IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN      ! add z-tilde and/or vvl corrections 
     131            DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 ) 
     132               zuu(ji,jj,jk) = zuu(ji,jj,jk) + un_td(ji,jj,jk) 
     133               zvv(ji,jj,jk) = zvv(ji,jj,jk) + vn_td(ji,jj,jk) 
     134            END_3D 
    117135         ENDIF 
    118136         ! 
     
    127145      ! 
    128146      CASE ( np_CEN )                                 ! Centered : 2nd / 4th order 
    129          IF (nn_hls.EQ.2) CALL lbc_lnk( 'trcadv', ptr(:,:,:,:,Kmm), 'T', 1.) 
     147         IF( nn_hls == 2 ) CALL lbc_lnk( 'trcadv', ptr(:,:,:,:,Kmm), 'T', 1._wp ) 
    130148         CALL tra_adv_cen( kt, nittrc000,'TRC',          zuu, zvv, zww,      Kmm, ptr, jptra, Krhs, nn_cen_h, nn_cen_v ) 
    131149      CASE ( np_FCT )                                 ! FCT      : 2nd / 4th order 
    132          IF (nn_hls.EQ.2) THEN 
    133             CALL lbc_lnk_multi( 'trcadv', ptr(:,:,:,:,Kbb), 'T', 1., ptr(:,:,:,:,Kmm), 'T', 1.) 
    134             CALL lbc_lnk_multi( 'traadv', zuu(:,:,:), 'U', -1., zvv(:,:,:), 'V', -1., zww(:,:,:), 'W', 1.) 
     150         IF( nn_hls == 2 ) THEN 
     151            CALL lbc_lnk_multi( 'trcadv', ptr(:,:,:,:,Kbb), 'T', 1._wp, ptr(:,:,:,:,Kmm), 'T', 1._wp) 
     152            CALL lbc_lnk_multi( 'traadv', zuu(:,:,:), 'U', -1._wp, zvv(:,:,:), 'V', -1._wp, zww(:,:,:), 'W', 1._wp) 
    135153#if defined key_loop_fusion 
    136154            CALL tra_adv_fct_lf( kt, nittrc000,'TRC', rDt_trc, zuu, zvv, zww, Kbb, Kmm, ptr, jptra, Krhs, nn_fct_h, nn_fct_v ) 
     
    142160         END IF 
    143161      CASE ( np_MUS )                                 ! MUSCL 
    144          IF (nn_hls.EQ.2) THEN 
    145             IF (nn_hls.EQ.2) CALL lbc_lnk( 'trcadv', ptr(:,:,:,:,Kbb), 'T', 1.) 
     162         IF( nn_hls == 2 ) THEN 
     163            CALL lbc_lnk( 'trcadv', ptr(:,:,:,:,Kbb), 'T', 1._wp) 
    146164#if defined key_loop_fusion 
    147165            CALL tra_adv_mus_lf( kt, nittrc000,'TRC', rDt_trc, zuu, zvv, zww, Kbb, Kmm, ptr, jptra, Krhs, ln_mus_ups )  
     
    153171         END IF 
    154172      CASE ( np_UBS )                                 ! UBS 
    155          IF (nn_hls.EQ.2) CALL lbc_lnk( 'trcadv', ptr(:,:,:,:,Kbb), 'T', 1.) 
     173         IF( nn_hls == 2 )   CALL lbc_lnk( 'trcadv', ptr(:,:,:,:,Kbb), 'T', 1._wp) 
    156174         CALL tra_adv_ubs( kt, nittrc000,'TRC', rDt_trc, zuu, zvv, zww, Kbb, Kmm, ptr, jptra, Krhs, nn_ubs_v           ) 
    157175      CASE ( np_QCK )                                 ! QUICKEST 
    158          IF (nn_hls.EQ.2) THEN 
    159             CALL lbc_lnk_multi( 'trcadv', zuu(:,:,:), 'U', -1., zvv(:,:,:), 'V', -1.) 
    160             CALL lbc_lnk( 'traadv', ptr(:,:,:,:,Kbb), 'T', 1.) 
     176         IF( nn_hls == 2 ) THEN 
     177            CALL lbc_lnk_multi( 'trcadv', zuu(:,:,:), 'U', -1._wp, zvv(:,:,:), 'V', -1._wp) 
     178            CALL lbc_lnk( 'traadv', ptr(:,:,:,:,Kbb), 'T', 1._wp) 
    161179         END IF 
    162180         CALL tra_adv_qck( kt, nittrc000,'TRC', rDt_trc, zuu, zvv, zww, Kbb, Kmm, ptr, jptra, Krhs                     ) 
Note: See TracChangeset for help on using the changeset viewer.