MODULE trcadv_crs !!============================================================================== !! *** MODULE trcadv *** !! Ocean passive tracers: advection trend !!============================================================================== !! History : 2.0 ! 05-11 (G. Madec) Original code !! 3.0 ! 10-06 (C. Ethe) Adapted to passive tracers !!---------------------------------------------------------------------- #if defined key_top !!---------------------------------------------------------------------- !! 'key_top' TOP models !!---------------------------------------------------------------------- !! trc_adv : compute ocean tracer advection trend !! trc_adv_ctl : control the different options of advection scheme !!---------------------------------------------------------------------- USE oce_trc ! ocean dynamics and active tracers !???USE oce_trc, ONLY: un,vn,wn USE trc ! ocean passive tracers variables USE trcnam_trp ! passive tracers transport namelist variables USE traadv_tvd_crs ! TVD scheme (tra_adv_tvd routine) USE traadv_ubs_crs ! TVD scheme (tra_adv_tvd routine) USE ldftra_oce ! lateral diffusion coefficient on tracers USE prtctl_trc ! Print control USE crs , ONLY : e2e3u_msk , e1e3v_msk , e1e2w_msk,jpi_crs,jpj_crs USE timing USE iom, ONLY: iom_put,iom_swap IMPLICIT NONE PRIVATE PUBLIC trc_adv_crs ! routine called by step module PUBLIC trc_adv_alloc_crs ! routine called by nemogcm module INTEGER :: nadv ! choice of the type of advection scheme REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: r2dt ! vertical profile time-step, = 2 rdttra ! ! except at nitrrc000 (=rdttra) if neuler=0 !! * Substitutions # include "domzgr_substitute.h90" # include "vectopt_loop_substitute.h90" !!---------------------------------------------------------------------- !! NEMO/TOP 3.3 , NEMO Consortium (2010) !! $Id$ !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) !!---------------------------------------------------------------------- CONTAINS INTEGER FUNCTION trc_adv_alloc_crs() !!---------------------------------------------------------------------- !! *** ROUTINE trc_adv_alloc *** !!---------------------------------------------------------------------- ALLOCATE( r2dt(jpk), STAT=trc_adv_alloc_crs ) IF( trc_adv_alloc_crs /= 0 ) CALL ctl_warn('trc_adv_alloc : failed to allocate array.') END FUNCTION trc_adv_alloc_crs SUBROUTINE trc_adv_crs( kt ) !!---------------------------------------------------------------------- !! *** ROUTINE trc_adv *** !! !! ** Purpose : compute the ocean tracer advection trend. !! !! ** Method : - Update the tracer with the advection term following nadv !!---------------------------------------------------------------------- !! INTEGER, INTENT(in) :: kt ! ocean time-step index ! INTEGER :: jk INTEGER :: ji,jj CHARACTER (len=22) :: charout REAL(wp), POINTER, DIMENSION(:,:,:) :: zun, zvn, zwn ! effective velocity !!---------------------------------------------------------------------- ! IF( nn_timing == 1 ) CALL timing_start('trc_adv') ! CALL wrk_alloc( jpi, jpj, jpk, zun, zvn, zwn ) ! IF( kt == nittrc000 ) CALL trc_adv_ctl_crs ! initialisation & control of options #if ! defined key_pisces IF( neuler == 0 .AND. kt == nittrc000 ) THEN ! at nittrc000 r2dt(:) = rdttrc(:) ! = rdttrc (restarting with Euler time stepping) ELSEIF( kt <= nittrc000 + 1 ) THEN ! at nittrc000 or nittrc000+1 r2dt(:) = 2. * rdttrc(:) ! = 2 rdttrc (leapfrog) ENDIF #else r2dt(:) = rdttrc(:) ! = rdttrc (for PISCES use Euler time stepping) #endif DO jk = 1, jpkm1 ! ! eulerian transport only zun(:,:,jk) = e2e3u_msk(:,:,jk) * un(:,:,jk) zvn(:,:,jk) = e1e3v_msk(:,:,jk) * vn(:,:,jk) zwn(:,:,jk) = e1e2w_msk(:,:,jk) * wn(:,:,jk) ! END DO zwn(:,:,jpk) = 0.e0 ! no transport trough the bottom ! SELECT CASE ( nadv ) !== compute advection trend and add it to general trend ==! CASE ( 2 ) ; CALL tra_adv_tvd_crs ( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb, trn, tra, jptra ) ! TVD CASE ( 5 ) ; CALL tra_adv_ubs_crs ( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb, trn, tra, jptra ) ! UBS ! CASE (-1 ) !== esopa: test all possibility with control print ==! CALL tra_adv_tvd_crs ( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb, trn, tra, jptra ) WRITE(charout, FMT="('adv2')") ; CALL prt_ctl_trc_info(charout) CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd') ! END SELECT ! ! print mean trends (used for debugging) IF( ln_ctl ) THEN WRITE(charout, FMT="('adv ')") ; CALL prt_ctl_trc_info(charout) CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) END IF ! CALL wrk_dealloc( jpi, jpj, jpk, zun, zvn, zwn ) ! IF( nn_timing == 1 ) CALL timing_stop('trc_adv') ! END SUBROUTINE trc_adv_crs SUBROUTINE trc_adv_ctl_crs !!--------------------------------------------------------------------- !! *** ROUTINE trc_adv_ctl *** !! !! ** Purpose : Control the consistency between namelist options for !! passive tracer advection schemes and set nadv !!---------------------------------------------------------------------- INTEGER :: ioptio !!---------------------------------------------------------------------- ioptio = 0 ! Parameter control IF( ln_trcadv_cen2 ) ioptio = ioptio + 1 IF( ln_trcadv_tvd ) ioptio = ioptio + 1 IF( ln_trcadv_muscl ) ioptio = ioptio + 1 IF( ln_trcadv_muscl2 ) ioptio = ioptio + 1 IF( ln_trcadv_ubs ) ioptio = ioptio + 1 IF( ln_trcadv_qck ) ioptio = ioptio + 1 IF( lk_esopa ) ioptio = 1 IF( ioptio /= 1 ) CALL ctl_stop( 'Choose ONE advection scheme in namelist namtrc_adv' ) ! ! Set nadv IF( ln_trcadv_cen2 ) nadv = 1 IF( ln_trcadv_tvd ) nadv = 2 IF( ln_trcadv_muscl ) nadv = 3 IF( ln_trcadv_muscl2 ) nadv = 4 IF( ln_trcadv_ubs ) nadv = 5 IF( ln_trcadv_qck ) nadv = 6 IF( lk_esopa ) nadv = -1 IF(lwp) THEN ! Print the choice WRITE(numout,*) IF( nadv == 1 ) WRITE(numout,*) ' 2nd order scheme is used' IF( nadv == 2 ) WRITE(numout,*) ' TVD scheme is used' IF( nadv == 3 ) WRITE(numout,*) ' MUSCL scheme is used' IF( nadv == 4 ) WRITE(numout,*) ' MUSCL2 scheme is used' IF( nadv == 5 ) WRITE(numout,*) ' UBS scheme is used' IF( nadv == 6 ) WRITE(numout,*) ' QUICKEST scheme is used' IF( nadv == -1 ) WRITE(numout,*) ' esopa test: use all advection scheme' ENDIF ! END SUBROUTINE trc_adv_ctl_crs #else !!---------------------------------------------------------------------- !! Default option Empty module !!---------------------------------------------------------------------- CONTAINS SUBROUTINE trc_adv_crs( kt ) INTEGER, INTENT(in) :: kt WRITE(*,*) 'trc_adv: You should not have seen this print! error?', kt END SUBROUTINE trc_adv_crs #endif !!====================================================================== END MODULE trcadv_crs