Index: branches/2013/dev_r3856_MERCATOR3_QSRMEAN24H/NEMOGCM/NEMO/TOP_SRC/TRP/trcdmp.F90
===================================================================
--- branches/2013/dev_r3856_MERCATOR3_QSRMEAN24H/NEMOGCM/NEMO/TOP_SRC/TRP/trcdmp.F90	(revision 3857)
+++ branches/2013/dev_r3856_MERCATOR3_QSRMEAN24H/NEMOGCM/NEMO/TOP_SRC/TRP/trcdmp.F90	(revision 4276)
@@ -11,7 +11,5 @@
    !!            3.3  !  2010-06  (C. Ethe, G. Madec) merge TRA-TRC 
    !!----------------------------------------------------------------------
-#if  defined key_top && defined key_trcdmp 
-   !!----------------------------------------------------------------------
-   !!   key_trcdmp                                         internal damping
+#if  defined key_top 
    !!----------------------------------------------------------------------
    !!   trc_dmp      : update the tracer trend with the internal damping
@@ -25,4 +23,5 @@
    USE prtctl_trc      ! Print control for debbuging
    USE trdtra
+   USE trdmod_oce
 
    IMPLICIT NONE
@@ -30,17 +29,13 @@
 
    PUBLIC trc_dmp            ! routine called by step.F90
+   PUBLIC trc_dmp_clo        ! routine called by step.F90
    PUBLIC trc_dmp_alloc      ! routine called by nemogcm.F90
 
-   LOGICAL , PUBLIC, PARAMETER ::   lk_trcdmp = .TRUE.   !: internal damping flag
-
    !                                !!* Namelist namtrc_dmp : passive tracer newtonian damping *
-   INTEGER  ::   nn_hdmp_tr =   -1   ! = 0/-1/'latitude' for damping over passive tracer
-   INTEGER  ::   nn_zdmp_tr =    0   ! = 0/1/2 flag for damping in the mixed layer
-   REAL(wp) ::   rn_surf_tr =   50.  ! surface time scale for internal damping        [days]
-   REAL(wp) ::   rn_bot_tr  =  360.  ! bottom time scale for internal damping         [days]
-   REAL(wp) ::   rn_dep_tr  =  800.  ! depth of transition between rn_surf and rn_bot [meters]
-   INTEGER  ::   nn_file_tr =    2   ! = 1 create a damping.coeff NetCDF file 
-
    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   restotr   ! restoring coeff. on tracers (s-1)
+
+   INTEGER, PARAMETER           ::   npncts   = 5        ! number of closed sea
+   INTEGER, DIMENSION(npncts)   ::   nctsi1, nctsj1      ! south-west closed sea limits (i,j)
+   INTEGER, DIMENSION(npncts)   ::   nctsi2, nctsj2      ! north-east closed sea limits (i,j)
 
    !! * Substitutions
@@ -86,8 +81,9 @@
       INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
       !!
-      INTEGER  ::   ji, jj, jk, jn       ! dummy loop indices
+      INTEGER  ::   ji, jj, jk, jn, jl       ! dummy loop indices
       REAL(wp) ::   ztra                 ! temporary scalars
       CHARACTER (len=22) :: charout
       REAL(wp), POINTER, DIMENSION(:,:,:) ::   ztrtrd
+      REAL(wp), POINTER, DIMENSION(:,:,:,:) ::  ztrcdta   ! 4D  workspace
       !!----------------------------------------------------------------------
       !
@@ -99,63 +95,68 @@
 
       IF( l_trdtrc )   CALL wrk_alloc( jpi, jpj, jpk, ztrtrd )   ! temporary save of trends
-
-      ! 1. Newtonian damping trends on tracer fields
-      ! --------------------------------------------
-      ! Initialize the input fields for newtonian damping
-      CALL trc_dta( kt )
-      !                                                          ! ===========
-      DO jn = 1, jptra                                           ! tracer loop
-         !                                                       ! ===========
-         IF( l_trdtrc ) ztrtrd(:,:,:) = tra(:,:,:,jn)    ! save trends 
-
-         IF( lutini(jn) ) THEN
+      !
+      IF( nb_trcdta > 0 ) THEN  ! Initialisation of tracer from a file that may also be used for damping
+         !
+         CALL wrk_alloc( jpi, jpj, jpk, nb_trcdta, ztrcdta )    ! Memory allocation
+         CALL trc_dta( kt, ztrcdta )   ! read tracer data at nit000
+         !                                                          ! ===========
+         DO jn = 1, jptra                                           ! tracer loop
+            !                                                       ! ===========
+            IF( l_trdtrc ) ztrtrd(:,:,:) = tra(:,:,:,jn)    ! save trends 
             !
-            SELECT CASE ( nn_zdmp_trc )
-            !
-            CASE( 0 )                !==  newtonian damping throughout the water column  ==!
-               DO jk = 1, jpkm1
-                  DO jj = 2, jpjm1
-                     DO ji = fs_2, fs_jpim1   ! vector opt.
-                        ztra = restotr(ji,jj,jk) * ( trdta(ji,jj,jk,jn) - trb(ji,jj,jk,jn) )
-                        tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + ztra
+            IF( ln_trc_ini(jn) ) THEN      ! update passive tracers arrays with input data read from file
+               
+               jl = n_trc_index(jn) 
+
+               SELECT CASE ( nn_zdmp_tr )
+               !
+               CASE( 0 )                !==  newtonian damping throughout the water column  ==!
+                  DO jk = 1, jpkm1
+                     DO jj = 2, jpjm1
+                        DO ji = fs_2, fs_jpim1   ! vector opt.
+                           ztra = restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk,jl) - trb(ji,jj,jk,jn) )
+                           tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + ztra
+                        END DO
                      END DO
                   END DO
-               END DO
-            !
-            CASE ( 1 )                !==  no damping in the turbocline (avt > 5 cm2/s)  ==!
-               DO jk = 1, jpkm1
-                  DO jj = 2, jpjm1
-                     DO ji = fs_2, fs_jpim1   ! vector opt.
-                        IF( avt(ji,jj,jk) <= 5.e-4 )  THEN 
-                           ztra = restotr(ji,jj,jk) * ( trdta(ji,jj,jk,jn) - trb(ji,jj,jk,jn) )
-                           tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + ztra
-                        ENDIF
+               !
+               CASE ( 1 )                !==  no damping in the turbocline (avt > 5 cm2/s)  ==!
+                  DO jk = 1, jpkm1
+                     DO jj = 2, jpjm1
+                        DO ji = fs_2, fs_jpim1   ! vector opt.
+                           IF( avt(ji,jj,jk) <= 5.e-4 )  THEN 
+                              ztra = restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk,jl) - trb(ji,jj,jk,jn) )
+                              tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + ztra
+                           ENDIF
+                        END DO
                      END DO
                   END DO
-               END DO
-            !
-            CASE ( 2 )               !==  no damping in the mixed layer   ==! 
-               DO jk = 1, jpkm1
-                  DO jj = 2, jpjm1
-                     DO ji = fs_2, fs_jpim1   ! vector opt.
-                        IF( fsdept(ji,jj,jk) >= hmlp (ji,jj) ) THEN
-                           ztra = restotr(ji,jj,jk,jn) * ( trdta(ji,jj,jk,jn) - trb(ji,jj,jk,jn) )
-                           tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + ztra
-                        END IF
+               !
+               CASE ( 2 )               !==  no damping in the mixed layer   ==! 
+                  DO jk = 1, jpkm1
+                     DO jj = 2, jpjm1
+                        DO ji = fs_2, fs_jpim1   ! vector opt.
+                           IF( fsdept(ji,jj,jk) >= hmlp (ji,jj) ) THEN
+                              ztra = restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk,jl) - trb(ji,jj,jk,jn) )
+                              tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + ztra
+                           END IF
+                        END DO
                      END DO
                   END DO
-               END DO
-            !  
-            END SELECT
-            ! 
-         ENDIF
-         !
-         IF( l_trdtrc ) THEN
-            ztrtrd(:,:,:) = tra(:,:,:,jn) -  ztrtrd(:,:,:)
-            CALL trd_tra( kt, 'TRC', jn, jptra_trd_dmp, ztrtrd )
-         END IF
-         !                                                       ! ===========
-      END DO                                                     ! tracer loop
-      !                                                          ! ===========
+               !  
+               END SELECT
+               ! 
+            ENDIF
+            !
+            IF( l_trdtrc ) THEN
+               ztrtrd(:,:,:) = tra(:,:,:,jn) -  ztrtrd(:,:,:)
+               CALL trd_tra( kt, 'TRC', jn, jptra_trd_dmp, ztrtrd )
+            END IF
+            !                                                       ! ===========
+         END DO                                                     ! tracer loop
+         !                                                          ! ===========
+         CALL wrk_dealloc( jpi, jpj, jpk, nb_trcdta, ztrcdta )
+      ENDIF
+      !
       IF( l_trdtrc )  CALL wrk_dealloc( jpi, jpj, jpk, ztrtrd )
       !                                          ! print mean trends (used for debugging)
@@ -168,4 +169,125 @@
       !
    END SUBROUTINE trc_dmp
+
+   SUBROUTINE trc_dmp_clo( kt )
+      !!---------------------------------------------------------------------
+      !!                  ***  ROUTINE trc_dmp_clo  ***
+      !!
+      !! ** Purpose :   Closed sea domain initialization
+      !!
+      !! ** Method  :   if a closed sea is located only in a model grid point
+      !!                we restore to initial data
+      !!
+      !! ** Action  :   nctsi1(), nctsj1() : south-west closed sea limits (i,j)
+      !!                nctsi2(), nctsj2() : north-east Closed sea limits (i,j)
+      !!----------------------------------------------------------------------
+      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
+      !
+      INTEGER :: ji, jj, jk, jn, jl, jc                     ! dummy loop indicesa
+      REAL(wp), POINTER, DIMENSION(:,:,:,:) ::  ztrcdta     ! 4D  workspace
+
+      !!----------------------------------------------------------------------
+
+      IF( kt == nit000 ) THEN
+         ! initial values
+         nctsi1(:) = 1  ;  nctsi2(:) = 1
+         nctsj1(:) = 1  ;  nctsj2(:) = 1
+
+         ! set the closed seas (in data domain indices)
+         ! -------------------
+
+         IF( cp_cfg == "orca" ) THEN
+            !
+            SELECT CASE ( jp_cfg )
+            !                                           ! =======================
+            CASE ( 2 )                                  !  ORCA_R2 configuration
+               !                                        ! =======================
+               !                                            ! Caspian Sea
+               nctsi1(1)   =  11  ;  nctsj1(1)   = 103
+               nctsi2(1)   =  17  ;  nctsj2(1)   = 112
+               !                                            ! Great North American Lakes
+               nctsi1(2)   =  97  ;  nctsj1(2)   = 107
+               nctsi2(2)   = 103  ;  nctsj2(2)   = 111
+               !                                            ! Black Sea 1 : west part of the Black Sea
+               nctsi1(3)   = 174  ;  nctsj1(3)   = 107
+               nctsi2(3)   = 181  ;  nctsj2(3)   = 112
+              !                                            ! Black Sea 2 : est part of the Black Sea
+               nctsi1(4)   =   2  ;  nctsj1(4)   = 107
+               nctsi2(4)   =   6  ;  nctsj2(4)   = 112
+               !                                            ! Baltic Sea
+               nctsi1(5)   =  145 ;  nctsj1(5)   = 116
+               nctsi2(5)   =  150 ;  nctsj2(5)   = 126
+               !                                        ! =======================
+            CASE ( 4 )                                  !  ORCA_R4 configuration
+               !                                        ! =======================
+               !                                            ! Caspian Sea
+               nctsi1(1)   =  4  ;  nctsj1(1)   = 53
+               nctsi2(1)   =  4  ;  nctsj2(1)   = 56
+               !                                            ! Great North American Lakes
+               nctsi1(2)   = 49  ;  nctsj1(2)   = 55
+               nctsi2(2)   = 51  ;  nctsj2(2)   = 56
+               !                                            ! Black Sea
+               nctsi1(3)   = 88  ;  nctsj1(3)   = 55
+               nctsi2(3)   = 91  ;  nctsj2(3)   = 56
+               !                                            ! Baltic Sea
+               nctsi1(4)   = 75  ;  nctsj1(4)   = 59
+               nctsi2(4)   = 76  ;  nctsj2(4)   = 61
+               !                                        ! =======================
+            CASE ( 025 )                                ! ORCA_R025 configuration
+               !                                        ! =======================
+                                                     ! Caspian + Aral sea
+               nctsi1(1)   = 1330 ; nctsj1(1)   = 645
+               nctsi2(1)   = 1400 ; nctsj2(1)   = 795
+               !                                        ! Azov Sea
+               nctsi1(2)   = 1284 ; nctsj1(2)   = 722
+               nctsi2(2)   = 1304 ; nctsj2(2)   = 747
+               !
+            END SELECT
+            !
+         ENDIF
+         !
+
+         ! convert the position in local domain indices
+         ! --------------------------------------------
+         DO jc = 1, npncts
+            nctsi1(jc)   = mi0( nctsi1(jc) )
+            nctsj1(jc)   = mj0( nctsj1(jc) )
+
+            nctsi2(jc)   = mi1( nctsi2(jc) )
+            nctsj2(jc)   = mj1( nctsj2(jc) )
+         END DO
+         !
+      ENDIF
+
+      ! Restore close seas values to initial data
+      IF( ln_trcdta .AND. nb_trcdta > 0 )  THEN   ! Initialisation of tracer from a file that may also be used for damping
+         !
+         IF(lwp)  WRITE(numout,*)
+         IF(lwp)  WRITE(numout,*) ' trc_dmp_clo : Restoring of nutrients on close seas at time-step kt = ', kt
+         IF(lwp)  WRITE(numout,*)
+         !
+         CALL wrk_alloc( jpi, jpj, jpk, nb_trcdta, ztrcdta )   ! Memory allocation
+         !
+         CALL trc_dta( kt , ztrcdta )   ! read tracer data at nittrc000
+         !
+         DO jn = 1, jptra
+            IF( ln_trc_ini(jn) ) THEN      ! update passive tracers arrays with input data read from file
+                jl = n_trc_index(jn)
+                DO jc = 1, npncts
+                   DO jk = 1, jpkm1
+                      DO jj = nctsj1(jc), nctsj2(jc)
+                         DO ji = nctsi1(jc), nctsi2(jc)
+                            trn(ji,jj,jk,jn) = ztrcdta(ji,jj,jk,jl) * tmask(ji,jj,jk)
+                            trb(ji,jj,jk,jn) = trn(ji,jj,jk,jn)
+                         ENDDO
+                      ENDDO
+                   ENDDO
+                ENDDO
+             ENDIF
+          ENDDO
+          CALL wrk_dealloc( jpi, jpj, jpk, nb_trcdta, ztrcdta )
+      ENDIF
+      !
+   END SUBROUTINE trc_dmp_clo
 
 
@@ -199,8 +321,5 @@
       END SELECT
 
-      IF( .NOT. lk_dtatrc )   &
-         &   CALL ctl_stop( 'no passive tracer data define key_dtatrc' )
-
-      IF( .NOT. lk_tradmp )   &
+      IF( .NOT. ln_tradmp )   &
          &   CALL ctl_stop( 'passive trace damping need key_tradmp to compute damping coef.' )
       !
@@ -214,9 +333,9 @@
       !
    END SUBROUTINE trc_dmp_init
+
 #else
    !!----------------------------------------------------------------------
-   !!   Default key                                     NO internal damping
-   !!----------------------------------------------------------------------
-   LOGICAL , PUBLIC, PARAMETER ::   lk_trcdmp = .FALSE.    !: internal damping flag
+   !!  Dummy module :                                     No passive tracer
+   !!----------------------------------------------------------------------
 CONTAINS
    SUBROUTINE trc_dmp( kt )        ! Empty routine
@@ -225,4 +344,6 @@
    END SUBROUTINE trc_dmp
 #endif
+
+
    !!======================================================================
 END MODULE trcdmp
