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 6225 for branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/TOP_SRC/TRP – NEMO

Ignore:
Timestamp:
2016-01-08T10:35:19+01:00 (8 years ago)
Author:
jamesharle
Message:

Update MPP_BDY_UPDATE branch to be consistent with head of trunk

Location:
branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/TOP_SRC/TRP
Files:
5 deleted
9 edited
4 copied

Legend:

Unmodified
Added
Removed

  • branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/TOP_SRC/TRP/trcadv.F90

    r4610 r6225  
    44   !! Ocean passive tracers:  advection trend  
    55   !!============================================================================== 
    6    !! History :  2.0  !  05-11  (G. Madec)  Original code 
    7    !!            3.0  !  10-06  (C. Ethe)   Adapted to passive tracers 
     6   !! History :  2.0  !  2005-11  (G. Madec)  Original code 
     7   !!            3.0  !  2010-06  (C. Ethe)   Adapted to passive tracers 
     8   !!            3.7  !  2014-05  (G. Madec, C. Ethe)  Add 2nd/4th order cases for CEN and FCT schemes  
    89   !!---------------------------------------------------------------------- 
    910#if defined key_top 
    … …  
    1112   !!   'key_top'                                                TOP models 
    1213   !!---------------------------------------------------------------------- 
    13    !!   trc_adv      : compute ocean tracer advection trend 
    14    !!   trc_adv_ctl  : control the different options of advection scheme 
    15    !!---------------------------------------------------------------------- 
    16    USE oce_trc         ! ocean dynamics and active tracers 
    17    USE trc             ! ocean passive tracers variables 
    18    USE trcnam_trp      ! passive tracers transport namelist variables 
    19    USE traadv_cen2     ! 2nd order centered scheme (tra_adv_cen2   routine) 
    20    USE traadv_tvd      ! TVD      scheme           (tra_adv_tvd    routine) 
    21    USE traadv_muscl    ! MUSCL    scheme           (tra_adv_muscl  routine) 
    22    USE traadv_muscl2   ! MUSCL2   scheme           (tra_adv_muscl2 routine) 
    23    USE traadv_ubs      ! UBS      scheme           (tra_adv_ubs    routine) 
    24    USE traadv_qck      ! QUICKEST scheme           (tra_adv_qck    routine) 
    25    USE traadv_eiv      ! eddy induced velocity     (tra_adv_eiv    routine) 
    26    USE traadv_mle      ! ML eddy induced velocity  (tra_adv_mle    routine) 
    27    USE ldftra_oce      ! lateral diffusion coefficient on tracers 
    28    USE prtctl_trc      ! Print control 
     14   !!   trc_adv       : compute ocean tracer advection trend 
     15   !!   trc_adv_ini   : control the different options of advection scheme 
     16   !!---------------------------------------------------------------------- 
     17   USE oce_trc        ! ocean dynamics and active tracers 
     18   USE trc            ! ocean passive tracers variables 
     19   USE traadv_cen     ! centered scheme           (tra_adv_cen  routine) 
     20   USE traadv_fct     ! FCT      scheme           (tra_adv_fct  routine) 
     21   USE traadv_mus     ! MUSCL    scheme           (tra_adv_mus  routine) 
     22   USE traadv_ubs     ! UBS      scheme           (tra_adv_ubs  routine) 
     23   USE traadv_qck     ! QUICKEST scheme           (tra_adv_qck  routine) 
     24   USE traadv_mle     ! ML eddy induced velocity  (tra_adv_mle  routine) 
     25   USE ldftra         ! lateral diffusion coefficient on tracers 
     26   USE ldfslp         ! Lateral diffusion: slopes of neutral surfaces 
     27   ! 
     28   USE prtctl_trc     ! Print control 
    2929 
    3030   IMPLICIT NONE 
    3131   PRIVATE 
    3232 
    33    PUBLIC   trc_adv          ! routine called by step module 
    34    PUBLIC   trc_adv_alloc    ! routine called by nemogcm module 
    35  
    36    INTEGER ::   nadv   ! choice of the type of advection scheme 
    37    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) ::   r2dt  ! vertical profile time-step, = 2 rdttra 
    38    !                                                    ! except at nitrrc000 (=rdttra) if neuler=0 
     33   PUBLIC   trc_adv        
     34   PUBLIC   trc_adv_ini   
     35 
     36   !                            !!* Namelist namtrc_adv * 
     37   LOGICAL ::   ln_trcadv_cen    ! centered scheme flag 
     38   INTEGER ::      nn_cen_h, nn_cen_v   ! =2/4 : horizontal and vertical choices of the order of CEN scheme 
     39   LOGICAL ::   ln_trcadv_fct    ! FCT scheme flag 
     40   INTEGER ::      nn_fct_h, nn_fct_v   ! =2/4 : horizontal and vertical choices of the order of FCT scheme 
     41   INTEGER ::      nn_fct_zts           ! >=1 : 2nd order FCT with vertical sub-timestepping 
     42   LOGICAL ::   ln_trcadv_mus    ! MUSCL scheme flag 
     43   LOGICAL ::      ln_mus_ups           ! use upstream scheme in vivcinity of river mouths 
     44   LOGICAL ::   ln_trcadv_ubs    ! UBS scheme flag 
     45   INTEGER ::      nn_ubs_v             ! =2/4 : vertical choice of the order of UBS scheme 
     46   LOGICAL ::   ln_trcadv_qck    ! QUICKEST scheme flag 
     47 
     48   !                                        ! choices of advection scheme: 
     49   INTEGER, PARAMETER ::   np_NO_adv  = 0   ! no T-S advection 
     50   INTEGER, PARAMETER ::   np_CEN     = 1   ! 2nd/4th order centered scheme 
     51   INTEGER, PARAMETER ::   np_FCT     = 2   ! 2nd/4th order Flux Corrected Transport scheme 
     52   INTEGER, PARAMETER ::   np_FCT_zts = 3   ! 2nd order FCT scheme with vertical sub-timestepping 
     53   INTEGER, PARAMETER ::   np_MUS     = 4   ! MUSCL scheme 
     54   INTEGER, PARAMETER ::   np_UBS     = 5   ! 3rd order Upstream Biased Scheme 
     55   INTEGER, PARAMETER ::   np_QCK     = 6   ! QUICK scheme 
     56 
     57   INTEGER ::              nadv             ! chosen advection scheme 
     58   ! 
     59   REAL(wp) ::   r2dttrc  ! vertical profile time-step, = 2 rdt 
     60   !                                                    ! except at nitrrc000 (=rdt) if neuler=0 
    3961 
    4062   !! * Substitutions 
    41 #  include "domzgr_substitute.h90" 
    4263#  include "vectopt_loop_substitute.h90" 
    4364   !!---------------------------------------------------------------------- 
    44    !! NEMO/TOP 3.3 , NEMO Consortium (2010) 
     65   !! NEMO/TOP 3.7 , NEMO Consortium (2015) 
    4566   !! $Id$  
    4667   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    4768   !!---------------------------------------------------------------------- 
    4869CONTAINS 
    49  
    50    INTEGER FUNCTION trc_adv_alloc() 
    51       !!---------------------------------------------------------------------- 
    52       !!                  ***  ROUTINE trc_adv_alloc  *** 
    53       !!---------------------------------------------------------------------- 
    54  
    55       ALLOCATE( r2dt(jpk), STAT=trc_adv_alloc ) 
    56  
    57       IF( trc_adv_alloc /= 0 ) CALL ctl_warn('trc_adv_alloc : failed to allocate array.') 
    58  
    59    END FUNCTION trc_adv_alloc 
    60  
    6170 
    6271   SUBROUTINE trc_adv( kt ) 
    … …  
    6877      !! ** Method  : - Update the tracer with the advection term following nadv 
    6978      !!---------------------------------------------------------------------- 
    70       !! 
    7179      INTEGER, INTENT(in) ::   kt   ! ocean time-step index 
    7280      ! 
    73       INTEGER ::   jk  
     81      INTEGER ::   jk   ! dummy loop index 
    7482      CHARACTER (len=22) ::   charout 
    7583      REAL(wp), POINTER, DIMENSION(:,:,:) :: zun, zvn, zwn  ! effective velocity 
    7684      !!---------------------------------------------------------------------- 
    7785      ! 
    78       IF( nn_timing == 1 )  CALL timing_start('trc_adv') 
    79       ! 
    80       CALL wrk_alloc( jpi, jpj, jpk, zun, zvn, zwn ) 
    81       ! 
    82  
    83       IF( kt == nittrc000 )   CALL trc_adv_ctl          ! initialisation & control of options 
    84  
    85       IF( ln_top_euler) THEN 
    86          r2dt(:) =  rdttrc(:)              ! = rdttrc (use Euler time stepping) 
    87       ELSE 
    88          IF( neuler == 0 .AND. kt == nittrc000 ) THEN     ! at nittrc000 
    89             r2dt(:) =  rdttrc(:)           ! = rdttrc (restarting with Euler time stepping) 
    90          ELSEIF( kt <= nittrc000 + 1 ) THEN          ! at nittrc000 or nittrc000+1 
    91             r2dt(:) = 2. * rdttrc(:)       ! = 2 rdttrc (leapfrog) 
    92          ENDIF 
    93       ENDIF 
    94  
    95       !                                                   ! effective transport 
     86      IF( nn_timing == 1 )   CALL timing_start('trc_adv') 
     87      ! 
     88      CALL wrk_alloc( jpi,jpj,jpk,   zun, zvn, zwn ) 
     89      ! 
     90      IF( ( neuler == 0 .AND. kt == nittrc000 ) .OR. ln_top_euler ) THEN     ! at nittrc000 
     91         r2dttrc =  rdttrc           ! = rdttrc (use or restarting with Euler time stepping) 
     92      ELSEIF( kt <= nittrc000 + nn_dttrc ) THEN          ! at nittrc000 or nittrc000+1 
     93         r2dttrc = 2. * rdttrc       ! = 2 rdttrc (leapfrog) 
     94      ENDIF 
     95      !                                               !==  effective transport  ==! 
    9696      DO jk = 1, jpkm1 
    97          !                                                ! eulerian transport only 
    98          zun(:,:,jk) = e2u  (:,:) * fse3u(:,:,jk) * un(:,:,jk) 
    99          zvn(:,:,jk) = e1v  (:,:) * fse3v(:,:,jk) * vn(:,:,jk) 
     97         zun(:,:,jk) = e2u  (:,:) * e3u_n(:,:,jk) * un(:,:,jk)                   ! eulerian transport 
     98         zvn(:,:,jk) = e1v  (:,:) * e3v_n(:,:,jk) * vn(:,:,jk) 
    10099         zwn(:,:,jk) = e1e2t(:,:)                 * wn(:,:,jk) 
    101          ! 
    102100      END DO 
    103101      ! 
    104       IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN 
     102      IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN                                 ! add z-tilde and/or vvl corrections 
    105103         zun(:,:,:) = zun(:,:,:) + un_td(:,:,:) 
    106104         zvn(:,:,:) = zvn(:,:,:) + vn_td(:,:,:) 
    107105      ENDIF 
    108106      ! 
    109       zun(:,:,jpk) = 0._wp                                                     ! no transport trough the bottom 
    110       zvn(:,:,jpk) = 0._wp                                                     ! no transport trough the bottom 
    111       zwn(:,:,jpk) = 0._wp                                                     ! no transport trough the bottom 
    112  
    113       IF( lk_traldf_eiv .AND. .NOT. ln_traldf_grif )   &  ! add the eiv transport (if necessary) 
    114          &              CALL tra_adv_eiv( kt, nittrc000, zun, zvn, zwn, 'TRC' ) 
    115       ! 
    116       IF( ln_mle    )   CALL tra_adv_mle( kt, nittrc000, zun, zvn, zwn, 'TRC' )    ! add the mle transport (if necessary) 
    117       ! 
    118       SELECT CASE ( nadv )                            !==  compute advection trend and add it to general trend  ==! 
    119       CASE ( 1 )   ;    CALL tra_adv_cen2  ( kt, nittrc000, 'TRC',       zun, zvn, zwn, trb, trn, tra, jptra )   !  2nd order centered 
    120       CASE ( 2 )   ;    CALL tra_adv_tvd   ( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb, trn, tra, jptra )   !  TVD  
    121       CASE ( 3 )   ;    CALL tra_adv_muscl ( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb,      tra, jptra, ln_trcadv_msc_ups )   !  MUSCL  
    122       CASE ( 4 )   ;    CALL tra_adv_muscl2( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb, trn, tra, jptra )   !  MUSCL2  
    123       CASE ( 5 )   ;    CALL tra_adv_ubs   ( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb, trn, tra, jptra )   !  UBS  
    124       CASE ( 6 )   ;    CALL tra_adv_qck   ( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb, trn, tra, jptra )   !  QUICKEST  
    125       ! 
    126       CASE (-1 )                                      !==  esopa: test all possibility with control print  ==! 
    127          CALL tra_adv_cen2  ( kt, nittrc000, 'TRC',       zun, zvn, zwn, trb, trn, tra, jptra )           
    128          WRITE(charout, FMT="('adv1')")  ; CALL prt_ctl_trc_info(charout) 
    129                                            CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd') 
    130          CALL tra_adv_tvd   ( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb, trn, tra, jptra )           
    131          WRITE(charout, FMT="('adv2')")  ; CALL prt_ctl_trc_info(charout) 
    132                                            CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd') 
    133          CALL tra_adv_muscl ( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb,      tra, jptra, ln_trcadv_msc_ups  )           
    134          WRITE(charout, FMT="('adv3')")  ; CALL prt_ctl_trc_info(charout) 
    135                                            CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd') 
    136          CALL tra_adv_muscl2( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb, trn, tra, jptra )           
    137          WRITE(charout, FMT="('adv4')")  ; CALL prt_ctl_trc_info(charout) 
    138                                            CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd') 
    139          CALL tra_adv_ubs   ( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb, trn, tra, jptra )           
    140          WRITE(charout, FMT="('adv5')")  ; CALL prt_ctl_trc_info(charout) 
    141                                            CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd') 
    142          CALL tra_adv_qck   ( kt, nittrc000, 'TRC', r2dt, zun, zvn, zwn, trb, trn, tra, jptra )           
    143          WRITE(charout, FMT="('adv6')")  ; CALL prt_ctl_trc_info(charout) 
    144                                            CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm,clinfo2='trd') 
    145          ! 
     107      IF( ln_ldfeiv .AND. .NOT. ln_traldf_triad )   &  
     108         &              CALL ldf_eiv_trp( kt, nittrc000, zun, zvn, zwn, 'TRC' )  ! add the eiv transport 
     109      ! 
     110      IF( ln_mle    )   CALL tra_adv_mle( kt, nittrc000, zun, zvn, zwn, 'TRC' )  ! add the mle transport 
     111      ! 
     112      zun(:,:,jpk) = 0._wp                                                       ! no transport trough the bottom 
     113      zvn(:,:,jpk) = 0._wp 
     114      zwn(:,:,jpk) = 0._wp 
     115      ! 
     116      ! 
     117      SELECT CASE ( nadv )                      !==  compute advection trend and add it to general trend  ==! 
     118      ! 
     119      CASE ( np_CEN )                                    ! Centered : 2nd / 4th order 
     120         CALL tra_adv_cen    ( kt, nittrc000,'TRC',       zun, zvn, zwn     , trn, tra, jptra, nn_cen_h, nn_cen_v ) 
     121      CASE ( np_FCT )                                    ! FCT      : 2nd / 4th order 
     122         CALL tra_adv_fct    ( kt, nittrc000,'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra, nn_fct_h, nn_fct_v ) 
     123      CASE ( np_FCT_zts )                                ! 2nd order FCT with vertical time-splitting 
     124         CALL tra_adv_fct_zts( kt, nittrc000,'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra        , nn_fct_zts ) 
     125      CASE ( np_MUS )                                    ! MUSCL 
     126         CALL tra_adv_mus    ( kt, nittrc000,'TRC', r2dttrc, zun, zvn, zwn, trb,      tra, jptra        , ln_mus_ups )  
     127      CASE ( np_UBS )                                    ! UBS 
     128         CALL tra_adv_ubs    ( kt, nittrc000,'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra        , nn_ubs_v   ) 
     129      CASE ( np_QCK )                                    ! QUICKEST 
     130         CALL tra_adv_qck    ( kt, nittrc000,'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra                     ) 
     131      ! 
    146132      END SELECT 
    147  
    148       !                                              ! print mean trends (used for debugging) 
    149       IF( ln_ctl )   THEN 
     133      !                   
     134      IF( ln_ctl )   THEN                             !== print mean trends (used for debugging) 
    150135         WRITE(charout, FMT="('adv ')")  ;  CALL prt_ctl_trc_info(charout) 
    151136                                            CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 
    152137      END IF 
    153138      ! 
    154       CALL wrk_dealloc( jpi, jpj, jpk, zun, zvn, zwn ) 
     139      CALL wrk_dealloc( jpi,jpj,jpk,  zun, zvn, zwn ) 
    155140      ! 
    156141      IF( nn_timing == 1 )  CALL timing_stop('trc_adv') 
    … …  
    159144 
    160145 
    161    SUBROUTINE trc_adv_ctl 
     146   SUBROUTINE trc_adv_ini 
    162147      !!--------------------------------------------------------------------- 
    163       !!                  ***  ROUTINE trc_adv_ctl  *** 
     148      !!                  ***  ROUTINE trc_adv_ini  *** 
    164149      !!                 
    165150      !! ** Purpose : Control the consistency between namelist options for  
    … …  
    167152      !!---------------------------------------------------------------------- 
    168153      INTEGER ::   ioptio 
    169       !!---------------------------------------------------------------------- 
    170  
    171       ioptio = 0                      ! Parameter control 
    172       IF( ln_trcadv_cen2   )   ioptio = ioptio + 1 
    173       IF( ln_trcadv_tvd    )   ioptio = ioptio + 1 
    174       IF( ln_trcadv_muscl  )   ioptio = ioptio + 1 
    175       IF( ln_trcadv_muscl2 )   ioptio = ioptio + 1 
    176       IF( ln_trcadv_ubs    )   ioptio = ioptio + 1 
    177       IF( ln_trcadv_qck    )   ioptio = ioptio + 1 
    178       IF( lk_esopa         )   ioptio =          1 
    179  
     154      INTEGER ::  ios                 ! Local integer output status for namelist read 
     155      !! 
     156      NAMELIST/namtrc_adv/ ln_trcadv_cen, nn_cen_h, nn_cen_v,               &   ! CEN 
     157         &                 ln_trcadv_fct, nn_fct_h, nn_fct_v, nn_fct_zts,   &   ! FCT 
     158         &                 ln_trcadv_mus,                     ln_mus_ups,   &   ! MUSCL 
     159         &                 ln_trcadv_ubs,           nn_ubs_v,               &   ! UBS 
     160         &                 ln_trcadv_qck                                        ! QCK 
     161      !!---------------------------------------------------------------------- 
     162      ! 
     163      REWIND( numnat_ref )              !  namtrc_adv in reference namelist  
     164      READ  ( numnat_ref, namtrc_adv, IOSTAT = ios, ERR = 901) 
     165901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrc_adv in reference namelist', lwp ) 
     166 
     167      REWIND( numnat_cfg )              ! namtrc_adv in configuration namelist 
     168      READ  ( numnat_cfg, namtrc_adv, IOSTAT = ios, ERR = 902 ) 
     169902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrc_adv in configuration namelist', lwp ) 
     170      IF(lwm) WRITE ( numont, namtrc_adv ) 
     171 
     172      IF(lwp) THEN                    ! Namelist print 
     173         WRITE(numout,*) 
     174         WRITE(numout,*) 'trc_adv_ini : choice/control of the tracer advection scheme' 
     175         WRITE(numout,*) '~~~~~~~~~~~' 
     176         WRITE(numout,*) '   Namelist namtrc_adv : chose a advection scheme for tracers' 
     177         WRITE(numout,*) '      centered scheme                           ln_trcadv_cen = ', ln_trcadv_cen 
     178         WRITE(numout,*) '            horizontal 2nd/4th order               nn_cen_h   = ', nn_fct_h 
     179         WRITE(numout,*) '            vertical   2nd/4th order               nn_cen_v   = ', nn_fct_v 
     180         WRITE(numout,*) '      Flux Corrected Transport scheme           ln_trcadv_fct = ', ln_trcadv_fct 
     181         WRITE(numout,*) '            horizontal 2nd/4th order               nn_fct_h   = ', nn_fct_h 
     182         WRITE(numout,*) '            vertical   2nd/4th order               nn_fct_v   = ', nn_fct_v 
     183         WRITE(numout,*) '            2nd order + vertical sub-timestepping  nn_fct_zts = ', nn_fct_zts 
     184         WRITE(numout,*) '      MUSCL scheme                              ln_trcadv_mus = ', ln_trcadv_mus 
     185         WRITE(numout,*) '            + upstream scheme near river mouths    ln_mus_ups = ', ln_mus_ups 
     186         WRITE(numout,*) '      UBS scheme                                ln_trcadv_ubs = ', ln_trcadv_ubs 
     187         WRITE(numout,*) '            vertical   2nd/4th order               nn_ubs_v   = ', nn_ubs_v 
     188         WRITE(numout,*) '      QUICKEST scheme                           ln_trcadv_qck = ', ln_trcadv_qck 
     189      ENDIF 
     190      ! 
     191 
     192      ioptio = 0                       !==  Parameter control  ==! 
     193      IF( ln_trcadv_cen )   ioptio = ioptio + 1 
     194      IF( ln_trcadv_fct )   ioptio = ioptio + 1 
     195      IF( ln_trcadv_mus )   ioptio = ioptio + 1 
     196      IF( ln_trcadv_ubs )   ioptio = ioptio + 1 
     197      IF( ln_trcadv_qck )   ioptio = ioptio + 1 
     198 
     199      ! 
     200      IF( ioptio == 0 ) THEN 
     201         nadv = np_NO_adv 
     202         CALL ctl_warn( 'trc_adv_init: You are running without tracer advection.' ) 
     203      ENDIF 
    180204      IF( ioptio /= 1 )   CALL ctl_stop( 'Choose ONE advection scheme in namelist namtrc_adv' ) 
    181  
    182       !                              ! Set nadv 
    183       IF( ln_trcadv_cen2   )   nadv =  1 
    184       IF( ln_trcadv_tvd    )   nadv =  2 
    185       IF( ln_trcadv_muscl  )   nadv =  3 
    186       IF( ln_trcadv_muscl2 )   nadv =  4 
    187       IF( ln_trcadv_ubs    )   nadv =  5 
    188       IF( ln_trcadv_qck    )   nadv =  6 
    189       IF( lk_esopa         )   nadv = -1 
    190  
     205      ! 
     206      IF( ln_trcadv_cen .AND. ( nn_cen_h /= 2 .AND. nn_cen_h /= 4 )   & 
     207                        .AND. ( nn_cen_v /= 2 .AND. nn_cen_v /= 4 )   ) THEN 
     208        CALL ctl_stop( 'trc_adv_init: CEN scheme, choose 2nd or 4th order' ) 
     209      ENDIF 
     210      IF( ln_trcadv_fct .AND. ( nn_fct_h /= 2 .AND. nn_fct_h /= 4 )   & 
     211                        .AND. ( nn_fct_v /= 2 .AND. nn_fct_v /= 4 )   ) THEN 
     212        CALL ctl_stop( 'trc_adv_init: FCT scheme, choose 2nd or 4th order' ) 
     213      ENDIF 
     214      IF( ln_trcadv_fct .AND. nn_fct_zts > 0 ) THEN 
     215         IF( nn_fct_h == 4 ) THEN 
     216            nn_fct_h = 2 
     217            CALL ctl_stop( 'trc_adv_init: force 2nd order FCT scheme, 4th order does not exist with sub-timestepping' ) 
     218         ENDIF 
     219         IF( .NOT.ln_linssh ) THEN 
     220            CALL ctl_stop( 'trc_adv_init: vertical sub-timestepping not allow in non-linear free surface' ) 
     221         ENDIF 
     222         IF( nn_fct_zts == 1 )   CALL ctl_warn( 'trc_adv_init: FCT with ONE sub-timestep = FCT without sub-timestep' ) 
     223      ENDIF 
     224      IF( ln_trcadv_ubs .AND. ( nn_ubs_v /= 2 .AND. nn_ubs_v /= 4 )   ) THEN 
     225        CALL ctl_stop( 'trc_adv_init: UBS scheme, choose 2nd or 4th order' ) 
     226      ENDIF 
     227      IF( ln_trcadv_ubs .AND. nn_ubs_v == 4 ) THEN 
     228         CALL ctl_warn( 'trc_adv_init: UBS scheme, only 2nd FCT scheme available on the vertical. It will be used' ) 
     229      ENDIF 
     230      IF( ln_isfcav ) THEN                                                       ! ice-shelf cavities 
     231         IF(  ln_trcadv_cen .AND. nn_cen_v /= 4    .OR.   &                            ! NO 4th order with ISF 
     232            & ln_trcadv_fct .AND. nn_fct_v /= 4   )   CALL ctl_stop( 'tra_adv_init: 4th order COMPACT scheme not allowed with ISF' ) 
     233      ENDIF 
     234      ! 
     235      !                                !==  used advection scheme  ==! 
     236      !                                      ! set nadv 
     237      IF( ln_trcadv_cen                      )   nadv = np_CEN 
     238      IF( ln_trcadv_fct                      )   nadv = np_FCT 
     239      IF( ln_trcadv_fct .AND. nn_fct_zts > 0 )   nadv = np_FCT_zts 
     240      IF( ln_trcadv_mus                      )   nadv = np_MUS 
     241      IF( ln_trcadv_ubs                      )   nadv = np_UBS 
     242      IF( ln_trcadv_qck                      )   nadv = np_QCK 
     243      ! 
    191244      IF(lwp) THEN                   ! Print the choice 
    192245         WRITE(numout,*) 
    193          IF( nadv ==  1 )   WRITE(numout,*) '         2nd order scheme is used' 
    194          IF( nadv ==  2 )   WRITE(numout,*) '         TVD       scheme is used' 
    195          IF( nadv ==  3 )   WRITE(numout,*) '         MUSCL     scheme is used' 
    196          IF( nadv ==  4 )   WRITE(numout,*) '         MUSCL2    scheme is used' 
    197          IF( nadv ==  5 )   WRITE(numout,*) '         UBS       scheme is used' 
    198          IF( nadv ==  6 )   WRITE(numout,*) '         QUICKEST  scheme is used' 
    199          IF( nadv == -1 )   WRITE(numout,*) '         esopa test: use all advection scheme' 
    200       ENDIF 
    201       ! 
    202    END SUBROUTINE trc_adv_ctl 
     246         IF( nadv == np_NO_adv  )   WRITE(numout,*) '         NO passive tracer advection' 
     247         IF( nadv == np_CEN     )   WRITE(numout,*) '         CEN      scheme is used. Horizontal order: ', nn_cen_h,   & 
     248            &                                                                        ' Vertical   order: ', nn_cen_v 
     249         IF( nadv == np_FCT     )   WRITE(numout,*) '         FCT      scheme is used. Horizontal order: ', nn_fct_h,   & 
     250            &                                                                       ' Vertical   order: ', nn_fct_v 
     251         IF( nadv == np_FCT_zts )   WRITE(numout,*) '         use 2nd order FCT with ', nn_fct_zts,'vertical sub-timestepping' 
     252         IF( nadv == np_MUS     )   WRITE(numout,*) '         MUSCL    scheme is used' 
     253         IF( nadv == np_UBS     )   WRITE(numout,*) '         UBS      scheme is used' 
     254         IF( nadv == np_QCK     )   WRITE(numout,*) '         QUICKEST scheme is used' 
     255      ENDIF 
     256      ! 
     257   END SUBROUTINE trc_adv_ini 
    203258    
    204259#else 

  • branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/TOP_SRC/TRP/trcbbl.F90

    r4513 r6225  
    2222   USE oce_trc             ! ocean dynamics and active tracers variables 
    2323   USE trc                 ! ocean passive tracers variables 
    24    USE trcnam_trp      ! passive tracers transport namelist variables 
    2524   USE trabbl              !  
    2625   USE prtctl_trc          ! Print control for debbuging 
    27    USE trdmod_oce 
     26   USE trd_oce 
    2827   USE trdtra 
    2928 
    3029   PUBLIC   trc_bbl       !  routine called by step.F90 
    3130 
    32  
    33    !! * Substitutions 
    34 #  include "top_substitute.h90" 
    3531   !!---------------------------------------------------------------------- 
    3632   !! NEMO/TOP 3.3 , NEMO Consortium (2010) 
    … …  
    9389        DO jn = 1, jptra 
    9490           ztrtrd(:,:,:,jn) = tra(:,:,:,jn) - ztrtrd(:,:,:,jn) 
    95            CALL trd_tra( kt, 'TRC', jn, jptra_trd_bbl, ztrtrd(:,:,:,jn) ) 
     91           CALL trd_tra( kt, 'TRC', jn, jptra_bbl, ztrtrd(:,:,:,jn) ) 
    9692        END DO 
    9793        CALL wrk_dealloc( jpi, jpj, jpk, jptra, ztrtrd ) ! temporary save of trends 

  • branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/TOP_SRC/TRP/trcdmp.F90

    r4359 r6225  
    1818   USE oce_trc         ! ocean dynamics and tracers variables 
    1919   USE trc             ! ocean passive tracers variables 
    20    USE trcnam_trp      ! passive tracers transport namelist variables 
    2120   USE trcdta 
    2221   USE tradmp 
    2322   USE prtctl_trc      ! Print control for debbuging 
    2423   USE trdtra 
    25    USE trdmod_oce 
     24   USE trd_oce 
     25   USE iom 
    2626 
    2727   IMPLICIT NONE 
    2828   PRIVATE 
    2929 
    30    PUBLIC trc_dmp            ! routine called by step.F90 
    31    PUBLIC trc_dmp_clo        ! routine called by step.F90 
    32    PUBLIC trc_dmp_alloc      ! routine called by nemogcm.F90 
     30   PUBLIC trc_dmp       
     31   PUBLIC trc_dmp_clo    
     32   PUBLIC trc_dmp_alloc   
     33   PUBLIC trc_dmp_ini     
     34 
     35   INTEGER , PUBLIC ::   nn_zdmp_tr    ! = 0/1/2 flag for damping in the mixed layer 
     36   CHARACTER(LEN=200) , PUBLIC :: cn_resto_tr    !File containing restoration coefficient 
    3337 
    3438   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   restotr   ! restoring coeff. on tracers (s-1) 
    … …  
    3943 
    4044   !! * Substitutions 
    41 #  include "top_substitute.h90" 
     45#  include "vectopt_loop_substitute.h90" 
    4246   !!---------------------------------------------------------------------- 
    4347   !! NEMO/TOP 3.3 , NEMO Consortium (2010) 
    44    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/TOP_SRC/TRP/trcdmp.F90,v 1.11 2006/09/01 14:03:49 opalod Exp $  
     48   !! $Id$  
    4549   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    4650   !!---------------------------------------------------------------------- 
    … …  
    7579      !! ** Action  : - update the tracer trends tra with the newtonian  
    7680      !!                damping trends. 
    77       !!              - save the trends ('key_trdmld_trc') 
    78       !!---------------------------------------------------------------------- 
    79       !! 
    80       INTEGER, INTENT( in ) ::   kt      ! ocean time-step index 
    81       !! 
    82       INTEGER  ::   ji, jj, jk, jn, jl       ! dummy loop indices 
    83       REAL(wp) ::   ztra                 ! temporary scalars 
    84       CHARACTER (len=22) :: charout 
     81      !!              - save the trends ('key_trdmxl_trc') 
     82      !!---------------------------------------------------------------------- 
     83      INTEGER, INTENT(in) ::   kt   ! ocean time-step index 
     84      ! 
     85      INTEGER ::   ji, jj, jk, jn, jl   ! dummy loop indices 
     86      CHARACTER (len=22) ::   charout 
    8587      REAL(wp), POINTER, DIMENSION(:,:,:) ::   ztrtrd 
    86       REAL(wp), POINTER, DIMENSION(:,:,:) ::  ztrcdta   ! 3D  workspace 
     88      REAL(wp), POINTER, DIMENSION(:,:,:) ::   ztrcdta   ! 3D  workspace 
    8789      !!---------------------------------------------------------------------- 
    8890      ! 
    8991      IF( nn_timing == 1 )  CALL timing_start('trc_dmp') 
    9092      ! 
    91       ! 0. Initialization (first time-step only) 
    92       !    -------------- 
    93       IF( kt == nittrc000 ) CALL trc_dmp_init 
    94  
    9593      IF( l_trdtrc )   CALL wrk_alloc( jpi, jpj, jpk, ztrtrd )   ! temporary save of trends 
    9694      ! 
    … …  
    104102            ! 
    105103            IF( ln_trc_ini(jn) ) THEN      ! update passive tracers arrays with input data read from file 
    106                 
     104               ! 
    107105               jl = n_trc_index(jn)  
    108106               CALL trc_dta( kt, sf_trcdta(jl),rf_trfac(jl) )   ! read tracer data at nit000 
    109107               ztrcdta(:,:,:) = sf_trcdta(jl)%fnow(:,:,:) 
    110  
     108               ! 
    111109               SELECT CASE ( nn_zdmp_tr ) 
    112110               ! 
    … …  
    115113                     DO jj = 2, jpjm1 
    116114                        DO ji = fs_2, fs_jpim1   ! vector opt. 
    117                            ztra = restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk) - trb(ji,jj,jk,jn) ) 
    118                            tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + ztra 
     115                           tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk) - trb(ji,jj,jk,jn) ) 
    119116                        END DO 
    120117                     END DO 
    121118                  END DO 
    122                ! 
     119                  ! 
    123120               CASE ( 1 )                !==  no damping in the turbocline (avt > 5 cm2/s)  ==! 
    124121                  DO jk = 1, jpkm1 
    125122                     DO jj = 2, jpjm1 
    126123                        DO ji = fs_2, fs_jpim1   ! vector opt. 
    127                            IF( avt(ji,jj,jk) <= 5.e-4 )  THEN  
    128                               ztra = restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk) - trb(ji,jj,jk,jn) ) 
    129                               tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + ztra 
     124                           IF( avt(ji,jj,jk) <= 5.e-4_wp )  THEN  
     125                              tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk) - trb(ji,jj,jk,jn) ) 
    130126                           ENDIF 
    131127                        END DO 
    132128                     END DO 
    133129                  END DO 
    134                ! 
     130                  ! 
    135131               CASE ( 2 )               !==  no damping in the mixed layer   ==!  
    136132                  DO jk = 1, jpkm1 
    137133                     DO jj = 2, jpjm1 
    138134                        DO ji = fs_2, fs_jpim1   ! vector opt. 
    139                            IF( fsdept(ji,jj,jk) >= hmlp (ji,jj) ) THEN 
    140                               ztra = restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk) - trb(ji,jj,jk,jn) ) 
    141                               tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + ztra 
     135                           IF( gdept_n(ji,jj,jk) >= hmlp (ji,jj) ) THEN 
     136                              tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk) - trb(ji,jj,jk,jn) ) 
    142137                           END IF 
    143138                        END DO 
    144139                     END DO 
    145140                  END DO 
    146                !   
     141                  !   
    147142               END SELECT 
    148143               !  
    … …  
    151146            IF( l_trdtrc ) THEN 
    152147               ztrtrd(:,:,:) = tra(:,:,:,jn) -  ztrtrd(:,:,:) 
    153                CALL trd_tra( kt, 'TRC', jn, jptra_trd_dmp, ztrtrd ) 
     148               CALL trd_tra( kt, 'TRC', jn, jptra_dmp, ztrtrd ) 
    154149            END IF 
    155150            !                                                       ! =========== 
    … …  
    161156      IF( l_trdtrc )  CALL wrk_dealloc( jpi, jpj, jpk, ztrtrd ) 
    162157      !                                          ! print mean trends (used for debugging) 
    163       IF( ln_ctl )   THEN 
    164          WRITE(charout, FMT="('dmp ')") ;  CALL prt_ctl_trc_info(charout) 
    165                                            CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 
     158      IF( ln_ctl ) THEN 
     159         WRITE(charout, FMT="('dmp ')") 
     160         CALL prt_ctl_trc_info(charout) 
     161         CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 
    166162      ENDIF 
    167163      ! 
    … …  
    169165      ! 
    170166   END SUBROUTINE trc_dmp 
     167 
     168 
     169   SUBROUTINE trc_dmp_ini 
     170      !!---------------------------------------------------------------------- 
     171      !!                  ***  ROUTINE trc_dmp_ini  *** 
     172      !!  
     173      !! ** Purpose :   Initialization for the newtonian damping  
     174      !! 
     175      !! ** Method  :   read the nammbf namelist and check the parameters 
     176      !!              called by trc_dmp at the first timestep (nittrc000) 
     177      !!---------------------------------------------------------------------- 
     178      INTEGER ::   ios, imask  ! local integers 
     179      !! 
     180      NAMELIST/namtrc_dmp/ nn_zdmp_tr , cn_resto_tr 
     181      !!---------------------------------------------------------------------- 
     182      ! 
     183      IF( nn_timing == 1 )  CALL timing_start('trc_dmp_init') 
     184      ! 
     185      REWIND( numnat_ref )              ! Namelist namtrc_dmp in reference namelist : Passive tracers newtonian damping 
     186      READ  ( numnat_ref, namtrc_dmp, IOSTAT = ios, ERR = 909) 
     187909   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrc_dmp in reference namelist', lwp ) 
     188 
     189      REWIND( numnat_cfg )              ! Namelist namtrc_dmp in configuration namelist : Passive tracers newtonian damping 
     190      READ  ( numnat_cfg, namtrc_dmp, IOSTAT = ios, ERR = 910) 
     191910   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrc_dmp in configuration namelist', lwp ) 
     192      IF(lwm) WRITE ( numont, namtrc_dmp ) 
     193 
     194      IF(lwp) THEN                       ! Namelist print 
     195         WRITE(numout,*) 
     196         WRITE(numout,*) 'trc_dmp : Passive tracers newtonian damping' 
     197         WRITE(numout,*) '~~~~~~~' 
     198         WRITE(numout,*) '   Namelist namtrc_dmp : set damping parameter' 
     199         WRITE(numout,*) '      mixed layer damping option     nn_zdmp_tr = ', nn_zdmp_tr, '(zoom: forced to 0)' 
     200         WRITE(numout,*) '      Restoration coeff file    cn_resto_tr = ', cn_resto_tr 
     201      ENDIF 
     202      ! 
     203      IF( lzoom .AND. .NOT.lk_c1d )   nn_zdmp_tr = 0           ! restoring to climatology at closed north or south boundaries 
     204      SELECT CASE ( nn_zdmp_tr ) 
     205      CASE ( 0 )   ;   IF(lwp) WRITE(numout,*) '   tracer damping throughout the water column' 
     206      CASE ( 1 )   ;   IF(lwp) WRITE(numout,*) '   no tracer damping in the turbocline (avt > 5 cm2/s)' 
     207      CASE ( 2 )   ;   IF(lwp) WRITE(numout,*) '   no tracer damping in the mixed layer' 
     208      CASE DEFAULT 
     209         WRITE(ctmp1,*) 'bad flag value for nn_zdmp_tr = ', nn_zdmp_tr 
     210         CALL ctl_stop(ctmp1) 
     211      END SELECT 
     212 
     213      IF( .NOT.lk_c1d ) THEN 
     214         IF( .NOT. ln_tradmp )   & 
     215            &   CALL ctl_stop( 'passive trace damping need ln_tradmp to compute damping coef.' ) 
     216         ! 
     217         !                          ! Read damping coefficients from file 
     218         !Read in mask from file 
     219         CALL iom_open ( cn_resto_tr, imask) 
     220         CALL iom_get  ( imask, jpdom_autoglo, 'resto', restotr) 
     221         CALL iom_close( imask ) 
     222         ! 
     223      ENDIF 
     224      IF( nn_timing == 1 )  CALL timing_stop('trc_dmp_init') 
     225      ! 
     226   END SUBROUTINE trc_dmp_ini 
     227 
    171228 
    172229   SUBROUTINE trc_dmp_clo( kt ) 
    … …  
    182239      !!                nctsi2(), nctsj2() : north-east Closed sea limits (i,j) 
    183240      !!---------------------------------------------------------------------- 
    184       INTEGER, INTENT( in ) ::   kt      ! ocean time-step index 
    185       ! 
    186       INTEGER :: ji, jj, jk, jn, jl, jc                     ! dummy loop indicesa 
    187       REAL(wp), POINTER, DIMENSION(:,:,:) ::  ztrcdta       ! 3D  workspace 
    188  
    189       !!---------------------------------------------------------------------- 
    190  
     241      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index 
     242      ! 
     243      INTEGER ::   ji , jj, jk, jn, jl, jc   ! dummy loop indicesa 
     244      INTEGER ::   isrow                     ! local index 
     245      REAL(wp), POINTER, DIMENSION(:,:,:) ::   ztrcdta   ! 3D  workspace 
     246      !!---------------------------------------------------------------------- 
     247      ! 
    191248      IF( kt == nit000 ) THEN 
    192249         ! initial values 
    … …  
    200257            ! 
    201258            SELECT CASE ( jp_cfg ) 
     259            !                                           ! ======================= 
     260            CASE ( 1 )                                  ! eORCA_R1 configuration 
     261            !                                           ! ======================= 
     262            isrow = 332 - jpjglo 
     263            ! 
     264                                                        ! Caspian Sea 
     265            nctsi1(1)   = 332  ; nctsj1(1)   = 243 - isrow 
     266            nctsi2(1)   = 344  ; nctsj2(1)   = 275 - isrow 
     267            !                                         
    202268            !                                           ! ======================= 
    203269            CASE ( 2 )                                  !  ORCA_R2 configuration 
    … …  
    291357   END SUBROUTINE trc_dmp_clo 
    292358 
    293  
    294    SUBROUTINE trc_dmp_init 
    295       !!---------------------------------------------------------------------- 
    296       !!                  ***  ROUTINE trc_dmp_init  *** 
    297       !!  
    298       !! ** Purpose :   Initialization for the newtonian damping  
    299       !! 
    300       !! ** Method  :   read the nammbf namelist and check the parameters 
    301       !!              called by trc_dmp at the first timestep (nittrc000) 
    302       !!---------------------------------------------------------------------- 
    303       ! 
    304       IF( nn_timing == 1 )  CALL timing_start('trc_dmp_init') 
    305       ! 
    306       SELECT CASE ( nn_hdmp_tr ) 
    307       CASE (  -1  )   ;   IF(lwp) WRITE(numout,*) '   tracer damping in the Med & Red seas only' 
    308       CASE ( 1:90 )   ;   IF(lwp) WRITE(numout,*) '   tracer damping poleward of', nn_hdmp_tr, ' degrees' 
    309       CASE DEFAULT 
    310          WRITE(ctmp1,*) '          bad flag value for nn_hdmp_tr = ', nn_hdmp_tr 
    311          CALL ctl_stop(ctmp1) 
    312       END SELECT 
    313  
    314       IF( lzoom )   nn_zdmp_tr = 0           ! restoring to climatology at closed north or south boundaries 
    315       SELECT CASE ( nn_zdmp_tr ) 
    316       CASE ( 0 )   ;   IF(lwp) WRITE(numout,*) '   tracer damping throughout the water column' 
    317       CASE ( 1 )   ;   IF(lwp) WRITE(numout,*) '   no tracer damping in the turbocline (avt > 5 cm2/s)' 
    318       CASE ( 2 )   ;   IF(lwp) WRITE(numout,*) '   no tracer damping in the mixed layer' 
    319       CASE DEFAULT 
    320          WRITE(ctmp1,*) 'bad flag value for nn_zdmp_tr = ', nn_zdmp_tr 
    321          CALL ctl_stop(ctmp1) 
    322       END SELECT 
    323  
    324       IF( .NOT. ln_tradmp )   & 
    325          &   CALL ctl_stop( 'passive trace damping need key_tradmp to compute damping coef.' ) 
    326       ! 
    327       !                          ! Damping coefficients initialization 
    328       IF( lzoom ) THEN   ;   CALL dtacof_zoom( restotr ) 
    329       ELSE               ;   CALL dtacof( nn_hdmp_tr, rn_surf_tr, rn_bot_tr, rn_dep_tr,  & 
    330                              &            nn_file_tr, 'TRC'     , restotr                ) 
    331       ENDIF 
    332       ! 
    333       IF( nn_timing == 1 )  CALL timing_stop('trc_dmp_init') 
    334       ! 
    335    END SUBROUTINE trc_dmp_init 
    336  
    337359#else 
    338360   !!---------------------------------------------------------------------- 
    … …  
    346368#endif 
    347369 
    348  
    349370   !!====================================================================== 
    350371END MODULE trcdmp 

  • branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/TOP_SRC/TRP/trcldf.F90

    r3294 r6225  
    44   !! Ocean Passive tracers : lateral diffusive trends 
    55   !!===================================================================== 
    6    !! History :  9.0  ! 2005-11 (G. Madec)  Original code 
    7    !!       NEMO 3.0  ! 2008-01  (C. Ethe, G. Madec)  merge TRC-TRA 
     6   !! History :  1.0  ! 2005-11  (G. Madec)  Original code 
     7   !!            3.0  ! 2008-01  (C. Ethe, G. Madec)  merge TRC-TRA 
     8   !!            3.7  ! 2014-03  (G. Madec)  LDF simplification 
    89   !!---------------------------------------------------------------------- 
    910#if defined key_top 
    … …  
    1112   !!   'key_top'                                                TOP models 
    1213   !!---------------------------------------------------------------------- 
    13    !!---------------------------------------------------------------------- 
    14    !!   trc_ldf     : update the tracer trend with the lateral diffusion 
    15    !!       ldf_ctl : initialization, namelist read, and parameters control 
    16    !!---------------------------------------------------------------------- 
    17    USE oce_trc         ! ocean dynamics and active tracers 
    18    USE trc             ! ocean passive tracers variables 
    19    USE trcnam_trp      ! passive tracers transport namelist variables 
    20    USE ldftra_oce      ! lateral diffusion coefficient on tracers 
    21    USE ldfslp          ! ??? 
    22    USE traldf_bilapg   ! lateral mixing            (tra_ldf_bilapg routine) 
    23    USE traldf_bilap    ! lateral mixing            (tra_ldf_bilap routine) 
    24    USE traldf_iso      ! lateral mixing            (tra_ldf_iso routine) 
    25    USE traldf_iso_grif ! lateral mixing          (tra_ldf_iso_grif routine) 
    26    USE traldf_lap      ! lateral mixing            (tra_ldf_lap routine) 
    27    USE trdmod_oce 
    28    USE trdtra 
    29    USE prtctl_trc      ! Print control 
     14   !!   trc_ldf       : update the tracer trend with the lateral diffusion 
     15   !!   trc_ldf_ini   : initialization, namelist read, and parameters control 
     16   !!---------------------------------------------------------------------- 
     17   USE trc            ! ocean passive tracers variables 
     18   USE oce_trc        ! ocean dynamics and active tracers 
     19   USE ldfslp         ! lateral diffusion: iso-neutral slope 
     20   USE traldf_lap_blp ! lateral diffusion: lap/bilaplacian iso-level      operator  (tra_ldf_lap/_blp   routine) 
     21   USE traldf_iso     ! lateral diffusion: laplacian iso-neutral standard operator  (tra_ldf_iso        routine) 
     22   USE traldf_triad   ! lateral diffusion: laplacian iso-neutral triad    operator  (tra_ldf_     triad routine) 
     23   USE trd_oce        ! trends: ocean variables 
     24   USE trdtra         ! trends manager: tracers 
     25   ! 
     26   USE prtctl_trc     ! Print control 
    3027 
    3128   IMPLICIT NONE 
    3229   PRIVATE 
    3330 
    34    PUBLIC   trc_ldf    ! called by step.F90 
    35    !                                                 !!: ** lateral mixing namelist (nam_trcldf) ** 
    36    REAL(wp) ::  rldf_rat    ! ratio between active and passive tracers diffusive coefficient 
     31   PUBLIC   trc_ldf     
     32   PUBLIC   trc_ldf_ini    
     33   ! 
     34   LOGICAL , PUBLIC ::   ln_trcldf_lap       !:   laplacian operator 
     35   LOGICAL , PUBLIC ::   ln_trcldf_blp       !: bilaplacian operator 
     36   LOGICAL , PUBLIC ::   ln_trcldf_lev       !: iso-level   direction 
     37   LOGICAL , PUBLIC ::   ln_trcldf_hor       !: horizontal  direction (rotation to geopotential) 
     38   LOGICAL , PUBLIC ::   ln_trcldf_iso       !: iso-neutral direction (standard) 
     39   LOGICAL , PUBLIC ::   ln_trcldf_triad     !: iso-neutral direction (triad) 
     40   REAL(wp), PUBLIC ::   rn_ahtrc_0          !:   laplacian diffusivity coefficient for passive tracer [m2/s] 
     41   REAL(wp), PUBLIC ::   rn_bhtrc_0          !: bilaplacian      -          --     -       -   [m4/s] 
     42   ! 
     43   !                      !!: ** lateral mixing namelist (nam_trcldf) ** 
     44   REAL(wp) ::  rldf       ! ratio between active and passive tracers diffusive coefficient 
     45    
    3746   INTEGER  ::  nldf = 0   ! type of lateral diffusion used defined from ln_trcldf_... namlist logicals) 
     47    
    3848   !! * Substitutions 
    39 #  include "domzgr_substitute.h90" 
    4049#  include "vectopt_loop_substitute.h90" 
    4150   !!---------------------------------------------------------------------- 
    42    !! NEMO/TOP 3.3 , NEMO Consortium (2010) 
     51   !! NEMO/TOP 3.7 , NEMO Consortium (2014) 
    4352   !! $Id$ 
    4453   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    4554   !!---------------------------------------------------------------------- 
    46  
    4755CONTAINS 
    4856 
    … …  
    5563      !!---------------------------------------------------------------------- 
    5664      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index 
    57       !! 
     65      ! 
    5866      INTEGER            :: jn 
    5967      CHARACTER (len=22) :: charout 
     68      REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zahu, zahv 
    6069      REAL(wp), POINTER, DIMENSION(:,:,:,:) ::   ztrtrd 
    6170      !!---------------------------------------------------------------------- 
    … …  
    6372      IF( nn_timing == 1 )   CALL timing_start('trc_ldf') 
    6473      ! 
    65       IF( kt == nittrc000 )   CALL ldf_ctl          ! initialisation & control of options 
    66  
    67       rldf = rldf_rat 
    68  
    6974      IF( l_trdtrc )  THEN 
    70          CALL wrk_alloc( jpi, jpj, jpk, jptra, ztrtrd ) 
     75         CALL wrk_alloc( jpi,jpj,jpk,jptra,  ztrtrd ) 
    7176         ztrtrd(:,:,:,:)  = tra(:,:,:,:) 
    7277      ENDIF 
    73  
    74       SELECT CASE ( nldf )                       ! compute lateral mixing trend and add it to the general trend 
    75       CASE ( 0 )   ;   CALL tra_ldf_lap   ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra            )  ! iso-level laplacian 
    76       CASE ( 1 )                                                                                            ! rotated laplacian 
    77                        IF( ln_traldf_grif ) THEN 
    78                           CALL tra_ldf_iso_grif( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 ) 
    79                        ELSE 
    80                           CALL tra_ldf_iso     ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 ) 
    81                        ENDIF 
    82       CASE ( 2 )   ;   CALL tra_ldf_bilap ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra            )  ! iso-level bilaplacian 
    83       CASE ( 3 )   ;   CALL tra_ldf_bilapg( kt, nittrc000, 'TRC',             trb, tra, jptra            )  ! s-coord. horizontal bilaplacian 
    84          ! 
    85       CASE ( -1 )                                     ! esopa: test all possibility with control print 
    86          CALL tra_ldf_lap   ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra            ) 
    87          WRITE(charout, FMT="('ldf0 ')") ;  CALL prt_ctl_trc_info(charout) 
    88                                             CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 
    89          IF( ln_traldf_grif ) THEN 
    90             CALL tra_ldf_iso_grif( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 ) 
    91          ELSE 
    92             CALL tra_ldf_iso     ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 ) 
    93          ENDIF 
    94          WRITE(charout, FMT="('ldf1 ')") ;  CALL prt_ctl_trc_info(charout) 
    95                                             CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 
    96          CALL tra_ldf_bilap ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra            ) 
    97          WRITE(charout, FMT="('ldf2 ')") ;  CALL prt_ctl_trc_info(charout) 
    98                                             CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 
    99          CALL tra_ldf_bilapg( kt, nittrc000, 'TRC',             trb, tra, jptra            ) 
    100          WRITE(charout, FMT="('ldf3 ')") ;  CALL prt_ctl_trc_info(charout) 
    101                                             CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 
     78      ! 
     79      !                                        !* set the lateral diffusivity coef. for passive tracer       
     80      CALL wrk_alloc( jpi,jpj,jpk,   zahu, zahv ) 
     81      zahu(:,:,:) = rldf * ahtu(:,:,:) 
     82      zahv(:,:,:) = rldf * ahtv(:,:,:) 
     83 
     84      SELECT CASE ( nldf )                     !* compute lateral mixing trend and add it to the general trend 
     85      ! 
     86      CASE ( np_lap   )                               ! iso-level laplacian 
     87         CALL tra_ldf_lap  ( kt, nittrc000,'TRC', zahu, zahv, gtru, gtrv, gtrui, gtrvi, trb,      tra, jptra,  1   ) 
     88         ! 
     89      CASE ( np_lap_i )                               ! laplacian : standard iso-neutral operator (Madec) 
     90         CALL tra_ldf_iso  ( kt, nittrc000,'TRC', zahu, zahv, gtru, gtrv, gtrui, gtrvi, trb, trb, tra, jptra,  1   ) 
     91         ! 
     92      CASE ( np_lap_it )                              ! laplacian : triad iso-neutral operator (griffies) 
     93         CALL tra_ldf_triad( kt, nittrc000,'TRC', zahu, zahv, gtru, gtrv, gtrui, gtrvi, trb, trb, tra, jptra,  1   ) 
     94         ! 
     95      CASE ( np_blp , np_blp_i , np_blp_it )          ! bilaplacian: all operator (iso-level, -neutral) 
     96         CALL tra_ldf_blp  ( kt, nittrc000,'TRC', zahu, zahv, gtru, gtrv, gtrui, gtrvi, trb     , tra, jptra, nldf ) 
     97         ! 
    10298      END SELECT 
    10399      ! 
    104       IF( l_trdtrc )   THEN                      ! save the horizontal diffusive trends for further diagnostics 
     100      IF( l_trdtrc )   THEN                    ! send the trends for further diagnostics 
    105101        DO jn = 1, jptra 
    106102           ztrtrd(:,:,:,jn) = tra(:,:,:,jn) - ztrtrd(:,:,:,jn) 
    107            CALL trd_tra( kt, 'TRC', jn, jptra_trd_ldf, ztrtrd(:,:,:,jn) ) 
     103           CALL trd_tra( kt, 'TRC', jn, jptra_ldf, ztrtrd(:,:,:,jn) ) 
    108104        END DO 
    109105        CALL wrk_dealloc( jpi, jpj, jpk, jptra, ztrtrd ) 
    110106      ENDIF 
    111       !                                          ! print mean trends (used for debugging) 
    112       IF( ln_ctl )   THEN 
    113          WRITE(charout, FMT="('ldf ')") ;  CALL prt_ctl_trc_info(charout) 
    114                                            CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 
    115       ENDIF 
     107      !                 
     108      IF( ln_ctl ) THEN                        ! print mean trends (used for debugging) 
     109         WRITE(charout, FMT="('ldf ')") 
     110         CALL prt_ctl_trc_info(charout) 
     111         CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 
     112      ENDIF 
     113      ! 
     114      CALL wrk_dealloc( jpi,jpj,jpk,   zahu, zahv ) 
    116115      ! 
    117116      IF( nn_timing == 1 )   CALL timing_stop('trc_ldf') 
    … …  
    120119 
    121120 
    122    SUBROUTINE ldf_ctl 
     121   SUBROUTINE trc_ldf_ini 
    123122      !!---------------------------------------------------------------------- 
    124123      !!                  ***  ROUTINE ldf_ctl  *** 
    125124      !! 
    126       !! ** Purpose :   Choice of the operator for the lateral tracer diffusion 
     125      !! ** Purpose :   Define the operator for the lateral diffusion 
    127126      !! 
    128127      !! ** Method  :   set nldf from the namtra_ldf logicals 
    129       !!      nldf == -2   No lateral diffusion 
    130       !!      nldf == -1   ESOPA test: ALL operators are used 
    131128      !!      nldf ==  0   laplacian operator 
    132129      !!      nldf ==  1   Rotated laplacian operator 
    … …  
    134131      !!      nldf ==  3   Rotated bilaplacian 
    135132      !!---------------------------------------------------------------------- 
    136       INTEGER ::   ioptio, ierr         ! temporary integers 
    137       !!---------------------------------------------------------------------- 
    138  
    139       IF (ABS(rn_aht_0) < 2._wp*TINY(1.e0)) THEN 
    140          IF (ABS(rn_ahtrc_0) < 2._wp*TINY(1.e0)) THEN 
    141             rldf_rat = 1.0_wp 
    142          ELSE 
    143             CALL ctl_stop( 'STOP', 'ldf_ctl : cannot define rldf_rat, rn_aht_0==0, rn_ahtrc_0 /=0' ) 
    144          END IF 
    145       ELSE 
    146          rldf_rat = rn_ahtrc_0 / rn_aht_0 
    147       END IF 
    148       !  Define the lateral mixing oparator for tracers 
    149       ! =============================================== 
    150  
    151       !                               ! control the input 
     133      INTEGER ::   ioptio, ierr   ! temporary integers 
     134      INTEGER ::   ios            ! Local integer output status for namelist read 
     135      !! 
     136      NAMELIST/namtrc_ldf/ ln_trcldf_lap, ln_trcldf_blp,                                  & 
     137         &                 ln_trcldf_lev, ln_trcldf_hor, ln_trcldf_iso, ln_trcldf_triad,  & 
     138         &                 rn_ahtrc_0   , rn_bhtrc_0 
     139      !!---------------------------------------------------------------------- 
     140      ! 
     141      REWIND( numnat_ref )             !  namtrc_ldf in reference namelist  
     142      READ  ( numnat_ref, namtrc_ldf, IOSTAT = ios, ERR = 903) 
     143903   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namtrc_ldf in reference namelist', lwp ) 
     144      ! 
     145      REWIND( numnat_cfg )             !  namtrc_ldf in configuration namelist  
     146      READ  ( numnat_cfg, namtrc_ldf, IOSTAT = ios, ERR = 904 ) 
     147904   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namtrc_ldf in configuration namelist', lwp ) 
     148      IF(lwm) WRITE ( numont, namtrc_ldf ) 
     149      ! 
     150      IF(lwp) THEN                     ! Namelist print 
     151         WRITE(numout,*) 
     152         WRITE(numout,*) 'trc_ldf_ini : lateral tracer diffusive operator' 
     153         WRITE(numout,*) '~~~~~~~~~~~' 
     154         WRITE(numout,*) '   Namelist namtrc_ldf : set lateral mixing parameters (type, direction, coefficients)' 
     155         WRITE(numout,*) '      operator' 
     156         WRITE(numout,*) '           laplacian                 ln_trcldf_lap   = ', ln_trcldf_lap 
     157         WRITE(numout,*) '         bilaplacian                 ln_trcldf_blp   = ', ln_trcldf_blp 
     158         WRITE(numout,*) '      direction of action' 
     159         WRITE(numout,*) '         iso-level                   ln_trcldf_lev   = ', ln_trcldf_lev 
     160         WRITE(numout,*) '         horizontal (geopotential)   ln_trcldf_hor   = ', ln_trcldf_hor 
     161         WRITE(numout,*) '         iso-neutral (standard)      ln_trcldf_iso   = ', ln_trcldf_iso 
     162         WRITE(numout,*) '         iso-neutral (triad)         ln_trcldf_triad = ', ln_trcldf_triad 
     163         WRITE(numout,*) '      diffusivity coefficient' 
     164         WRITE(numout,*) '           laplacian                 rn_ahtrc_0      = ', rn_ahtrc_0 
     165         WRITE(numout,*) '         bilaplacian                 rn_bhtrc_0      = ', rn_bhtrc_0 
     166      ENDIF 
     167      !       
     168      !                                ! control the namelist parameters 
    152169      ioptio = 0 
    153       IF( ln_trcldf_lap   )   ioptio = ioptio + 1 
    154       IF( ln_trcldf_bilap )   ioptio = ioptio + 1 
    155       IF( ioptio >  1 )   CALL ctl_stop( '          use ONE or NONE of the 2 lap/bilap operator type on tracer' ) 
    156       IF( ioptio == 0 )   nldf = -2   ! No lateral diffusion 
     170      IF( ln_trcldf_lap )   ioptio = ioptio + 1 
     171      IF( ln_trcldf_blp )   ioptio = ioptio + 1 
     172      IF( ioptio >  1   )   CALL ctl_stop( 'trc_ldf_ctl: use ONE or NONE of the 2 lap/bilap operator type on tracer' ) 
     173      IF( ioptio == 0   )   nldf = np_no_ldf   ! No lateral diffusion 
     174       
     175      IF( ln_trcldf_lap .AND. ln_trcldf_blp )   CALL ctl_stop( 'trc_ldf_ctl: bilaplacian should be used on both TRC and TRA' ) 
     176      IF( ln_trcldf_blp .AND. ln_trcldf_lap )   CALL ctl_stop( 'trc_ldf_ctl:   laplacian should be used on both TRC and TRA' ) 
     177      ! 
    157178      ioptio = 0 
    158       IF( ln_trcldf_level )   ioptio = ioptio + 1 
    159       IF( ln_trcldf_hor   )   ioptio = ioptio + 1 
    160       IF( ln_trcldf_iso   )   ioptio = ioptio + 1 
    161       IF( ioptio /= 1 )   CALL ctl_stop( '          use only ONE direction (level/hor/iso)' ) 
    162  
     179      IF( ln_trcldf_lev )   ioptio = ioptio + 1 
     180      IF( ln_trcldf_hor )   ioptio = ioptio + 1 
     181      IF( ln_trcldf_iso )   ioptio = ioptio + 1 
     182      IF( ioptio /= 1   )   CALL ctl_stop( 'trc_ldf_ctl: use only ONE direction (level/hor/iso)' ) 
     183      ! 
    163184      ! defined the type of lateral diffusion from ln_trcldf_... logicals 
    164185      ! CAUTION : nldf = 1 is used in trazdf_imp, change it carefully 
    165186      ierr = 0 
    166       IF( ln_trcldf_lap ) THEN       ! laplacian operator 
     187      IF( ln_trcldf_lap ) THEN      !==  laplacian operator  ==! 
    167188         IF ( ln_zco ) THEN                ! z-coordinate 
    168             IF ( ln_trcldf_level )   nldf = 0      ! iso-level  (no rotation) 
    169             IF ( ln_trcldf_hor   )   nldf = 0      ! horizontal (no rotation) 
    170             IF ( ln_trcldf_iso   )   nldf = 1      ! isoneutral (   rotation) 
    171          ENDIF 
    172          IF ( ln_zps ) THEN             ! z-coordinate 
    173             IF ( ln_trcldf_level )   ierr = 1      ! iso-level not allowed 
    174             IF ( ln_trcldf_hor   )   nldf = 0      ! horizontal (no rotation) 
    175             IF ( ln_trcldf_iso   )   nldf = 1      ! isoneutral (   rotation) 
    176          ENDIF 
    177          IF ( ln_sco ) THEN             ! z-coordinate 
    178             IF ( ln_trcldf_level )   nldf = 0      ! iso-level  (no rotation) 
    179             IF ( ln_trcldf_hor   )   nldf = 1      ! horizontal (   rotation) 
    180             IF ( ln_trcldf_iso   )   nldf = 1      ! isoneutral (   rotation) 
    181          ENDIF 
    182       ENDIF 
    183  
    184       IF( ln_trcldf_bilap ) THEN      ! bilaplacian operator 
     189            IF ( ln_trcldf_lev   )   nldf = np_lap     ! iso-level = horizontal (no rotation) 
     190            IF ( ln_trcldf_hor   )   nldf = np_lap     ! iso-level = horizontal (no rotation) 
     191            IF ( ln_trcldf_iso   )   nldf = np_lap_i   ! iso-neutral: standard  (   rotation) 
     192            IF ( ln_trcldf_triad )   nldf = np_lap_it  ! iso-neutral: triad     (   rotation) 
     193         ENDIF 
     194         IF ( ln_zps ) THEN             ! z-coordinate with partial step 
     195            IF ( ln_trcldf_lev   )   ierr = 1         ! iso-level not allowed  
     196            IF ( ln_trcldf_hor   )   nldf = np_lap     ! horizontal (no rotation) 
     197            IF ( ln_trcldf_iso   )   nldf = np_lap_i   ! iso-neutral: standard (rotation) 
     198            IF ( ln_trcldf_triad )   nldf = np_lap_it  ! iso-neutral: triad    (rotation) 
     199         ENDIF 
     200         IF ( ln_sco ) THEN             ! s-coordinate 
     201            IF ( ln_trcldf_lev   )   nldf = np_lap     ! iso-level  (no rotation) 
     202            IF ( ln_trcldf_hor   )   nldf = np_lap_it  ! horizontal (   rotation)       !!gm   a checker.... 
     203            IF ( ln_trcldf_iso   )   nldf = np_lap_i   ! iso-neutral: standard (rotation) 
     204            IF ( ln_trcldf_triad )   nldf = np_lap_it  ! iso-neutral: triad    (rotation) 
     205         ENDIF 
     206         !                                ! diffusivity ratio: passive / active tracers  
     207         IF( ABS(rn_aht_0) < 2._wp*TINY(1._wp) ) THEN 
     208            IF( ABS(rn_ahtrc_0) < 2._wp*TINY(1._wp) ) THEN 
     209               rldf = 1.0_wp 
     210            ELSE 
     211               CALL ctl_stop( 'trc_ldf_ctl : cannot define rldf, rn_aht_0==0, rn_ahtrc_0 /=0' ) 
     212            ENDIF 
     213         ELSE 
     214            rldf = rn_ahtrc_0 / rn_aht_0 
     215         ENDIF 
     216      ENDIF 
     217      ! 
     218      IF( ln_trcldf_blp ) THEN      !==  bilaplacian operator  ==! 
    185219         IF ( ln_zco ) THEN                ! z-coordinate 
    186             IF ( ln_trcldf_level )   nldf = 2      ! iso-level  (no rotation) 
    187             IF ( ln_trcldf_hor   )   nldf = 2      ! horizontal (no rotation) 
    188             IF ( ln_trcldf_iso   )   ierr = 2      ! isoneutral (   rotation) 
    189          ENDIF 
    190          IF ( ln_zps ) THEN             ! z-coordinate 
    191             IF ( ln_trcldf_level )   ierr = 1      ! iso-level not allowed 
    192             IF ( ln_trcldf_hor   )   nldf = 2      ! horizontal (no rotation) 
    193             IF ( ln_trcldf_iso   )   ierr = 2      ! isoneutral (   rotation) 
    194          ENDIF 
    195          IF ( ln_sco ) THEN             ! z-coordinate 
    196             IF ( ln_trcldf_level )   nldf = 2      ! iso-level  (no rotation) 
    197             IF ( ln_trcldf_hor   )   nldf = 3      ! horizontal (   rotation) 
    198             IF ( ln_trcldf_iso   )   ierr = 2      ! isoneutral (   rotation) 
    199          ENDIF 
    200       ENDIF 
    201  
    202       IF( ierr == 1 )   CALL ctl_stop( ' iso-level in z-coordinate - partial step, not allowed' ) 
    203       IF( ierr == 2 )   CALL ctl_stop( ' isoneutral bilaplacian operator does not exist' ) 
    204       IF( lk_traldf_eiv .AND. .NOT.ln_trcldf_iso )   & 
    205            CALL ctl_stop( '          eddy induced velocity on tracers',   & 
    206            &              ' the eddy induced velocity on tracers requires isopycnal laplacian diffusion' ) 
    207       IF( nldf == 1 .OR. nldf == 3 ) THEN      ! rotation 
    208          IF( .NOT.lk_ldfslp )   CALL ctl_stop( '          the rotation of the diffusive tensor require key_ldfslp' ) 
    209 #if defined key_offline 
    210          l_traldf_rot = .TRUE.                 ! needed for trazdf_imp 
    211 #endif 
    212       ENDIF 
    213  
    214       IF( lk_esopa ) THEN 
    215          IF(lwp) WRITE(numout,*) '          esopa control: use all lateral physics options' 
    216          nldf = -1 
    217       ENDIF 
    218  
    219       IF( .NOT. ln_trcldf_diff ) THEN 
    220          IF(lwp) WRITE(numout,*) '          No lateral diffusion on passive tracers' 
    221          nldf = -2 
    222       ENDIF 
    223  
     220            IF ( ln_trcldf_lev   )   nldf = np_blp     ! iso-level = horizontal (no rotation) 
     221            IF ( ln_trcldf_hor   )   nldf = np_blp     ! iso-level = horizontal (no rotation) 
     222            IF ( ln_trcldf_iso   )   nldf = np_blp_i   ! iso-neutral: standard (rotation) 
     223            IF ( ln_trcldf_triad )   nldf = np_blp_it  ! iso-neutral: triad    (rotation) 
     224         ENDIF 
     225         IF ( ln_zps ) THEN             ! z-coordinate with partial step 
     226            IF ( ln_trcldf_lev   )   ierr = 1         ! iso-level not allowed  
     227            IF ( ln_trcldf_hor   )   nldf = np_blp     ! horizontal (no rotation) 
     228            IF ( ln_trcldf_iso   )   nldf = np_blp_i   ! iso-neutral: standard (rotation) 
     229            IF ( ln_trcldf_triad )   nldf = np_blp_it  ! iso-neutral: triad    (rotation) 
     230         ENDIF 
     231         IF ( ln_sco ) THEN             ! s-coordinate 
     232            IF ( ln_trcldf_lev   )   nldf = np_blp     ! iso-level  (no rotation) 
     233            IF ( ln_trcldf_hor   )   nldf = np_blp_it  ! horizontal (   rotation)       !!gm   a checker.... 
     234            IF ( ln_trcldf_iso   )   nldf = np_blp_i   ! iso-neutral: standard (rotation) 
     235            IF ( ln_trcldf_triad )   nldf = np_blp_it  ! iso-neutral: triad    (rotation) 
     236         ENDIF 
     237         !                                ! diffusivity ratio: passive / active tracers  
     238         IF( ABS(rn_bht_0) < 2._wp*TINY(1._wp) ) THEN 
     239            IF( ABS(rn_bhtrc_0) < 2._wp*TINY(1._wp) ) THEN 
     240               rldf = 1.0_wp 
     241            ELSE 
     242               CALL ctl_stop( 'trc_ldf_ctl : cannot define rldf, rn_aht_0==0, rn_ahtrc_0 /=0' ) 
     243            ENDIF 
     244         ELSE 
     245            rldf = SQRT(  ABS( rn_bhtrc_0 / rn_bht_0 )  ) 
     246         ENDIF 
     247      ENDIF 
     248      ! 
     249      IF( ierr == 1 )   CALL ctl_stop( 'trc_ldf_ctl: iso-level in z-partial step, not allowed' ) 
     250      IF( ln_ldfeiv .AND. .NOT.ln_trcldf_iso )   CALL ctl_stop( 'trc_ldf_ctl: eiv requires isopycnal laplacian diffusion' ) 
     251      IF( nldf == 1 .OR. nldf == 3 )   l_ldfslp = .TRUE.    ! slope of neutral surfaces required  
     252      ! 
    224253      IF(lwp) THEN 
    225254         WRITE(numout,*) 
    226          IF( nldf == -2 )   WRITE(numout,*) '          NO lateral diffusion' 
    227          IF( nldf == -1 )   WRITE(numout,*) '          ESOPA test All scheme used' 
    228          IF( nldf ==  0 )   WRITE(numout,*) '          laplacian operator' 
    229          IF( nldf ==  1 )   WRITE(numout,*) '          Rotated laplacian operator' 
    230          IF( nldf ==  2 )   WRITE(numout,*) '          bilaplacian operator' 
    231          IF( nldf ==  3 )   WRITE(numout,*) '          Rotated bilaplacian' 
    232       ENDIF 
    233  
    234       IF( ln_trcldf_bilap ) THEN 
    235          IF(lwp) WRITE(numout,*) '          biharmonic tracer diffusion' 
    236          IF( rn_ahtrc_0 > 0 .AND. .NOT. lk_esopa )   CALL ctl_stop( 'The horizontal diffusivity coef. rn_ahtrc_0 must be negative' ) 
    237       ELSE 
    238          IF(lwp) WRITE(numout,*) '          harmonic tracer diffusion (default)' 
    239          IF( rn_ahtrc_0 < 0 .AND. .NOT. lk_esopa )   CALL ctl_stop('The horizontal diffusivity coef. rn_ahtrc_0 must be positive' ) 
    240       ENDIF 
    241  
    242       ! ratio between active and passive tracers diffusive coef. 
    243       IF (ABS(rn_aht_0) < 2._wp*TINY(1.e0)) THEN 
    244          IF (ABS(rn_ahtrc_0) < 2._wp*TINY(1.e0)) THEN 
    245             rldf_rat = 1.0_wp 
    246          ELSE 
    247             CALL ctl_stop( 'STOP', 'ldf_ctl : cannot define rldf_rat, rn_aht_0==0, rn_ahtrc_0 /=0' ) 
    248          END IF 
    249       ELSE 
    250          rldf_rat = rn_ahtrc_0 / rn_aht_0 
    251       END IF 
    252       IF( rldf_rat < 0 ) THEN 
    253          IF( .NOT.lk_offline ) THEN  
    254             CALL ctl_stop( 'Choose the same type of diffusive scheme both for active & passive tracers' ) 
    255          ELSE 
    256             CALL ctl_stop( 'Change the sign of rn_aht_0 in namelist to -/+1' ) 
    257          ENDIF  
    258       ENDIF 
    259       ! 
    260    END SUBROUTINE ldf_ctl 
     255         SELECT CASE( nldf ) 
     256         CASE( np_no_ldf )   ;   WRITE(numout,*) '          NO lateral diffusion' 
     257         CASE( np_lap    )   ;   WRITE(numout,*) '          laplacian iso-level operator' 
     258         CASE( np_lap_i  )   ;   WRITE(numout,*) '          Rotated laplacian operator (standard)' 
     259         CASE( np_lap_it )   ;   WRITE(numout,*) '          Rotated laplacian operator (triad)' 
     260         CASE( np_blp    )   ;   WRITE(numout,*) '          bilaplacian iso-level operator' 
     261         CASE( np_blp_i  )   ;   WRITE(numout,*) '          Rotated bilaplacian operator (standard)' 
     262         CASE( np_blp_it )   ;   WRITE(numout,*) '          Rotated bilaplacian operator (triad)' 
     263         END SELECT 
     264      ENDIF 
     265      ! 
     266   END SUBROUTINE trc_ldf_ini 
    261267#else 
    262268   !!---------------------------------------------------------------------- 

  • branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/TOP_SRC/TRP/trcnxt.F90

    r4611 r6225  
    3030   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    3131   USE prtctl_trc      ! Print control for debbuging 
    32    USE trdmod_oce 
     32   USE trd_oce 
    3333   USE trdtra 
    3434   USE tranxt 
     35   USE trcbdy          ! BDY open boundaries 
     36   USE bdy_par, only: lk_bdy 
    3537# if defined key_agrif 
    3638   USE agrif_top_interp 
    … …  
    4143 
    4244   PUBLIC   trc_nxt          ! routine called by step.F90 
    43    PUBLIC   trc_nxt_alloc    ! routine called by nemogcm.F90 
    4445 
    45    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) ::   r2dt 
     46   REAL(wp) ::   r2dttrc 
    4647 
    4748   !!---------------------------------------------------------------------- 
    … …  
    5152   !!---------------------------------------------------------------------- 
    5253CONTAINS 
    53  
    54    INTEGER FUNCTION trc_nxt_alloc() 
    55       !!---------------------------------------------------------------------- 
    56       !!                   ***  ROUTINE trc_nxt_alloc  *** 
    57       !!---------------------------------------------------------------------- 
    58       ALLOCATE( r2dt(jpk), STAT=trc_nxt_alloc ) 
    59       ! 
    60       IF( trc_nxt_alloc /= 0 )   CALL ctl_warn('trc_nxt_alloc : failed to allocate array') 
    61       ! 
    62    END FUNCTION trc_nxt_alloc 
    63  
    6454 
    6555   SUBROUTINE trc_nxt( kt ) 
    … …  
    10191         WRITE(numout,*) 'trc_nxt : time stepping on passive tracers' 
    10292      ENDIF 
    103  
    104       ! Update after tracer on domain lateral boundaries 
    105       DO jn = 1, jptra 
     93      ! 
     94#if defined key_agrif 
     95      CALL Agrif_trc                   ! AGRIF zoom boundaries 
     96#endif 
     97      DO jn = 1, jptra                 ! Update after tracer on domain lateral boundaries 
    10698         CALL lbc_lnk( tra(:,:,:,jn), 'T', 1. )    
    10799      END DO 
    108100 
     101      IF( lk_bdy )  CALL trc_bdy( kt ) 
    109102 
    110 #if defined key_bdy 
    111 !!      CALL bdy_trc( kt )               ! BDY open boundaries 
    112 #endif 
    113 #if defined key_agrif 
    114       CALL Agrif_trc                   ! AGRIF zoom boundaries 
    115 #endif 
    116  
    117  
    118       ! set time step size (Euler/Leapfrog) 
    119       IF( neuler == 0 .AND. kt ==  nittrc000 ) THEN  ;  r2dt(:) =     rdttrc(:)   ! at nittrc000             (Euler) 
    120       ELSEIF( kt <= nittrc000 + 1 )            THEN  ;  r2dt(:) = 2.* rdttrc(:)   ! at nit000 or nit000+1 (Leapfrog) 
     103      !                                ! set time step size (Euler/Leapfrog) 
     104      IF( neuler == 0 .AND. kt ==  nittrc000 ) THEN  ;  r2dttrc =     rdttrc   ! at nittrc000             (Euler) 
     105      ELSEIF( kt <= nittrc000 + nn_dttrc )     THEN  ;  r2dttrc = 2.* rdttrc   ! at nit000 or nit000+1 (Leapfrog) 
    121106      ENDIF 
    122107 
    123       ! trends computation initialisation 
    124       IF( l_trdtrc )  THEN 
    125          CALL wrk_alloc( jpi, jpj, jpk, jptra, ztrdt )  !* store now fields before applying the Asselin filter 
     108      IF( l_trdtrc )  THEN             ! trends: store now fields before the Asselin filter application 
     109         CALL wrk_alloc( jpi, jpj, jpk, jptra, ztrdt ) 
    126110         ztrdt(:,:,:,:)  = trn(:,:,:,:) 
    127111      ENDIF 
    128       ! Leap-Frog + Asselin filter time stepping 
    129       IF( neuler == 0 .AND. kt == nittrc000 ) THEN        ! Euler time-stepping at first time-step 
    130          !                                                ! (only swap) 
     112      !                                ! Leap-Frog + Asselin filter time stepping 
     113      IF( neuler == 0 .AND. kt == nittrc000 ) THEN    ! Euler time-stepping at first time-step (only swap) 
    131114         DO jn = 1, jptra 
    132115            DO jk = 1, jpkm1 
    … …  
    134117            END DO 
    135118         END DO 
    136          !                                               
    137       ELSE 
    138          ! Leap-Frog + Asselin filter time stepping 
    139          IF( lk_vvl ) THEN   ;   CALL tra_nxt_vvl( kt, nittrc000, 'TRC', trb, trn, tra, jptra )      ! variable volume level (vvl)  
    140          ELSE                ;   CALL tra_nxt_fix( kt, nittrc000, 'TRC', trb, trn, tra, jptra )      ! fixed    volume level  
     119      ELSE                                            ! Asselin filter + swap 
     120         IF( ln_linssh ) THEN   ;   CALL tra_nxt_fix( kt, nittrc000,         'TRC', trb, trn, tra, jptra )  !     linear ssh 
     121         ELSE                   ;   CALL tra_nxt_vvl( kt, nittrc000, rdttrc, 'TRC', trb, trn, tra,      & 
     122           &                                                                   sbc_trc, sbc_trc_b, jptra )  ! non-linear ssh 
    141123         ENDIF 
     124         ! 
     125         DO jn = 1, jptra 
     126            CALL lbc_lnk( trb(:,:,:,jn), 'T', 1._wp )  
     127            CALL lbc_lnk( trn(:,:,:,jn), 'T', 1._wp ) 
     128            CALL lbc_lnk( tra(:,:,:,jn), 'T', 1._wp ) 
     129         END DO 
    142130      ENDIF 
    143  
    144       ! trends computation 
    145       IF( l_trdtrc ) THEN                                      ! trends 
     131      ! 
     132      IF( l_trdtrc ) THEN              ! trends: send Asselin filter trends to trdtra manager for further diagnostics 
    146133         DO jn = 1, jptra 
    147134            DO jk = 1, jpkm1 
    148                zfact = 1.e0 / r2dt(jk)   
     135               zfact = 1._wp / r2dttrc   
    149136               ztrdt(:,:,jk,jn) = ( trb(:,:,jk,jn) - ztrdt(:,:,jk,jn) ) * zfact  
    150                CALL trd_tra( kt, 'TRC', jn, jptra_trd_atf, ztrdt ) 
     137               CALL trd_tra( kt, 'TRC', jn, jptra_atf, ztrdt ) 
    151138            END DO 
    152139         END DO 

  • branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/TOP_SRC/TRP/trcrad.F90

    r3680 r6225  
    1515   USE oce_trc             ! ocean dynamics and tracers variables 
    1616   USE trc                 ! ocean passive tracers variables 
    17    USE trdmod_oce 
     17   USE trd_oce 
    1818   USE trdtra 
    1919   USE prtctl_trc          ! Print control for debbuging 
    … …  
    2222   PRIVATE 
    2323 
    24    PUBLIC trc_rad         ! routine called by trcstp.F90 
    25  
    26    !! * Substitutions 
    27 #  include "top_substitute.h90" 
     24   PUBLIC trc_rad      
     25   PUBLIC trc_rad_ini   
     26 
     27   LOGICAL , PUBLIC ::   ln_trcrad           !: flag to artificially correct negative concentrations 
     28 
    2829   !!---------------------------------------------------------------------- 
    2930   !! NEMO/TOP 3.3 , NEMO Consortium (2010) 
    … …  
    7677      ! 
    7778   END SUBROUTINE trc_rad 
     79 
     80   SUBROUTINE trc_rad_ini 
     81      !!--------------------------------------------------------------------- 
     82      !!                  ***  ROUTINE trc _rad_ini *** 
     83      !! 
     84      !! ** Purpose : read  namelist options  
     85      !!---------------------------------------------------------------------- 
     86      INTEGER ::  ios                 ! Local integer output status for namelist read 
     87      NAMELIST/namtrc_rad/ ln_trcrad 
     88      !!---------------------------------------------------------------------- 
     89 
     90      ! 
     91      REWIND( numnat_ref )              ! namtrc_rad in reference namelist  
     92      READ  ( numnat_ref, namtrc_rad, IOSTAT = ios, ERR = 907) 
     93907   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrc_rad in reference namelist', lwp ) 
     94 
     95      REWIND( numnat_cfg )              ! namtrc_rad in configuration namelist  
     96      READ  ( numnat_cfg, namtrc_rad, IOSTAT = ios, ERR = 908 ) 
     97908   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrc_rad in configuration namelist', lwp ) 
     98      IF(lwm) WRITE ( numont, namtrc_rad ) 
     99 
     100      IF(lwp) THEN                     !   ! Control print 
     101         WRITE(numout,*) 
     102         WRITE(numout,*) '   Namelist namtrc_rad : treatment of negative concentrations' 
     103         WRITE(numout,*) '      correct artificially negative concen. or not ln_trcrad = ', ln_trcrad 
     104      ENDIF 
     105      ! 
     106   END SUBROUTINE trc_rad_ini 
    78107 
    79108   SUBROUTINE trc_rad_sms( kt, ptrb, ptrn, jp_sms0, jp_sms1, cpreserv ) 
    … …  
    156185               ztrtrdb(:,:,:) = ( ptrb(:,:,:,jn) - ztrtrdb(:,:,:) ) * zs2rdt 
    157186               ztrtrdn(:,:,:) = ( ptrn(:,:,:,jn) - ztrtrdn(:,:,:) ) * zs2rdt  
    158                CALL trd_tra( kt, 'TRC', jn, jptra_trd_radb, ztrtrdb )       ! Asselin-like trend handling 
    159                CALL trd_tra( kt, 'TRC', jn, jptra_trd_radn, ztrtrdn )       ! standard     trend handling 
     187               CALL trd_tra( kt, 'TRC', jn, jptra_radb, ztrtrdb )       ! Asselin-like trend handling 
     188               CALL trd_tra( kt, 'TRC', jn, jptra_radn, ztrtrdn )       ! standard     trend handling 
    160189              ! 
    161190            ENDIF 
    … …  
    187216               ztrtrdb(:,:,:) = ( ptrb(:,:,:,jn) - ztrtrdb(:,:,:) ) * zs2rdt 
    188217               ztrtrdn(:,:,:) = ( ptrn(:,:,:,jn) - ztrtrdn(:,:,:) ) * zs2rdt  
    189                CALL trd_tra( kt, 'TRC', jn, jptra_trd_radb, ztrtrdb )       ! Asselin-like trend handling 
    190                CALL trd_tra( kt, 'TRC', jn, jptra_trd_radn, ztrtrdn )       ! standard     trend handling 
     218               CALL trd_tra( kt, 'TRC', jn, jptra_radb, ztrtrdb )       ! Asselin-like trend handling 
     219               CALL trd_tra( kt, 'TRC', jn, jptra_radn, ztrtrdn )       ! standard     trend handling 
    191220              ! 
    192221            ENDIF 

  • branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/TOP_SRC/TRP/trcsbc.F90

    r3719 r6225  
    1919   USE trc             ! ocean  passive tracers variables 
    2020   USE prtctl_trc      ! Print control for debbuging 
    21    USE trdmod_oce 
     21   USE iom 
     22   USE trd_oce 
    2223   USE trdtra 
    2324 
    … …  
    2829 
    2930   !! * Substitutions 
    30 #  include "top_substitute.h90" 
     31#  include "vectopt_loop_substitute.h90" 
    3132   !!---------------------------------------------------------------------- 
    3233   !! NEMO/TOP 3.3 , NEMO Consortium (2010) 
    … …  
    6061      INTEGER, INTENT( in ) ::   kt          ! ocean time-step index 
    6162      ! 
    62       INTEGER  ::   ji, jj, jn           ! dummy loop indices 
    63       REAL(wp) ::   zsrau, zse3t   ! temporary scalars 
     63      INTEGER  ::   ji, jj, jn                                     ! dummy loop indices 
     64      REAL(wp) ::   zse3t, zrtrn, zratio, zfact                    ! temporary scalars 
     65      REAL(wp) ::   zswitch, zftra, zcd, zdtra, ztfx, ztra         ! temporary scalars 
    6466      CHARACTER (len=22) :: charout 
    6567      REAL(wp), POINTER, DIMENSION(:,:  ) :: zsfx 
    6668      REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrtrd 
     69 
    6770      !!--------------------------------------------------------------------- 
    6871      ! 
    … …  
    7073      ! 
    7174      ! Allocate temporary workspace 
    72                       CALL wrk_alloc( jpi, jpj,      zsfx   ) 
    73       IF( l_trdtrc )  CALL wrk_alloc( jpi, jpj, jpk, ztrtrd ) 
     75                      CALL wrk_alloc( jpi,jpj,       zsfx   ) 
     76      IF( l_trdtrc )  CALL wrk_alloc( jpi,jpj,jpk,   ztrtrd ) 
     77      ! 
     78      zrtrn = 1.e-15_wp 
     79 
     80      SELECT CASE( nn_ice_embd )         ! levitating or embedded sea-ice option 
     81         CASE( 0    )   ;   zswitch = 1  ! (0) standard levitating sea-ice : salt exchange only 
     82         CASE( 1, 2 )   ;   zswitch = 0  ! (1) levitating sea-ice: salt and volume exchange but no pressure effect                                 
     83      !                                  ! (2) embedded sea-ice : salt and volume fluxes and pressure 
     84      END SELECT 
     85 
     86      IF( ln_top_euler) THEN 
     87         r2dt =  rdttrc              ! = rdttrc (use Euler time stepping) 
     88      ELSE 
     89         IF( neuler == 0 .AND. kt == nittrc000 ) THEN     ! at nittrc000 
     90            r2dt = rdttrc            ! = rdttrc (restarting with Euler time stepping) 
     91         ELSEIF( kt <= nittrc000 + nn_dttrc ) THEN          ! at nittrc000 or nittrc000+1 
     92            r2dt = 2. * rdttrc       ! = 2 rdttrc (leapfrog) 
     93         ENDIF 
     94      ENDIF 
     95 
    7496 
    7597      IF( kt == nittrc000 ) THEN 
    … …  
    7799         IF(lwp) WRITE(numout,*) 'trc_sbc : Passive tracers surface boundary condition' 
    78100         IF(lwp) WRITE(numout,*) '~~~~~~~ ' 
     101 
     102         IF( ln_rsttr .AND.    &                     ! Restart: read in restart  file 
     103            iom_varid( numrtr, 'sbc_'//TRIM(ctrcnm(1))//'_b', ldstop = .FALSE. ) > 0 ) THEN 
     104            IF(lwp) WRITE(numout,*) '          nittrc000-nn_dttrc surface tracer content forcing fields red in the restart file' 
     105            zfact = 0.5_wp 
     106            DO jn = 1, jptra 
     107               CALL iom_get( numrtr, jpdom_autoglo, 'sbc_'//TRIM(ctrcnm(jn))//'_b', sbc_trc_b(:,:,jn) )   ! before tracer content sbc 
     108            END DO 
     109         ELSE                                         ! No restart or restart not found: Euler forward time stepping 
     110           zfact = 1._wp 
     111           sbc_trc_b(:,:,:) = 0._wp 
     112         ENDIF 
     113      ELSE                                         ! Swap of forcing fields 
     114         IF( ln_top_euler ) THEN 
     115            zfact = 1._wp 
     116            sbc_trc_b(:,:,:) = 0._wp 
     117         ELSE 
     118            zfact = 0.5_wp 
     119            sbc_trc_b(:,:,:) = sbc_trc(:,:,:) 
     120         ENDIF 
     121         ! 
    79122      ENDIF 
    80123 
    … …  
    83126      ! Coupling offline : runoff are in emp which contains E-P-R 
    84127      ! 
    85       IF( .NOT. lk_offline .AND. lk_vvl ) THEN  ! online coupling with vvl 
     128      IF( .NOT. lk_offline .AND. .NOT.ln_linssh ) THEN  ! online coupling with vvl 
    86129         zsfx(:,:) = 0._wp 
    87130      ELSE                                      ! online coupling free surface or offline with free surface 
    … …  
    90133 
    91134      ! 0. initialization 
    92       zsrau = 1. / rau0 
    93135      DO jn = 1, jptra 
    94136         ! 
    95          IF( l_trdtrc ) ztrtrd(:,:,:) = tra(:,:,:,jn)  ! save trends 
    96          !                                             ! add the trend to the general tracer trend 
     137         IF( l_trdtrc )   ztrtrd(:,:,:) = tra(:,:,:,jn)  ! save trends 
     138 
     139         IF ( nn_ice_tr == -1 ) THEN  ! No tracers in sea ice (null concentration in sea ice) 
     140 
     141            DO jj = 2, jpj 
     142               DO ji = fs_2, fs_jpim1   ! vector opt. 
     143                  sbc_trc(ji,jj,jn) = zsfx(ji,jj) * r1_rau0 * trn(ji,jj,1,jn) 
     144               END DO 
     145            END DO 
     146 
     147         ELSE 
     148 
     149            DO jj = 2, jpj 
     150               DO ji = fs_2, fs_jpim1   ! vector opt. 
     151                  zse3t = 1. / e3t_n(ji,jj,1) 
     152                  ! tracer flux at the ice/ocean interface (tracer/m2/s) 
     153                  zftra = - trc_i(ji,jj,jn) * fmmflx(ji,jj) ! uptake of tracer in the sea ice 
     154                  zcd   =   trc_o(ji,jj,jn) * fmmflx(ji,jj) ! concentration dilution due to freezing-melting, 
     155                                                               ! only used in the levitating sea ice case 
     156                  ! tracer flux only       : add concentration dilution term in net tracer flux, no F-M in volume flux 
     157                  ! tracer and mass fluxes : no concentration dilution term in net tracer flux, F-M term in volume flux 
     158                  ztfx  = zftra + zswitch * zcd                ! net tracer flux (+C/D if no ice/ocean mass exchange) 
     159    
     160                  zdtra = r1_rau0 * ( ztfx + zsfx(ji,jj) * trn(ji,jj,1,jn) )  
     161                  IF ( zdtra < 0. ) THEN 
     162                     zratio = -zdtra * zse3t * r2dt / ( trn(ji,jj,1,jn) + zrtrn ) 
     163                     zdtra = MIN(1.0, zratio) * zdtra ! avoid negative concentrations to arise 
     164                  ENDIF 
     165                  sbc_trc(ji,jj,jn) =  zdtra  
     166               END DO 
     167            END DO 
     168         ENDIF 
     169         !                                       Concentration dilution effect on tracers due to evaporation & precipitation  
    97170         DO jj = 2, jpj 
    98171            DO ji = fs_2, fs_jpim1   ! vector opt. 
    99                zse3t = 1. / fse3t(ji,jj,1) 
    100                tra(ji,jj,1,jn) = tra(ji,jj,1,jn) + zsfx(ji,jj) *  zsrau * trn(ji,jj,1,jn) * zse3t 
     172               zse3t = zfact / e3t_n(ji,jj,1) 
     173               tra(ji,jj,1,jn) = tra(ji,jj,1,jn) + ( sbc_trc_b(ji,jj,jn) + sbc_trc(ji,jj,jn) ) * zse3t 
    101174            END DO 
    102175         END DO 
    103           
     176         ! 
    104177         IF( l_trdtrc ) THEN 
    105178            ztrtrd(:,:,:) = tra(:,:,:,jn) - ztrtrd(:,:,:) 
    106             CALL trd_tra( kt, 'TRC', jn, jptra_trd_nsr, ztrtrd ) 
     179            CALL trd_tra( kt, 'TRC', jn, jptra_nsr, ztrtrd ) 
    107180         END IF 
    108181         !                                                       ! =========== 
    109182      END DO                                                     ! tracer loop 
    110183      !                                                          ! =========== 
     184 
     185      !                                           Write in the tracer restar  file 
     186      !                                          ******************************* 
     187      IF( lrst_trc ) THEN 
     188         IF(lwp) WRITE(numout,*) 
     189         IF(lwp) WRITE(numout,*) 'sbc : ocean surface tracer content forcing fields written in tracer restart file ',   & 
     190            &                    'at it= ', kt,' date= ', ndastp 
     191         IF(lwp) WRITE(numout,*) '~~~~' 
     192         DO jn = 1, jptra 
     193            CALL iom_rstput( kt, nitrst, numrtw, 'sbc_'//TRIM(ctrcnm(jn))//'_b', sbc_trc(:,:,jn) ) 
     194         END DO 
     195      ENDIF 
     196      ! 
    111197      IF( ln_ctl )   THEN 
    112198         WRITE(charout, FMT="('sbc ')") ;  CALL prt_ctl_trc_info(charout) 
    113199                                           CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 
    114200      ENDIF 
    115                       CALL wrk_dealloc( jpi, jpj,      zsfx   ) 
    116       IF( l_trdtrc )  CALL wrk_dealloc( jpi, jpj, jpk, ztrtrd ) 
     201                      CALL wrk_dealloc( jpi,jpj,       zsfx   ) 
     202      IF( l_trdtrc )  CALL wrk_dealloc( jpi,jpj,jpk,  ztrtrd ) 
    117203      ! 
    118204      IF( nn_timing == 1 )  CALL timing_stop('trc_sbc') 

  • branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/TOP_SRC/TRP/trctrp.F90

    r4148 r6225  
    1515   USE oce_trc         ! ocean dynamics and active tracers variables 
    1616   USE trc             ! ocean passive tracers variables  
    17    USE trcnam_trp      ! passive tracers transport namelist variables 
    1817   USE trabbl          ! bottom boundary layer               (trc_bbl routine) 
    1918   USE trcbbl          ! bottom boundary layer               (trc_bbl routine) 
    20    USE zdfkpp          ! KPP non-local tracer fluxes         (trc_kpp routine) 
    2119   USE trcdmp          ! internal damping                    (trc_dmp routine) 
    2220   USE trcldf          ! lateral mixing                      (trc_ldf routine) 
    … …  
    2725   USE trcsbc          ! surface boundary condition          (trc_sbc routine) 
    2826   USE zpshde          ! partial step: hor. derivative       (zps_hde routine) 
     27   USE trcbdy          ! BDY open boundaries 
     28   USE bdy_par, only: lk_bdy 
    2929 
    3030#if defined key_agrif 
    … …  
    3838   PUBLIC   trc_trp    ! called by trc_stp 
    3939 
    40    !! * Substitutions 
    41 #  include "top_substitute.h90" 
    4240   !!---------------------------------------------------------------------- 
    4341   !! NEMO/TOP 3.3 , NEMO Consortium (2010) 
    … …  
    4846CONTAINS 
    4947 
    50    SUBROUTINE trc_trp( kstp ) 
     48   SUBROUTINE trc_trp( kt ) 
    5149      !!---------------------------------------------------------------------- 
    5250      !!                     ***  ROUTINE trc_trp  *** 
    … …  
    5755      !!              - Update the passive tracers 
    5856      !!---------------------------------------------------------------------- 
    59       INTEGER, INTENT( in ) ::  kstp  ! ocean time-step index 
     57      INTEGER, INTENT( in ) ::  kt  ! ocean time-step index 
    6058      !! --------------------------------------------------------------------- 
    6159      ! 
    … …  
    6462      IF( .NOT. lk_c1d ) THEN 
    6563         ! 
    66                                 CALL trc_sbc( kstp )            ! surface boundary condition 
    67          IF( lk_trabbl )        CALL trc_bbl( kstp )            ! advective (and/or diffusive) bottom boundary layer scheme 
    68          IF( ln_trcdmp )        CALL trc_dmp( kstp )            ! internal damping trends 
    69          IF( ln_trcdmp_clo )    CALL trc_dmp_clo( kstp )        ! internal damping trends on closed seas only 
    70                                 CALL trc_adv( kstp )            ! horizontal & vertical advection  
    71                                 CALL trc_ldf( kstp )            ! lateral mixing 
    72          IF( .NOT. lk_offline .AND. lk_zdfkpp )    & 
    73             &                   CALL trc_kpp( kstp )            ! KPP non-local tracer fluxes 
     64                                CALL trc_sbc    ( kt )      ! surface boundary condition 
     65         IF( lk_trabbl )        CALL trc_bbl    ( kt )      ! advective (and/or diffusive) bottom boundary layer scheme 
     66         IF( ln_trcdmp )        CALL trc_dmp    ( kt )      ! internal damping trends 
     67         IF( ln_trcdmp_clo )    CALL trc_dmp_clo( kt )      ! internal damping trends on closed seas only 
     68         IF( lk_bdy )           CALL trc_bdy_dmp( kt )      ! BDY damping trends 
     69                                CALL trc_adv    ( kt )      ! horizontal & vertical advection  
     70         !                                                         ! Partial top/bottom cell: GRADh( trb )   
     71         IF( ln_zps ) THEN 
     72           IF( ln_isfcav ) THEN ; CALL zps_hde_isf( kt, jptra, trb, pgtu=gtru, pgtv=gtrv, pgtui=gtrui, pgtvi=gtrvi )  ! both top & bottom 
     73           ELSE                 ; CALL zps_hde    ( kt, jptra, trb, gtru, gtrv )                                      !  only bottom 
     74           ENDIF 
     75         ENDIF 
     76         !                                                       
     77                                CALL trc_ldf    ( kt )      ! lateral mixing 
    7478#if defined key_agrif 
    75          IF(.NOT. Agrif_Root()) CALL Agrif_Sponge_trc           ! tracers sponge 
     79         IF(.NOT. Agrif_Root()) CALL Agrif_Sponge_trc       ! tracers sponge 
    7680#endif 
    77                                 CALL trc_zdf( kstp )            ! vertical mixing and after tracer fields 
    78                                 CALL trc_nxt( kstp )            ! tracer fields at next time step      
    79          IF( ln_trcrad )        CALL trc_rad( kstp )            ! Correct artificial negative concentrations 
     81                                CALL trc_zdf    ( kt )      ! vertical mixing and after tracer fields 
     82                                CALL trc_nxt    ( kt )      ! tracer fields at next time step      
     83         IF( ln_trcrad )        CALL trc_rad    ( kt )      ! Correct artificial negative concentrations 
    8084 
    8185#if defined key_agrif 
    82       IF( .NOT. Agrif_Root())   CALL Agrif_Update_Trc( kstp )  ! Update tracer at AGRIF zoom boundaries : children only 
     86         IF( .NOT.Agrif_Root()) CALL Agrif_Update_Trc( kt ) ! Update tracer at AGRIF zoom boundaries : children only 
    8387#endif 
    84          IF( ln_zps    )        CALL zps_hde( kstp, jptra, trn, gtru, gtrv )  ! Partial steps: now horizontal gradient of passive 
    85                                                                 ! tracers at the bottom ocean level 
    8688         ! 
    8789      ELSE                                               ! 1D vertical configuration 
    88                                 CALL trc_sbc( kstp )            ! surface boundary condition 
    89          IF( .NOT. lk_offline .AND. lk_zdfkpp )    & 
    90             &                   CALL trc_kpp( kstp )            ! KPP non-local tracer fluxes 
    91                                 CALL trc_zdf( kstp )            ! vertical mixing and after tracer fields 
    92                                 CALL trc_nxt( kstp )            ! tracer fields at next time step      
    93           IF( ln_trcrad )       CALL trc_rad( kstp )            ! Correct artificial negative concentrations 
     90                                CALL trc_sbc( kt )            ! surface boundary condition 
     91         IF( ln_trcdmp )        CALL trc_dmp( kt )            ! internal damping trends 
     92                                CALL trc_zdf( kt )            ! vertical mixing and after tracer fields 
     93                                CALL trc_nxt( kt )            ! tracer fields at next time step      
     94          IF( ln_trcrad )       CALL trc_rad( kt )            ! Correct artificial negative concentrations 
    9495         ! 
    9596      END IF 
    … …  
    104105   !!---------------------------------------------------------------------- 
    105106CONTAINS 
    106    SUBROUTINE trc_trp( kstp )              ! Empty routine 
    107       INTEGER, INTENT(in) ::   kstp 
    108       WRITE(*,*) 'trc_trp: You should not have seen this print! error?', kstp 
     107   SUBROUTINE trc_trp( kt )              ! Empty routine 
     108      INTEGER, INTENT(in) ::   kt 
     109      WRITE(*,*) 'trc_trp: You should not have seen this print! error?', kt 
    109110   END SUBROUTINE trc_trp 
    110111#endif 

  • branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/TOP_SRC/TRP/trczdf.F90

    r3680 r6225  
    1111   !!   'key_top'                                                TOP models 
    1212   !!---------------------------------------------------------------------- 
    13    !!   trc_ldf     : update the tracer trend with the lateral diffusion 
    14    !!       ldf_ctl : initialization, namelist read, and parameters control 
     13   !!   trc_zdf      : update the tracer trend with the lateral diffusion 
     14   !!   trc_zdf_ini : initialization, namelist read, and parameters control 
    1515   !!---------------------------------------------------------------------- 
    16    USE oce_trc         ! ocean dynamics and active tracers 
    17    USE trc             ! ocean passive tracers variables 
    18    USE trcnam_trp      ! passive tracers transport namelist variables 
    19    USE trazdf_exp      ! vertical diffusion: explicit (tra_zdf_exp     routine) 
    20    USE trazdf_imp      ! vertical diffusion: implicit (tra_zdf_imp     routine) 
    21    USE trdmod_oce 
    22    USE trdtra 
    23    USE prtctl_trc      ! Print control 
     16   USE trc           ! ocean passive tracers variables 
     17   USE oce_trc       ! ocean dynamics and active tracers 
     18   USE trd_oce       ! trends: ocean variables 
     19   USE trazdf_exp    ! vertical diffusion: explicit (tra_zdf_exp     routine) 
     20   USE trazdf_imp    ! vertical diffusion: implicit (tra_zdf_imp     routine) 
     21   USE trcldf        ! passive tracers: lateral diffusion 
     22   USE trdtra        ! trends manager: tracers  
     23   USE prtctl_trc    ! Print control 
    2424 
    2525   IMPLICIT NONE 
    2626   PRIVATE 
    2727 
    28    PUBLIC   trc_zdf          ! called by step.F90  
    29    PUBLIC   trc_zdf_alloc    ! called by nemogcm.F90  
     28   PUBLIC   trc_zdf         ! called by step.F90  
     29   PUBLIC   trc_zdf_ini     ! called by nemogcm.F90  
     30    
     31   !                                        !!** Vertical diffusion (nam_trczdf) ** 
     32   LOGICAL , PUBLIC ::   ln_trczdf_exp       !: explicit vertical diffusion scheme flag 
     33   INTEGER , PUBLIC ::   nn_trczdf_exp       !: number of sub-time step (explicit time stepping) 
    3034 
    3135   INTEGER ::   nzdf = 0               ! type vertical diffusion algorithm used 
    3236      !                                ! defined from ln_zdf...  namlist logicals) 
    33    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) ::  r2dt   ! vertical profile time-step, = 2 rdttra 
    34       !                                                 ! except at nittrc000 (=rdttra) if neuler=0 
     37   REAL(wp) ::  r2dttrc   ! vertical profile time-step, = 2 rdt 
     38      !                   ! except at nittrc000 (=rdt) if neuler=0 
    3539 
    3640   !! * Substitutions 
    37 #  include "domzgr_substitute.h90" 
    3841#  include "zdfddm_substitute.h90" 
    3942#  include "vectopt_loop_substitute.h90" 
    4043   !!---------------------------------------------------------------------- 
    41    !! NEMO/TOP 3.3 , NEMO Consortium (2010) 
     44   !! NEMO/TOP 3.7 , NEMO Consortium (2015) 
    4245   !! $Id$  
    4346   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    4447   !!---------------------------------------------------------------------- 
    4548CONTAINS 
    46     
    47    INTEGER FUNCTION trc_zdf_alloc() 
    48       !!---------------------------------------------------------------------- 
    49       !!                  ***  ROUTINE trc_zdf_alloc  *** 
    50       !!---------------------------------------------------------------------- 
    51       ALLOCATE( r2dt(jpk) , STAT=trc_zdf_alloc ) 
    52       ! 
    53       IF( trc_zdf_alloc /= 0 )   CALL ctl_warn('trc_zdf_alloc : failed to allocate array.') 
    54       ! 
    55    END FUNCTION trc_zdf_alloc 
    56  
    5749 
    5850   SUBROUTINE trc_zdf( kt ) 
    … …  
    7163      IF( nn_timing == 1 )  CALL timing_start('trc_zdf') 
    7264      ! 
    73       IF( kt == nittrc000 )   CALL zdf_ctl          ! initialisation & control of options 
    74  
    75       IF( ln_top_euler) THEN 
    76          r2dt(:) =  rdttrc(:)              ! = rdttrc (use Euler time stepping) 
    77       ELSE 
    78          IF( neuler == 0 .AND. kt == nittrc000 ) THEN     ! at nittrc000 
    79             r2dt(:) =  rdttrc(:)           ! = rdttrc (restarting with Euler time stepping) 
    80          ELSEIF( kt <= nittrc000 + 1 ) THEN          ! at nittrc000 or nittrc000+1 
    81             r2dt(:) = 2. * rdttrc(:)       ! = 2 rdttrc (leapfrog) 
    82          ENDIF 
     65      IF( ( neuler == 0 .AND. kt == nittrc000 ) .OR. ln_top_euler ) THEN     ! at nittrc000 
     66         r2dttrc =  rdttrc           ! = rdttrc (use or restarting with Euler time stepping) 
     67      ELSEIF( kt <= nittrc000 + nn_dttrc ) THEN          ! at nittrc000 or nittrc000+1 
     68         r2dttrc = 2. * rdttrc       ! = 2 rdttrc (leapfrog) 
    8369      ENDIF 
    8470 
    … …  
    8975 
    9076      SELECT CASE ( nzdf )                       ! compute lateral mixing trend and add it to the general trend 
    91       CASE ( -1 )                                       ! esopa: test all possibility with control print 
    92          CALL tra_zdf_exp( kt, nittrc000, 'TRC', r2dt, nn_trczdf_exp, trb, tra, jptra )  
    93          WRITE(charout, FMT="('zdf1 ')") ;  CALL prt_ctl_trc_info(charout) 
    94                                             CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 
    95          CALL tra_zdf_imp( kt, nittrc000, 'TRC', r2dt,                trb, tra, jptra )  
    96          WRITE(charout, FMT="('zdf2 ')") ;  CALL prt_ctl_trc_info(charout) 
    97                                             CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' ) 
    98       CASE ( 0 ) ;  CALL tra_zdf_exp( kt, nittrc000, 'TRC', r2dt, nn_trczdf_exp, trb, tra, jptra )    !   explicit scheme  
    99       CASE ( 1 ) ;  CALL tra_zdf_imp( kt, nittrc000, 'TRC', r2dt,                trb, tra, jptra )    !   implicit scheme           
    100  
     77      CASE ( 0 ) ;  CALL tra_zdf_exp( kt, nittrc000, 'TRC', r2dttrc, nn_trczdf_exp, trb, tra, jptra )    !   explicit scheme  
     78      CASE ( 1 ) ;  CALL tra_zdf_imp( kt, nittrc000, 'TRC', r2dttrc,                trb, tra, jptra )    !   implicit scheme           
    10179      END SELECT 
    10280 
    … …  
    10482         DO jn = 1, jptra 
    10583            DO jk = 1, jpkm1 
    106                ztrtrd(:,:,jk,jn) = ( ( tra(:,:,jk,jn) - trb(:,:,jk,jn) ) / r2dt(jk) ) - ztrtrd(:,:,jk,jn) 
     84               ztrtrd(:,:,jk,jn) = ( ( tra(:,:,jk,jn) - trb(:,:,jk,jn) ) / r2dttrc ) - ztrtrd(:,:,jk,jn) 
    10785            END DO 
    108             CALL trd_tra( kt, 'TRC', jn, jptra_trd_zdf, ztrtrd(:,:,:,jn) ) 
     86            CALL trd_tra( kt, 'TRC', jn, jptra_zdf, ztrtrd(:,:,:,jn) ) 
    10987         END DO 
    11088         CALL wrk_dealloc( jpi, jpj, jpk, jptra, ztrtrd ) 
    … …  
    12199 
    122100 
    123    SUBROUTINE zdf_ctl 
     101   SUBROUTINE trc_zdf_ini 
    124102      !!---------------------------------------------------------------------- 
    125       !!                 ***  ROUTINE zdf_ctl  *** 
     103      !!                 ***  ROUTINE trc_zdf_ini  *** 
    126104      !! 
    127105      !! ** Purpose :   Choose the vertical mixing scheme 
    … …  
    132110      !!      NB: The implicit scheme is required when using :  
    133111      !!             - rotated lateral mixing operator 
    134       !!             - TKE, GLS or KPP vertical mixing scheme 
     112      !!             - TKE, GLS vertical mixing scheme 
    135113      !!---------------------------------------------------------------------- 
    136  
    137       !  Define the vertical tracer physics scheme 
    138       ! ========================================== 
    139  
    140       ! Choice from ln_zdfexp already read in namelist in zdfini module 
    141       IF( ln_trczdf_exp ) THEN           ! use explicit scheme 
    142          nzdf = 0 
    143       ELSE                               ! use implicit scheme 
    144          nzdf = 1 
     114      INTEGER ::  ios                 ! Local integer output status for namelist read 
     115      !! 
     116      NAMELIST/namtrc_zdf/ ln_trczdf_exp  , nn_trczdf_exp 
     117      !!---------------------------------------------------------------------- 
     118      ! 
     119      REWIND( numnat_ref )             ! namtrc_zdf in reference namelist  
     120      READ  ( numnat_ref, namtrc_zdf, IOSTAT = ios, ERR = 905) 
     121905   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrc_zdf in reference namelist', lwp ) 
     122      ! 
     123      REWIND( numnat_cfg )             ! namtrc_zdf in configuration namelist  
     124      READ  ( numnat_cfg, namtrc_zdf, IOSTAT = ios, ERR = 906 ) 
     125906   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrc_zdf in configuration namelist', lwp ) 
     126      IF(lwm) WRITE ( numont, namtrc_zdf ) 
     127      ! 
     128      IF(lwp) THEN                     ! Control print 
     129         WRITE(numout,*) 
     130         WRITE(numout,*) '   Namelist namtrc_zdf : set vertical diffusion  parameters' 
     131         WRITE(numout,*) '      time splitting / backward scheme ln_trczdf_exp = ', ln_trczdf_exp 
     132         WRITE(numout,*) '      number of time step              nn_trczdf_exp = ', nn_trczdf_exp 
    145133      ENDIF 
    146134 
    147       ! Force implicit schemes 
    148       IF( ln_trcldf_iso                               )   nzdf = 1      ! iso-neutral lateral physics 
    149       IF( ln_trcldf_hor .AND. ln_sco                  )   nzdf = 1      ! horizontal lateral physics in s-coordinate 
    150 #if defined key_zdftke || defined key_zdfgls || defined key_zdfkpp 
    151                                                           nzdf = 1      ! TKE, GLS or KPP physics        
    152 #endif 
    153       IF( ln_trczdf_exp .AND. nzdf == 1 )   THEN 
    154          CALL ctl_stop( 'trc_zdf : If using the rotated lateral mixing operator or TKE, GLS or KPP vertical scheme ', & 
    155             &           '          the implicit scheme is required, set ln_trczdf_exp = .false.' ) 
     135      !                                ! Define the vertical tracer physics scheme 
     136      IF( ln_trczdf_exp ) THEN   ;   nzdf = 0     ! explicit scheme 
     137      ELSE                       ;   nzdf = 1     ! implicit scheme 
    156138      ENDIF 
    157139 
    158       ! Test: esopa 
    159       IF( lk_esopa )    nzdf = -1                      ! All schemes used 
     140      !                                ! Force implicit schemes 
     141      IF( ln_trcldf_iso              )   nzdf = 1      ! iso-neutral lateral physics 
     142      IF( ln_trcldf_hor .AND. ln_sco )   nzdf = 1      ! horizontal lateral physics in s-coordinate 
     143#if defined key_zdftke || defined key_zdfgls  
     144                                         nzdf = 1      ! TKE or GLS physics        
     145#endif 
     146      IF( ln_trczdf_exp .AND. nzdf == 1 )  &  
     147         CALL ctl_stop( 'trc_zdf : If using the rotated lateral mixing operator or TKE, GLS vertical scheme ', & 
     148            &           '          the implicit scheme is required, set ln_trczdf_exp = .false.' ) 
    160149 
    161150      IF(lwp) THEN 
    … …  
    163152         WRITE(numout,*) 'trc:zdf_ctl : vertical passive tracer physics scheme' 
    164153         WRITE(numout,*) '~~~~~~~~~~~' 
    165          IF( nzdf == -1 )   WRITE(numout,*) '              ESOPA test All scheme used' 
    166154         IF( nzdf ==  0 )   WRITE(numout,*) '              Explicit time-splitting scheme' 
    167155         IF( nzdf ==  1 )   WRITE(numout,*) '              Implicit (euler backward) scheme' 
    168156      ENDIF 
    169  
    170    END SUBROUTINE zdf_ctl 
     157      ! 
     158   END SUBROUTINE trc_zdf_ini 
     159    
    171160#else 
    172161   !!---------------------------------------------------------------------- 
Note: See TracChangeset for help on using the changeset viewer.