Index: /branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/DYN/dynkeg.F90
===================================================================
--- /branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/DYN/dynkeg.F90	(revision 11062)
+++ /branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/DYN/dynkeg.F90	(revision 11063)
@@ -132,4 +132,8 @@
          ENDIF 
       ENDDO  
+!CEOD ADD
+      CALL lbc_lnk( un,'U',-1. )
+      CALL lbc_lnk( vn,'V',-1. )
+!CEOD ADD
 #endif       
 
Index: /branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/TOP_SRC/TRP/trcnam_trp.F90
===================================================================
--- /branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/TOP_SRC/TRP/trcnam_trp.F90	(revision 11062)
+++ /branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/TOP_SRC/TRP/trcnam_trp.F90	(revision 11063)
@@ -16,4 +16,5 @@
    USE in_out_manager      ! ocean dynamics and active tracers variables
    USE lib_mpp           ! distributed memory computing library
+   USE fldread
 
    IMPLICIT NONE
@@ -52,4 +53,8 @@
    INTEGER , PUBLIC ::   nn_zdmp_tr    ! = 0/1/2 flag for damping in the mixed layer
    CHARACTER(LEN=200) , PUBLIC :: cn_resto_tr    !File containing restoration coefficient
+   CHARACTER(LEN=200) , PUBLIC :: cn_dir_chldmp = './'    !: Directory containing chlorophyll file
+   INTEGER , PUBLIC ::    nn_chldmp = 0    !: = 0/1/2 flag for surface chlorophyll damping
+   REAL(wp), PUBLIC ::    rn_chldmp = 0.0  !: chlorophyll damping coefficient
+   TYPE(FLD_N), PUBLIC :: sn_chldmp        !: informations about the fields to be read
 
    !!----------------------------------------------------------------------
@@ -77,5 +82,6 @@
       NAMELIST/namtrc_zdf/ ln_trczdf_exp  , nn_trczdf_exp
       NAMELIST/namtrc_rad/ ln_trcrad
-      NAMELIST/namtrc_dmp/ nn_zdmp_tr , cn_resto_tr
+      NAMELIST/namtrc_dmp/ nn_zdmp_tr , cn_resto_tr, cn_dir_chldmp, nn_chldmp, &
+         &                 sn_chldmp  , rn_chldmp
       !!----------------------------------------------------------------------
 
@@ -179,4 +185,7 @@
          WRITE(numout,*) '      mixed layer damping option     nn_zdmp_tr = ', nn_zdmp_tr, '(zoom: forced to 0)'
          WRITE(numout,*) '      Restoration coeff file    cn_resto_tr = ', cn_resto_tr
+         WRITE(numout,*) '      Surface chlorophyll damping     nn_chldmp = ', nn_chldmp
+         WRITE(numout,*) '      Damping coefficient             rn_chldmp = ', rn_chldmp
+         WRITE(numout,*) '      Chlorophyll directory       cn_dir_chldmp = ', cn_dir_chldmp
       ENDIF
       !
Index: /branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/TOP_SRC/TRP/trcsbcssr.F90
===================================================================
--- /branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/TOP_SRC/TRP/trcsbcssr.F90	(revision 11063)
+++ /branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/TOP_SRC/TRP/trcsbcssr.F90	(revision 11063)
@@ -0,0 +1,379 @@
+MODULE trcsbcssr
+   !!======================================================================
+   !!                       ***  MODULE  trcsbcssr  ***
+   !! Surface module :  restoring term towards surface chlorophyll climatology
+   !!======================================================================
+   !! History :  3.6  !  2017-06  (D. Ford)  Adapt from sbcssr.F90
+   !!                 !  2019-04  (D. Ford)  Add biomass balancing and error checking
+   !!----------------------------------------------------------------------
+#if defined key_top
+   !!----------------------------------------------------------------------
+   !!   trc_sbc_ssr       : add a restoring term toward chl climatology
+   !!   trc_sbc_ssr_init  : initialisation of surface restoring
+   !!----------------------------------------------------------------------
+   USE dom_oce        ! ocean space and time domain
+   USE oce_trc       !  shared variables between ocean and passive tracers
+   USE trc
+   USE trcnam_trp
+   !
+   USE fldread        ! read input fields
+   USE iom            ! I/O manager
+   USE in_out_manager ! I/O manager
+   USE lib_mpp        ! distribued memory computing library
+   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
+   USE timing         ! Timing
+   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)  
+#if defined key_fabm
+   USE par_fabm
+#endif
+
+   IMPLICIT NONE
+   PRIVATE
+
+   PUBLIC   trc_sbc_ssr        ! routine called in trctrp
+
+   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_chldmp   ! structure of input Chl (file informations, fields read)
+
+   !! * Substitutions
+#  include "top_substitute.h90"
+
+CONTAINS
+
+   SUBROUTINE trc_sbc_ssr( kt )
+      !!---------------------------------------------------------------------
+      !!                     ***  ROUTINE trc_sbc_ssr  ***
+      !!
+      !! ** Purpose :   Add to chlorophyll a damping term
+      !!                toward chlorophyll climatology
+      !!
+      !! ** Method  : - Read chlorophyll climatology
+      !!              - at each trc time step add term to each PFT
+      !!                   surface only    (nn_chldmp = 1)
+      !!                   mixed layer     (nn_chldmp = 2)
+      !!              - balance biomass variables
+      !!---------------------------------------------------------------------
+      INTEGER, INTENT(in   ) ::   kt   ! ocean time step
+      !!
+      INTEGER  ::   ji, jj, jk   ! dummy loop indices
+ 
+      REAL(wp), DIMENSION(jpi,jpj)           :: zdiff, zchl
+      REAL(wp), DIMENSION(jpi,jpj,jpk,jptra) :: ztra_save, ztra_new
+      REAL(wp)                               :: zpft, zrat, zinc
+      REAL(wp)                               :: zmax_chl, zmax_c, zmax_n, zmax_p, zmax_s
+      INTEGER                                :: inum
+      !!----------------------------------------------------------------------
+      !
+      IF( nn_timing == 1 )  CALL timing_start('trc_sbc_ssr')
+      !
+      ! Maximum permitted (absolute) updates to phytoplankton variables
+      zmax_chl = 0.01
+      zmax_c   = 0.1
+      zmax_n   = 0.01
+      zmax_p   = 0.001
+      zmax_s   = 0.01
+      
+      IF( kt == nittrc000 )  THEN
+         !
+         CALL trc_sbc_ssr_init
+         !
+         IF( nn_chldmp > 0 ) THEN
+            !
+            IF (lwp) WRITE(numout,*) 'Damping chlorophyll on timestep ', kt
+            
+            ztra_save(:,:,:,:) = tra(:,:,:,:)
+            !
+            CALL fld_read( kt, 1, sf_chldmp )   ! Read Chl data and provides it at kt
+            !
+#if defined key_fabm
+            zchl(:,:) = trb(:,:,1,jp_fabm_m1+jp_fabm_chl1) + &
+               &        trb(:,:,1,jp_fabm_m1+jp_fabm_chl2) + &
+               &        trb(:,:,1,jp_fabm_m1+jp_fabm_chl3) + &
+               &        trb(:,:,1,jp_fabm_m1+jp_fabm_chl4)
+            zdiff(:,:) = rn_chldmp * ( sf_chldmp(1)%fnow(:,:,1) - zchl(:,:) )
+            !
+            DO jj = 2, jpjm1
+               DO ji = fs_2, fs_jpim1   ! vector opt.
+                  IF ( ( sf_chldmp(1)%fnow(ji,jj,1) >   0.0 ) .AND. &
+                     & ( sf_chldmp(1)%fnow(ji,jj,1) < 100.0 ) .AND. &
+                     & ( zchl(ji,jj)                >   0.0 ) ) THEN
+                     ! P1
+                     IF ( ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) > 0.0 ) .AND. &
+                        & ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_p1c)  > 0.0 ) .AND. &
+                        & ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_p1n)  > 0.0 ) .AND. &
+                        & ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_p1p)  > 0.0 ) .AND. &
+                        & ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_p1s)  > 0.0 ) ) THEN
+                        ! Chlorophyll
+                        zpft = ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) / zchl(ji,jj) ) * zdiff(ji,jj)
+                        zpft = MAX( -1.0 * zmax_chl, MIN( zpft, zmax_chl ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_chl1) + zpft
+                        ! Carbon
+                        zrat = trb(ji,jj,1,jp_fabm_m1+jp_fabm_p1c) / trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl1)
+                        zinc = MAX( -1.0 * zmax_c, MIN( (zpft * zrat), zmax_c ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_p1c) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_p1c) + zinc
+                        ! Nitrogen
+                        zrat = trb(ji,jj,1,jp_fabm_m1+jp_fabm_p1n) / trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl1)
+                        zinc = MAX( -1.0 * zmax_n, MIN( (zpft * zrat), zmax_n ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_p1n) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_p1n) + zinc
+                        ! Phosphorus
+                        zrat = trb(ji,jj,1,jp_fabm_m1+jp_fabm_p1p) / trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl1)
+                        zinc = MAX( -1.0 * zmax_p, MIN( (zpft * zrat), zmax_p ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_p1p) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_p1p) + zinc
+                        ! Silicon
+                        zrat = trb(ji,jj,1,jp_fabm_m1+jp_fabm_p1s) / trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl1)
+                        zinc = MAX( -1.0 * zmax_s, MIN( (zpft * zrat), zmax_s ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_p1s) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_p1s) + zinc
+                        !
+                        IF( nn_chldmp == 2 ) THEN
+                           ! Mixed layer
+                           DO jk = 2, jpkm1
+                              IF( fsdept(ji,jj,jk) < hmlp (ji,jj) ) THEN
+                                 IF ( ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl1) > 0.0 ) .AND. &
+                                    & ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p1c)  > 0.0 ) .AND. &
+                                    & ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p1n)  > 0.0 ) .AND. &
+                                    & ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p1p)  > 0.0 ) .AND. &
+                                    & ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p1s)  > 0.0 ) ) THEN
+                                    ! Chlorophyll
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_chl1) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_chl1) + zpft
+                                    ! Carbon
+                                    zrat = trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p1c) / trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl1)
+                                    zinc = MAX( -1.0 * zmax_c, MIN( (zpft * zrat), zmax_c ) )
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p1c) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p1c) + zinc
+                                    ! Nitrogen
+                                    zrat = trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p1n) / trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl1)
+                                    zinc = MAX( -1.0 * zmax_n, MIN( (zpft * zrat), zmax_n ) )
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p1n) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p1n) + zinc
+                                    ! Phosphorus
+                                    zrat = trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p1p) / trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl1)
+                                    zinc = MAX( -1.0 * zmax_p, MIN( (zpft * zrat), zmax_p ) )
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p1p) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p1p) + zinc
+                                    ! Silicon
+                                    zrat = trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p1s) / trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl1)
+                                    zinc = MAX( -1.0 * zmax_s, MIN( (zpft * zrat), zmax_s ) )
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p1s) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p1s) + zinc
+                                 ENDIF
+                              ENDIF
+                           END DO
+                        ENDIF
+                     ENDIF
+                     ! P2
+                     IF ( ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) > 0.0 ) .AND. &
+                        & ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_p2c)  > 0.0 ) .AND. &
+                        & ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_p2n)  > 0.0 ) .AND. &
+                        & ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_p2p)  > 0.0 ) ) THEN
+                        ! Chlorophyll
+                        zpft = ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) / zchl(ji,jj) ) * zdiff(ji,jj)
+                        zpft = MAX( -1.0 * zmax_chl, MIN( zpft, zmax_chl ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_chl2) + zpft
+                        ! Carbon
+                        zrat = trb(ji,jj,1,jp_fabm_m1+jp_fabm_p2c) / trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl2)
+                        zinc = MAX( -1.0 * zmax_c, MIN( (zpft * zrat), zmax_c ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_p2c) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_p2c) + zinc
+                        ! Nitrogen
+                        zrat = trb(ji,jj,1,jp_fabm_m1+jp_fabm_p2n) / trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl2)
+                        zinc = MAX( -1.0 * zmax_n, MIN( (zpft * zrat), zmax_n ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_p2n) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_p2n) + zinc
+                        ! Phosphorus
+                        zrat = trb(ji,jj,1,jp_fabm_m1+jp_fabm_p2p) / trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl2)
+                        zinc = MAX( -1.0 * zmax_p, MIN( (zpft * zrat), zmax_p ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_p2p) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_p2p) + zinc
+                        !
+                        IF( nn_chldmp == 2 ) THEN
+                           ! Mixed layer
+                           DO jk = 2, jpkm1
+                              IF( fsdept(ji,jj,jk) < hmlp (ji,jj) ) THEN
+                                 IF ( ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl2) > 0.0 ) .AND. &
+                                    & ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p2c)  > 0.0 ) .AND. &
+                                    & ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p2n)  > 0.0 ) .AND. &
+                                    & ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p2p)  > 0.0 ) ) THEN
+                                    ! Chlorophyll
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_chl2) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_chl2) + zpft
+                                    ! Carbon
+                                    zrat = trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p2c) / trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl2)
+                                    zinc = MAX( -1.0 * zmax_c, MIN( (zpft * zrat), zmax_c ) )
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p2c) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p2c) + zinc
+                                    ! Nitrogen
+                                    zrat = trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p2n) / trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl2)
+                                    zinc = MAX( -1.0 * zmax_n, MIN( (zpft * zrat), zmax_n ) )
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p2n) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p2n) + zinc
+                                    ! Phosphorus
+                                    zrat = trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p2p) / trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl2)
+                                    zinc = MAX( -1.0 * zmax_p, MIN( (zpft * zrat), zmax_p ) )
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p2p) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p2p) + zinc
+                                 ENDIF
+                              ENDIF
+                           END DO
+                        ENDIF
+                     ENDIF
+                     ! P3
+                     IF ( ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) > 0.0 ) .AND. &
+                        & ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_p3c)  > 0.0 ) .AND. &
+                        & ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_p3n)  > 0.0 ) .AND. &
+                        & ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_p3p)  > 0.0 ) ) THEN
+                        ! Chlorophyll
+                        zpft = ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) / zchl(ji,jj) ) * zdiff(ji,jj)
+                        zpft = MAX( -1.0 * zmax_chl, MIN( zpft, zmax_chl ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_chl3) + zpft
+                        ! Carbon
+                        zrat = trb(ji,jj,1,jp_fabm_m1+jp_fabm_p3c) / trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl3)
+                        zinc = MAX( -1.0 * zmax_c, MIN( (zpft * zrat), zmax_c ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_p3c) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_p3c) + zinc
+                        ! Nitrogen
+                        zrat = trb(ji,jj,1,jp_fabm_m1+jp_fabm_p3n) / trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl3)
+                        zinc = MAX( -1.0 * zmax_n, MIN( (zpft * zrat), zmax_n ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_p3n) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_p3n) + zinc
+                        ! Phosphorus
+                        zrat = trb(ji,jj,1,jp_fabm_m1+jp_fabm_p3p) / trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl3)
+                        zinc = MAX( -1.0 * zmax_p, MIN( (zpft * zrat), zmax_p ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_p3p) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_p3p) + zinc
+                        !
+                        IF( nn_chldmp == 2 ) THEN
+                           ! Mixed layer
+                           DO jk = 2, jpkm1
+                              IF( fsdept(ji,jj,jk) < hmlp (ji,jj) ) THEN
+                                 IF ( ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl3) > 0.0 ) .AND. &
+                                    & ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p3c)  > 0.0 ) .AND. &
+                                    & ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p3n)  > 0.0 ) .AND. &
+                                    & ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p3p)  > 0.0 ) ) THEN
+                                    ! Chlorophyll
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_chl3) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_chl3) + zpft
+                                    ! Carbon
+                                    zrat = trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p3c) / trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl3)
+                                    zinc = MAX( -1.0 * zmax_c, MIN( (zpft * zrat), zmax_c ) )
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p3c) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p3c) + zinc
+                                    ! Nitrogen
+                                    zrat = trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p3n) / trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl3)
+                                    zinc = MAX( -1.0 * zmax_n, MIN( (zpft * zrat), zmax_n ) )
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p3n) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p3n) + zinc
+                                    ! Phosphorus
+                                    zrat = trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p3p) / trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl3)
+                                    zinc = MAX( -1.0 * zmax_p, MIN( (zpft * zrat), zmax_p ) )
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p3p) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p3p) + zinc
+                                 ENDIF
+                              ENDIF
+                           END DO
+                        ENDIF
+                     ENDIF
+                     ! P4
+                     IF ( ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) > 0.0 ) .AND. &
+                        & ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_p4c)  > 0.0 ) .AND. &
+                        & ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_p4n)  > 0.0 ) .AND. &
+                        & ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_p4p)  > 0.0 ) ) THEN
+                        ! Chlorophyll
+                        zpft = ( trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) / zchl(ji,jj) ) * zdiff(ji,jj)
+                        zpft = MAX( -1.0 * zmax_chl, MIN( zpft, zmax_chl ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_chl4) + zpft
+                        ! Carbon
+                        zrat = trb(ji,jj,1,jp_fabm_m1+jp_fabm_p4c) / trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl4)
+                        zinc = MAX( -1.0 * zmax_c, MIN( (zpft * zrat), zmax_c ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_p4c) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_p4c) + zinc
+                        ! Nitrogen
+                        zrat = trb(ji,jj,1,jp_fabm_m1+jp_fabm_p4n) / trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl4)
+                        zinc = MAX( -1.0 * zmax_n, MIN( (zpft * zrat), zmax_n ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_p4n) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_p4n) + zinc
+                        ! Phosphorus
+                        zrat = trb(ji,jj,1,jp_fabm_m1+jp_fabm_p4p) / trb(ji,jj,1,jp_fabm_m1+jp_fabm_chl4)
+                        zinc = MAX( -1.0 * zmax_p, MIN( (zpft * zrat), zmax_p ) )
+                        tra(ji,jj,1,jp_fabm_m1+jp_fabm_p4p) = tra(ji,jj,1,jp_fabm_m1+jp_fabm_p4p) + zinc
+                        !
+                        IF( nn_chldmp == 2 ) THEN
+                           ! Mixed layer
+                           DO jk = 2, jpkm1
+                              IF( fsdept(ji,jj,jk) < hmlp (ji,jj) ) THEN
+                                 IF ( ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl4) > 0.0 ) .AND. &
+                                    & ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p4c)  > 0.0 ) .AND. &
+                                    & ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p4n)  > 0.0 ) .AND. &
+                                    & ( trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p4p)  > 0.0 ) ) THEN
+                                    ! Chlorophyll
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_chl4) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_chl4) + zpft
+                                    ! Carbon
+                                    zrat = trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p4c) / trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl4)
+                                    zinc = MAX( -1.0 * zmax_c, MIN( (zpft * zrat), zmax_c ) )
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p4c) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p4c) + zinc
+                                    ! Nitrogen
+                                    zrat = trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p4n) / trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl4)
+                                    zinc = MAX( -1.0 * zmax_n, MIN( (zpft * zrat), zmax_n ) )
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p4n) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p4n) + zinc
+                                    ! Phosphorus
+                                    zrat = trb(ji,jj,jk,jp_fabm_m1+jp_fabm_p4p) / trb(ji,jj,jk,jp_fabm_m1+jp_fabm_chl4)
+                                    zinc = MAX( -1.0 * zmax_p, MIN( (zpft * zrat), zmax_p ) )
+                                    tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p4p) = tra(ji,jj,jk,jp_fabm_m1+jp_fabm_p4p) + zinc
+                                 ENDIF
+                              ENDIF
+                           END DO
+                        ENDIF
+                     ENDIF
+                  ENDIF
+               END DO
+            END DO
+            
+            ! Save updates
+            ztra_new(:,:,:,:) = tra(:,:,:,:) - ztra_save(:,:,:,:)
+            CALL iom_open( 'chldmp.balincs', inum, ldwrt = .TRUE., kiolib = jprstlib)
+            zdiff(:,:) = zdiff(:,:) * tmask(:,:,1)
+            CALL iom_rstput( kt, kt, inum, 'ChlTot', zdiff(:,:) )
+            CALL iom_rstput( kt, kt, inum, 'P1_Chl', ztra_new(:,:,:,jp_fabm_chl1) )
+            CALL iom_rstput( kt, kt, inum, 'P1_c',   ztra_new(:,:,:,jp_fabm_p1c)  )
+            CALL iom_rstput( kt, kt, inum, 'P1_n',   ztra_new(:,:,:,jp_fabm_p1n)  )
+            CALL iom_rstput( kt, kt, inum, 'P1_p',   ztra_new(:,:,:,jp_fabm_p1p)  )
+            CALL iom_rstput( kt, kt, inum, 'P1_s',   ztra_new(:,:,:,jp_fabm_p1s)  )
+            CALL iom_rstput( kt, kt, inum, 'P2_Chl', ztra_new(:,:,:,jp_fabm_chl2) )
+            CALL iom_rstput( kt, kt, inum, 'P2_c',   ztra_new(:,:,:,jp_fabm_p2c)  )
+            CALL iom_rstput( kt, kt, inum, 'P2_n',   ztra_new(:,:,:,jp_fabm_p2n)  )
+            CALL iom_rstput( kt, kt, inum, 'P2_p',   ztra_new(:,:,:,jp_fabm_p2p)  )
+            CALL iom_rstput( kt, kt, inum, 'P3_Chl', ztra_new(:,:,:,jp_fabm_chl3) )
+            CALL iom_rstput( kt, kt, inum, 'P3_c',   ztra_new(:,:,:,jp_fabm_p3c)  )
+            CALL iom_rstput( kt, kt, inum, 'P3_n',   ztra_new(:,:,:,jp_fabm_p3n)  )
+            CALL iom_rstput( kt, kt, inum, 'P3_p',   ztra_new(:,:,:,jp_fabm_p3p)  )
+            CALL iom_rstput( kt, kt, inum, 'P4_Chl', ztra_new(:,:,:,jp_fabm_chl4) )
+            CALL iom_rstput( kt, kt, inum, 'P4_c',   ztra_new(:,:,:,jp_fabm_p4c)  )
+            CALL iom_rstput( kt, kt, inum, 'P4_n',   ztra_new(:,:,:,jp_fabm_p4n)  )
+            CALL iom_rstput( kt, kt, inum, 'P4_p',   ztra_new(:,:,:,jp_fabm_p4p)  )
+            CALL iom_close( inum )
+#else
+            CALL ctl_stop( 'STOP', 'trc_sbc_ssr: only works with FABM-ERSEM' )
+#endif
+            !
+         ENDIF
+         !
+      ENDIF
+      !
+      IF( nn_timing == 1 )  CALL timing_stop('trc_sbc_ssr')
+      !
+   END SUBROUTINE trc_sbc_ssr
+
+ 
+   SUBROUTINE trc_sbc_ssr_init
+      !!---------------------------------------------------------------------
+      !!                  ***  ROUTINE trc_sbc_ssr_init  ***
+      !!
+      !! ** Purpose :   initialisation of surface damping term
+      !!
+      !! ** Method  : - Read chlorophyll
+      !!---------------------------------------------------------------------
+      INTEGER  ::   ierror   ! return error code
+      !!----------------------------------------------------------------------
+      !
+      IF( nn_chldmp > 0 ) THEN      !* set sf_sss structure & allocate arrays
+         !
+         ALLOCATE( sf_chldmp(1), STAT=ierror )
+         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'trc_sbc_ssr: unable to allocate sf_chldmp structure' )
+         ALLOCATE( sf_chldmp(1)%fnow(jpi,jpj,1), STAT=ierror )
+         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'trc_sbc_ssr: unable to allocate sf_chldmp now array' )
+         !
+         ! fill sf_sss with sn_sss and control print
+         CALL fld_fill( sf_chldmp, (/ sn_chldmp /), cn_dir_chldmp, 'trc_sbc_ssr', 'Chl restoring term', 'namtrc_dmp' )
+         IF( sf_chldmp(1)%ln_tint )   ALLOCATE( sf_chldmp(1)%fdta(jpi,jpj,1,2), STAT=ierror )
+         IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'trc_sbc_ssr: unable to allocate sf_chldmp data array' )
+         !
+      ENDIF
+      !
+   END SUBROUTINE trc_sbc_ssr_init
+
+#else
+   SUBROUTINE trc_sbc_ssr( kt )        ! Empty routine
+      INTEGER, INTENT(in) :: kt
+      WRITE(*,*) 'trc_sbc_ssr: You should not have seen this print! error?', kt
+   END SUBROUTINE trc_sbc_ssr
+#endif
+   !!======================================================================
+END MODULE trcsbcssr
Index: /branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/TOP_SRC/TRP/trctrp.F90
===================================================================
--- /branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/TOP_SRC/TRP/trctrp.F90	(revision 11062)
+++ /branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/TOP_SRC/TRP/trctrp.F90	(revision 11063)
@@ -29,4 +29,5 @@
    USE trcbdy          ! BDY open boundaries
    USE bdy_par, only: lk_bdy
+   USE trcsbcssr
 
 #if defined key_agrif
@@ -66,4 +67,5 @@
       IF( .NOT. lk_c1d ) THEN
          !
+         IF( nn_chldmp > 0 )    CALL trc_sbc_ssr( kstp )        ! add Chl damping term
                                 CALL trc_sbc( kstp )            ! surface boundary condition
          IF( lk_trabbl )        CALL trc_bbl( kstp )            ! advective (and/or diffusive) bottom boundary layer scheme
@@ -93,4 +95,5 @@
          !
       ELSE                                               ! 1D vertical configuration
+         IF( nn_chldmp > 0 )    CALL trc_sbc_ssr( kstp )        ! add Chl damping term
                                 CALL trc_sbc( kstp )            ! surface boundary condition
          IF( .NOT. lk_offline .AND. lk_zdfkpp )    &
