Index: trunk/NEMO/LIM_SRC_3/dom_ice.F90
===================================================================
--- trunk/NEMO/LIM_SRC_3/dom_ice.F90	(revision 1607)
+++ trunk/NEMO/LIM_SRC_3/dom_ice.F90	(revision 1608)
@@ -2,15 +2,8 @@
    !!======================================================================
    !!                   ***  MODULE  dom_ice  ***
-   !! LIM Sea Ice :   Domain  variables
+   !! LIM-3 Sea Ice :   Domain  variables
    !!======================================================================
-   !! History :
-   !!   2.0  !  03-08  (C. Ethe)  Free form and module
-   !!   3.0  !  03-08  (M. Vancoppenolle), Ice cats, salinity and EVP-C
+   !! History :  3.0  ! 2003-08  (M. Vancoppenolle)  LIM-3
    !!----------------------------------------------------------------------
-   !!   LIM 3.0, UCL-ASTR-LOCEAN-IPSL (2005)
-   !! $Id$
-   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
-   !!----------------------------------------------------------------------
-   !! * Modules used
    USE par_ice
 
@@ -18,29 +11,22 @@
    PRIVATE
 
-   !! * Share module variables
-   LOGICAL, PUBLIC ::   &  !:
-      l_jeq = .TRUE.       !: Equator inside the domain flag
+   LOGICAL, PUBLIC ::   l_jeq = .TRUE.       !: Equator inside the domain flag
 
-   INTEGER, PUBLIC ::   &  !:
-      njeq , njeqm1        !: j-index of the equator if it is inside the domain
-   !                    !  (otherwise = jpj+10 (SH) or -10 (SH) )
+   INTEGER, PUBLIC ::   njeq , njeqm1        !: j-index of the equator if it is inside the domain
 
-   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   &  !:
-      fs2cor ,          &  !: coriolis factor
-      fcor   ,          &  !: coriolis coefficient
-      covrai ,          &  !: sine of geographic latitude
-      area   ,          &  !: surface of grid cell 
-      tms , tmi ,       &  !: temperature mask, mask for stress
-      tmu , tmv ,       &  !: mask at u-v velocity points
-      tmf
+   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fs2cor     !: coriolis factor
+   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fcor       !: coriolis coefficient
+   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   covrai     !: sine of geographic latitude
+   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   area       !: surface of grid cell 
+   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   tms, tmi   !: temperature mask, mask for stress
+   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   tmu, tmv   !: mask at u and v velocity points
+   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   tmf        !: mask at f-point
 
-   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,2,2) ::   &  !:
-      wght   ,          &  !: weight of the 4 neighbours to compute averages
-      akappa ,          &  !: first group of metric coefficients
-      bkappa               !: third group of metric coefficients
+   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,2,2) ::   wght     !: weight of the 4 neighbours to compute averages
 
-   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,2,2,2,2) ::   &  !:
-      alambd               !: second group of metric coefficients
-
+   !!----------------------------------------------------------------------
+   !! NEMO/LIM 3.2, UCL-ASTR-LOCEAN-IPSL (2009)
+   !! $Id$
+   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
    !!======================================================================
 END MODULE dom_ice
Index: trunk/NEMO/LIM_SRC_3/limmsh.F90
===================================================================
--- trunk/NEMO/LIM_SRC_3/limmsh.F90	(revision 1607)
+++ trunk/NEMO/LIM_SRC_3/limmsh.F90	(revision 1608)
@@ -4,4 +4,6 @@
    !! LIM ice model :   definition of the ice mesh parameters
    !!======================================================================
+   !! History :  3.2  !  2008-01 (NEMO team)  LIM-3: adaptation from LIM-2
+   !!----------------------------------------------------------------------
 #if defined key_lim3
    !!----------------------------------------------------------------------
@@ -10,21 +12,19 @@
    !!   lim_msh   : definition of the ice mesh
    !!----------------------------------------------------------------------
-   !! * Modules used
-   USE phycst
-   USE dom_oce
-   USE dom_ice
-   USE lbclnk
-   USE in_out_manager
+   USE phycst         ! physical constants
+   USE dom_oce        ! ocean domain
+   USE dom_ice        ! sea-ice domain
+   USE in_out_manager ! I/O manager
+   USE lbclnk         ! 
 
    IMPLICIT NONE
    PRIVATE
 
-   !! * Accessibility
-   PUBLIC lim_msh      ! routine called by ice_ini.F90
+   PUBLIC   lim_msh   ! routine called by ice_ini.F90
 
    !!----------------------------------------------------------------------
-   !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008) 
+   !! NEMO/LIM 3.2,  UCL-ASTR-LOCEAN-IPSL (2009) 
    !! $Id$
-   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
+   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
    !!----------------------------------------------------------------------
 
@@ -42,37 +42,22 @@
       !!              - Initialization of the ice masks (tmsk, umsk)
       !! 
-      !! ** Refer.  : Deleersnijder et al. Ocean Modelling 100, 7-10 
-      !!
-      !! ** History :
-      !!         original    : 01-04 (LIM)
-      !!         addition    : 02-08 (C. Ethe, G. Madec)
+      !! Reference  : Deleersnijder et al. Ocean Modelling 100, 7-10 
       !!--------------------------------------------------------------------- 
-      !! * Local variables
-      INTEGER :: ji, jj      ! dummy loop indices
-
-      REAL(wp), DIMENSION(jpi,jpj) ::  &
-         zd2d1 , zd1d2       ! Derivative of zh2 (resp. zh1) in the x direction
-      !                   ! (resp. y direction) (defined at the center)
-      REAL(wp) ::         &
-         zh1p  , zh2p   , &  ! Idem zh1, zh2 for the bottom left corner of the grid
-         zd2d1p, zd1d2p , &  ! Idem zd2d1, zd1d2 for the bottom left corner of the grid
-         zusden, zusden2     ! temporary scalars
+      INTEGER  ::   ji, jj   ! dummy loop indices
+      REAL(wp) ::   zusden   ! temporary scalar
       !!---------------------------------------------------------------------
 
       IF(lwp) THEN
          WRITE(numout,*)
-         WRITE(numout,*) 'lim_msh : LIM sea-ice model, mesh initialization'
+         WRITE(numout,*) 'lim_msh : LIM-3 sea-ice model, mesh initialization'
          WRITE(numout,*) '~~~~~~~'
       ENDIF
 
-      !----------------------------------------------------------                          
-      !    Initialization of local and some global (common) variables 
-      !------------------------------------------------------------------ 
-
-      njeq   = INT( jpj / 2 )   !i bug mpp potentiel
+      !                           !==  coriolis factor & Equator position ==!
+      njeq   = INT( jpj / 2 ) 
       njeqm1 = njeq - 1 
-
-      fcor(:,:) = 2. * omega * SIN( gphit(:,:) * rad )   !  coriolis factor
-
+      !
+      fcor(:,:) = 2. * omega * SIN( gphit(:,:) * rad )   ! coriolis factor
+      !
       IF( fcor(1,1) * fcor(1,nlcj) < 0.e0 ) THEN   ! local domain include both hemisphere
          l_jeq = .TRUE.
@@ -91,49 +76,15 @@
          IF(lwp ) WRITE(numout,*) '          the model domain is entirely in the northern hemisphere: njeq = ', njeq
       ENDIF
-
+      !
       njeqm1 = njeq - 1
 
 
-      !   For each grid, definition of geometric tables 
-      !------------------------------------------------------------------
-
-      !-------------------
-      ! Conventions :    !
-      !-------------------
-      !  indices 1 \ 2 <-> localisation in the 2 direction x \ y
-      !  3rd indice <-> localisation on the mesh :
-      !  0 = Centre ;  1 = corner W x(i-1/2) ; 2 = corner S y(j-1/2) ;
-      !  3 = corner SW x(i-1/2),y(j-1/2)
-      !-------------------
-      !!ibug ???
-      akappa(:,:,:,:) = 0.e0
+      !                           !==  metric coefficients for sea ice dynamic  ==!
       wght(:,:,:,:) = 0.e0
-      alambd(:,:,:,:,:,:) = 0.e0
-      tmu(:,:) = 0.e0
-      tmv(:,:) = 0.0e0 ! CGrid EVP
-      !!i
-      ! metric coefficients for sea ice dynamic
-      !----------------------------------------
-      !                                                       ! akappa
-      DO jj = 2, jpj
-         zd1d2(:,jj) = e1v(:,jj) - e1v(:,jj-1)
-      END DO
-      CALL lbc_lnk( zd1d2, 'T', -1. )
-
-      DO ji = 2, jpi
-         zd2d1(ji,:) = e2u(ji,:) - e2u(ji-1,:)
-      END DO
-      CALL lbc_lnk( zd2d1, 'T', -1. )
-
-      akappa(:,:,1,1) =        1.0 / ( 2.0 * e1t(:,:) )
-      akappa(:,:,1,2) = zd1d2(:,:) / ( 4.0 * e1t(:,:) * e2t(:,:) )
-      akappa(:,:,2,1) = zd2d1(:,:) / ( 4.0 * e1t(:,:) * e2t(:,:) )
-      akappa(:,:,2,2) =        1.0 / ( 2.0 * e2t(:,:) )
-
-      !                                                      ! weights (wght)
+!!gm  Optimisation :  wght to be defined at F-point, not I-point  and change in limrhg
       DO jj = 2, jpj
          DO ji = 2, jpi
-            zusden = 1. / (  ( e1t(ji,jj) + e1t(ji-1,jj  ) )   &
-               &           * ( e2t(ji,jj) + e2t(ji  ,jj-1) ) )
+            zusden = 1.e0 / (  ( e1t(ji,jj) + e1t(ji-1,jj  ) )   &
+               &             * ( e2t(ji,jj) + e2t(ji  ,jj-1) ) )
             wght(ji,jj,1,1) = zusden * e1t(ji  ,jj) * e2t(ji,jj  )
             wght(ji,jj,1,2) = zusden * e1t(ji  ,jj) * e2t(ji,jj-1)
@@ -146,96 +97,20 @@
       CALL lbc_lnk( wght(:,:,2,1), 'I', 1. )      ! but it is never used
       CALL lbc_lnk( wght(:,:,2,2), 'I', 1. )
+!!gm end
 
-      ! Coefficients for divergence of the stress tensor
-      !-------------------------------------------------
-
-      DO jj = 2, jpj
-         DO ji = 2, jpi
-            zh1p  =  e1t(ji  ,jj  ) * wght(ji,jj,2,2)   &
-               &   + e1t(ji-1,jj  ) * wght(ji,jj,1,2)   &
-               &   + e1t(ji  ,jj-1) * wght(ji,jj,2,1)   &
-               &   + e1t(ji-1,jj-1) * wght(ji,jj,1,1)
-
-            zh2p  =  e2t(ji  ,jj  ) * wght(ji,jj,2,2)   &
-               &   + e2t(ji-1,jj  ) * wght(ji,jj,1,2)   &
-               &   + e2t(ji  ,jj-1) * wght(ji,jj,2,1)   &
-               &   + e2t(ji-1,jj-1) * wght(ji,jj,1,1)
-
-            zusden = 1.0 / MAX( zh1p * zh2p * 4.e0 , 1.e-20 )
-            zusden2 = zusden * 2.0 
-
-            zd1d2p = zusden * 0.5 * ( -e1t(ji-1,jj-1) + e1t(ji-1,jj  ) - e1t(ji,jj-1) + e1t(ji  ,jj)   )
-            zd2d1p = zusden * 0.5 * (  e2t(ji  ,jj-1) - e2t(ji-1,jj-1) + e2t(ji,jj  ) - e2t(ji-1,jj)   )
-
-            alambd(ji,jj,2,2,2,1) = zusden2 * e2t(ji  ,jj-1)
-            alambd(ji,jj,2,2,2,2) = zusden2 * e2t(ji  ,jj  )
-            alambd(ji,jj,2,2,1,1) = zusden2 * e2t(ji-1,jj-1)
-            alambd(ji,jj,2,2,1,2) = zusden2 * e2t(ji-1,jj  )
-
-            alambd(ji,jj,1,1,2,1) = zusden2 * e1t(ji  ,jj-1)
-            alambd(ji,jj,1,1,2,2) = zusden2 * e1t(ji  ,jj  )
-            alambd(ji,jj,1,1,1,1) = zusden2 * e1t(ji-1,jj-1)
-            alambd(ji,jj,1,1,1,2) = zusden2 * e1t(ji-1,jj  )
-
-            alambd(ji,jj,1,2,2,1) = zd1d2p
-            alambd(ji,jj,1,2,2,2) = zd1d2p
-            alambd(ji,jj,1,2,1,1) = zd1d2p
-            alambd(ji,jj,1,2,1,2) = zd1d2p
-
-            alambd(ji,jj,2,1,2,1) = zd2d1p
-            alambd(ji,jj,2,1,2,2) = zd2d1p
-            alambd(ji,jj,2,1,1,1) = zd2d1p
-            alambd(ji,jj,2,1,1,2) = zd2d1p
+      !                           !==  ice masks  ==!
+      tms(:,:) = tmask(:,:,1)             ! ice T-point  : use surface tmask
+      tmu(:,:) = umask(:,:,1)             ! ice U-point  : use surface umask  (C-grid EVP)
+      tmv(:,:) = vmask(:,:,1)             ! ice V-point  : use surface vmask  (C-grid EVP)
+      DO jj = 1, jpjm1                    ! ice F-point  : recompute fmask (due to nn_shlat)
+         DO ji = 1 , jpim1
+            tmf(ji,jj) =  tms(ji,jj) * tms(ji+1,jj) * tms(ji,jj+1) * tms(ji+1,jj+1)
          END DO
       END DO
+      CALL lbc_lnk( tmf(:,:), 'F', 1. )           ! lateral boundary conditions
 
-      CALL lbc_lnk( alambd(:,:,2,2,2,1), 'I', 1. )      ! CAUTION: even with the lbc_lnk at ice U-V point
-      CALL lbc_lnk( alambd(:,:,2,2,2,2), 'I', 1. )      ! the value of wght at jpj is wrong
-      CALL lbc_lnk( alambd(:,:,2,2,1,1), 'I', 1. )      ! but it is never used
-      CALL lbc_lnk( alambd(:,:,2,2,1,2), 'I', 1. )      ! 
-
-      CALL lbc_lnk( alambd(:,:,1,1,2,1), 'I', 1. )      ! CAUTION: idem
-      CALL lbc_lnk( alambd(:,:,1,1,2,2), 'I', 1. )      ! 
-      CALL lbc_lnk( alambd(:,:,1,1,1,1), 'I', 1. )      !
-      CALL lbc_lnk( alambd(:,:,1,1,1,2), 'I', 1. )      !
-
-      CALL lbc_lnk( alambd(:,:,1,2,2,1), 'I', 1. )      ! CAUTION: idem
-      CALL lbc_lnk( alambd(:,:,1,2,2,2), 'I', 1. )      !
-      CALL lbc_lnk( alambd(:,:,1,2,1,1), 'I', 1. )      !
-      CALL lbc_lnk( alambd(:,:,1,2,1,2), 'I', 1. )      !
-
-      CALL lbc_lnk( alambd(:,:,2,1,2,1), 'I', 1. )      ! CAUTION: idem
-      CALL lbc_lnk( alambd(:,:,2,1,2,2), 'I', 1. )      !
-      CALL lbc_lnk( alambd(:,:,2,1,1,1), 'I', 1. )      !
-      CALL lbc_lnk( alambd(:,:,2,1,1,2), 'I', 1. )      !
-
-
-      ! Initialization of ice masks
-      !----------------------------
-
-      tms(:,:) = tmask(:,:,1)      ! ice T-point  : use surface tmask
-
-      !      tmu(:,1) = 0.e0
-      !      tmu(1,:) = 0.e0
-      !      tmv(:,1) = 0.e0
-      !      tmv(1,:) = 0.e0
-
-      DO jj = 1, jpj - 1
-         DO ji = 1 , jpi - 1
-            tmu(ji,jj) =  tms(ji,jj) * tms(ji+1,jj)
-            tmv(ji,jj) =  tms(ji,jj) * tms(ji,jj+1)
-            tmf(ji,jj) =  tms(ji,jj) * tms(ji+1,jj) * tms(ji,jj+1) * &
-               tms(ji+1,jj+1)
-         END DO
-      END DO
-
-      !--lateral boundary conditions    
-      CALL lbc_lnk( tmu(:,:), 'U', 1. )
-      CALL lbc_lnk( tmv(:,:), 'V', 1. )
-      CALL lbc_lnk( tmf(:,:), 'F', 1. )
-
-      ! unmasked and masked area of T-grid cell
+      !                           !==  unmasked and masked area of T-grid cell
       area(:,:) = e1t(:,:) * e2t(:,:)
-
+      !
    END SUBROUTINE lim_msh
 
Index: trunk/NEMO/LIM_SRC_3/limthd.F90
===================================================================
--- trunk/NEMO/LIM_SRC_3/limthd.F90	(revision 1607)
+++ trunk/NEMO/LIM_SRC_3/limthd.F90	(revision 1608)
@@ -98,5 +98,5 @@
       ! 1) Initialization of diagnostic variables                                    !
       !------------------------------------------------------------------------------!
-      zeps = 1.0e-10
+      zeps = 1.e-10
 
       !--------------------
@@ -478,10 +478,9 @@
       !!
       INTEGER  ::   ji,jk   ! loop indices
-      REAL(wp) ::   zdes    ! snow heat content increment (dummy)
       REAL(wp) ::   zeps    ! very small value (1.e-10)
       !!-----------------------------------------------------------------------
       eti(:,:) = 0.e0
       ets(:,:) = 0.e0
-      zeps     = 1.0e-10
+      zeps     = 1.e-10
 
       DO jk = 1, nlay_i                ! total q over all layers, ice [J.m-2]
Index: trunk/NEMO/LIM_SRC_3/par_ice.F90
===================================================================
--- trunk/NEMO/LIM_SRC_3/par_ice.F90	(revision 1607)
+++ trunk/NEMO/LIM_SRC_3/par_ice.F90	(revision 1608)
@@ -2,27 +2,25 @@
    !!======================================================================
    !!                       ***  MODULE par_ice   ***
-   !! Sea-Ice model : definition of the parameters
+   !! LIM-3 Sea Ice :   definition of parameters
    !!======================================================================
+   !! History :  3.0  ! 2003-08  (M. Vancoppenolle)  LIM-3
    !!----------------------------------------------------------------------
-   !!  LIM 3.0, UCL-LOCEAN-IPSL (2008)
-   !! $Id$
-   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
-   !!----------------------------------------------------------------------
-   !! * Modules used
    USE par_oce
 
    IMPLICIT NONE
-   PUBLIC               ! allows par_oce and par_kind to be known in ice modules
+   PUBLIC          
 
-   !ICE THERMODYNAMICS
-   INTEGER , PARAMETER ::      &  !:
-      jkmax    = 6 ,           &  !: maximum number of ice layers
-      nlay_s   = 1                !: number of snow layers
+   !                                             !!! ice thermodynamics
+   INTEGER, PUBLIC, PARAMETER ::   jkmax    = 6   !: maximum number of ice layers
+   INTEGER, PUBLIC, PARAMETER ::   nlay_s   = 1   !: number of snow layers
 
-   !ICE MECHANICAL REDISTRIBUTION
-   INTEGER , PARAMETER ::      &  !:
-      jpl      = 5,            &  !: number of ice categories
-      jpm      = 1                !: number of ice types
+   !                                             !!! ice mechanical redistribution
+   INTEGER, PUBLIC, PARAMETER ::   jpl      = 5   !: number of ice categories
+   INTEGER, PUBLIC, PARAMETER ::   jpm      = 1   !: number of ice types
 
+   !!----------------------------------------------------------------------
+   !! NEMO/LIM 3.2, UCL-ASTR-LOCEAN-IPSL (2009)
+   !! $Id$
+   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
    !!======================================================================
 END MODULE par_ice
