Index: trunk/NEMOGCM/NEMO/TOP_SRC/TRP/trcldf.F90
===================================================================
--- trunk/NEMOGCM/NEMO/TOP_SRC/TRP/trcldf.F90	(revision 2715)
+++ trunk/NEMOGCM/NEMO/TOP_SRC/TRP/trcldf.F90	(revision 3294)
@@ -2,8 +2,8 @@
    !!======================================================================
    !!                       ***  MODULE  trcldf  ***
-   !! Ocean Passive tracers : lateral diffusive trends 
+   !! Ocean Passive tracers : lateral diffusive trends
    !!=====================================================================
    !! History :  9.0  ! 2005-11 (G. Madec)  Original code
-   !!       NEMO 3.0  ! 2008-01  (C. Ethe, G. Madec)  merge TRC-TRA 
+   !!       NEMO 3.0  ! 2008-01  (C. Ethe, G. Madec)  merge TRC-TRA
    !!----------------------------------------------------------------------
 #if defined key_top
@@ -23,4 +23,5 @@
    USE traldf_bilap    ! lateral mixing            (tra_ldf_bilap routine)
    USE traldf_iso      ! lateral mixing            (tra_ldf_iso routine)
+   USE traldf_iso_grif ! lateral mixing          (tra_ldf_iso_grif routine)
    USE traldf_lap      ! lateral mixing            (tra_ldf_lap routine)
    USE trdmod_oce
@@ -31,7 +32,8 @@
    PRIVATE
 
-   PUBLIC   trc_ldf    ! called by step.F90 
+   PUBLIC   trc_ldf    ! called by step.F90
    !                                                 !!: ** lateral mixing namelist (nam_trcldf) **
-   INTEGER ::   nldf = 0   ! type of lateral diffusion used defined from ln_trcldf_... namlist logicals)
+   REAL(wp) ::  rldf_rat    ! ratio between active and passive tracers diffusive coefficient
+   INTEGER  ::  nldf = 0   ! type of lateral diffusion used defined from ln_trcldf_... namlist logicals)
    !! * Substitutions
 #  include "domzgr_substitute.h90"
@@ -39,5 +41,5 @@
    !!----------------------------------------------------------------------
    !! NEMO/TOP 3.3 , NEMO Consortium (2010)
-   !! $Id$ 
+   !! $Id$
    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
    !!----------------------------------------------------------------------
@@ -48,5 +50,5 @@
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE tra_ldf  ***
-      !! 
+      !!
       !! ** Purpose :   compute the lateral ocean tracer physics.
       !!
@@ -56,31 +58,44 @@
       INTEGER            :: jn
       CHARACTER (len=22) :: charout
-      REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE ::   ztrtrd
-      !!----------------------------------------------------------------------
-
-      IF( kt == nit000 )   CALL ldf_ctl          ! initialisation & control of options
-
-      IF( l_trdtrc )  THEN 
-         ALLOCATE( ztrtrd(jpi,jpj,jpk,jptra) )  ! temporary save of trends
+      REAL(wp), POINTER, DIMENSION(:,:,:,:) ::   ztrtrd
+      !!----------------------------------------------------------------------
+      !
+      IF( nn_timing == 1 )   CALL timing_start('trc_ldf')
+      !
+      IF( kt == nittrc000 )   CALL ldf_ctl          ! initialisation & control of options
+
+      rldf = rldf_rat
+
+      IF( l_trdtrc )  THEN
+         CALL wrk_alloc( jpi, jpj, jpk, jptra, ztrtrd )
          ztrtrd(:,:,:,:)  = tra(:,:,:,:)
       ENDIF
 
       SELECT CASE ( nldf )                       ! compute lateral mixing trend and add it to the general trend
-      CASE ( 0 )   ;   CALL tra_ldf_lap   ( kt, 'TRC', gtru, gtrv, trb, tra, jptra            )  ! iso-level laplacian
-      CASE ( 1 )   ;   CALL tra_ldf_iso   ( kt, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 )  ! rotated laplacian 
-      CASE ( 2 )   ;   CALL tra_ldf_bilap ( kt, 'TRC', gtru, gtrv, trb, tra, jptra            )  ! iso-level bilaplacian
-      CASE ( 3 )   ;   CALL tra_ldf_bilapg( kt, 'TRC',             trb, tra, jptra            )  ! s-coord. horizontal bilaplacian
+      CASE ( 0 )   ;   CALL tra_ldf_lap   ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra            )  ! iso-level laplacian
+      CASE ( 1 )                                                                                            ! rotated laplacian
+                       IF( ln_traldf_grif ) THEN
+                          CALL tra_ldf_iso_grif( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 )
+                       ELSE
+                          CALL tra_ldf_iso     ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 )
+                       ENDIF
+      CASE ( 2 )   ;   CALL tra_ldf_bilap ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra            )  ! iso-level bilaplacian
+      CASE ( 3 )   ;   CALL tra_ldf_bilapg( kt, nittrc000, 'TRC',             trb, tra, jptra            )  ! s-coord. horizontal bilaplacian
          !
       CASE ( -1 )                                     ! esopa: test all possibility with control print
-         CALL tra_ldf_lap   ( kt, 'TRC', gtru, gtrv, trb, tra, jptra            )
+         CALL tra_ldf_lap   ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra            )
          WRITE(charout, FMT="('ldf0 ')") ;  CALL prt_ctl_trc_info(charout)
                                             CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
-         CALL tra_ldf_iso   ( kt, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 )
+         IF( ln_traldf_grif ) THEN
+            CALL tra_ldf_iso_grif( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 )
+         ELSE
+            CALL tra_ldf_iso     ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra, rn_ahtb_0 )
+         ENDIF
          WRITE(charout, FMT="('ldf1 ')") ;  CALL prt_ctl_trc_info(charout)
                                             CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
-         CALL tra_ldf_bilap ( kt, 'TRC', gtru, gtrv, trb, tra, jptra            )
+         CALL tra_ldf_bilap ( kt, nittrc000, 'TRC', gtru, gtrv, trb, tra, jptra            )
          WRITE(charout, FMT="('ldf2 ')") ;  CALL prt_ctl_trc_info(charout)
                                             CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
-         CALL tra_ldf_bilapg( kt, 'TRC',             trb, tra, jptra            )
+         CALL tra_ldf_bilapg( kt, nittrc000, 'TRC',             trb, tra, jptra            )
          WRITE(charout, FMT="('ldf3 ')") ;  CALL prt_ctl_trc_info(charout)
                                             CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
@@ -92,5 +107,5 @@
            CALL trd_tra( kt, 'TRC', jn, jptra_trd_ldf, ztrtrd(:,:,:,jn) )
         END DO
-        DEALLOCATE( ztrtrd ) 
+        CALL wrk_dealloc( jpi, jpj, jpk, jptra, ztrtrd )
       ENDIF
       !                                          ! print mean trends (used for debugging)
@@ -100,4 +115,6 @@
       ENDIF
       !
+      IF( nn_timing == 1 )   CALL timing_stop('trc_ldf')
+      !
    END SUBROUTINE trc_ldf
 
@@ -106,9 +123,9 @@
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE ldf_ctl  ***
-      !! 
+      !!
       !! ** Purpose :   Choice of the operator for the lateral tracer diffusion
       !!
       !! ** Method  :   set nldf from the namtra_ldf logicals
-      !!      nldf == -2   No lateral diffusion  
+      !!      nldf == -2   No lateral diffusion
       !!      nldf == -1   ESOPA test: ALL operators are used
       !!      nldf ==  0   laplacian operator
@@ -117,10 +134,19 @@
       !!      nldf ==  3   Rotated bilaplacian
       !!----------------------------------------------------------------------
-      INTEGER ::   ioptio, ierr         ! temporary integers 
-      !!----------------------------------------------------------------------
-
+      INTEGER ::   ioptio, ierr         ! temporary integers
+      !!----------------------------------------------------------------------
+
+      IF (ABS(rn_aht_0) < 2._wp*TINY(1.e0)) THEN
+         IF (ABS(rn_ahtrc_0) < 2._wp*TINY(1.e0)) THEN
+            rldf_rat = 1.0_wp
+         ELSE
+            CALL ctl_stop( 'STOP', 'ldf_ctl : cannot define rldf_rat, rn_aht_0==0, rn_ahtrc_0 /=0' )
+         END IF
+      ELSE
+         rldf_rat = rn_ahtrc_0 / rn_aht_0
+      END IF
       !  Define the lateral mixing oparator for tracers
       ! ===============================================
-    
+
       !                               ! control the input
       ioptio = 0
@@ -163,5 +189,5 @@
          ENDIF
          IF ( ln_zps ) THEN             ! z-coordinate
-            IF ( ln_trcldf_level )   ierr = 1      ! iso-level not allowed 
+            IF ( ln_trcldf_level )   ierr = 1      ! iso-level not allowed
             IF ( ln_trcldf_hor   )   nldf = 2      ! horizontal (no rotation)
             IF ( ln_trcldf_iso   )   ierr = 2      ! isoneutral (   rotation)
@@ -206,4 +232,29 @@
       ENDIF
 
+      IF( ln_trcldf_bilap ) THEN
+         IF(lwp) WRITE(numout,*) '          biharmonic tracer diffusion'
+         IF( rn_ahtrc_0 > 0 .AND. .NOT. lk_esopa )   CALL ctl_stop( 'The horizontal diffusivity coef. rn_ahtrc_0 must be negative' )
+      ELSE
+         IF(lwp) WRITE(numout,*) '          harmonic tracer diffusion (default)'
+         IF( rn_ahtrc_0 < 0 .AND. .NOT. lk_esopa )   CALL ctl_stop('The horizontal diffusivity coef. rn_ahtrc_0 must be positive' )
+      ENDIF
+
+      ! ratio between active and passive tracers diffusive coef.
+      IF (ABS(rn_aht_0) < 2._wp*TINY(1.e0)) THEN
+         IF (ABS(rn_ahtrc_0) < 2._wp*TINY(1.e0)) THEN
+            rldf_rat = 1.0_wp
+         ELSE
+            CALL ctl_stop( 'STOP', 'ldf_ctl : cannot define rldf_rat, rn_aht_0==0, rn_ahtrc_0 /=0' )
+         END IF
+      ELSE
+         rldf_rat = rn_ahtrc_0 / rn_aht_0
+      END IF
+      IF( rldf_rat < 0 ) THEN
+         IF( .NOT.lk_offline ) THEN 
+            CALL ctl_stop( 'Choose the same type of diffusive scheme both for active & passive tracers' )
+         ELSE
+            CALL ctl_stop( 'Change the sign of rn_aht_0 in namelist to -/+1' )
+         ENDIF 
+      ENDIF
       !
    END SUBROUTINE ldf_ctl
