Index: /branches/GRISLIv3/SOURCES/3D-physique-gen_mod.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/3D-physique-gen_mod.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/3D-physique-gen_mod.f90	(revision 467)
@@ -63,5 +63,4 @@
   real ::  V_limit                      !< vitesse maxi (limitateur de flux pour conserv masse)
   real ::  betamax                      !< (Pa) frottement maxi sous les streams
-  real ::  pvimin                       !< valeur de pvi pour les noeuds fictifs
 
 
@@ -77,6 +76,4 @@
 
   integer,dimension(nx,ny) :: FRONT     !< point d'un shelf situe sur le front
-  integer,dimension(nx,ny) :: FRONTFACEX!< type de front sur les faces x
-  integer,dimension(nx,ny) :: FRONTFACEY!< type de front sur les faces y
   integer,dimension(nx,ny) :: gr_line_schoof ! points ou on impose le flux de schoof (pour sorties)
   integer,dimension(nx,ny) :: gr_line   !< points grounding line pour les sorties
@@ -88,5 +85,4 @@
   integer,dimension(nx,ny) :: Mk_init      !< initial mask (with islands, outcrops, ...
   integer,dimension(nx,ny) :: MK           !< masks (ice sheet, max, above water, below water, 1)
-  integer,dimension(nx,ny) :: MNEG         !< masks (ice sheet, max, above water, below water, 1)
   integer,dimension(nx,ny) :: IBASE     !< type de base (froide, temperee)
 
@@ -95,14 +91,9 @@
   real,dimension(nx,ny) :: ABLBORD      !< Ablation dans les zones non englacees
   real,dimension(nx,ny) :: ablbord_dtt  !< Ablation sur le bord de la calotte (flux de glace inférieur à ablation sur le bord) sur dtt
-  real,dimension(nx,ny) :: ACQUA        !< Surface des surfaces en eau
   double precision,dimension(nx,ny) :: B            !< Altitude de la base de la glace  'o'
   real,dimension(nx,ny) :: BDOT         !< derivee de B / t
-  real,dimension(nx,ny) :: B1           !< 
   double precision,dimension(nx,ny) :: BSOC0        !< ice free bedrock
-  real,dimension(nx,ny) :: betamx       !< coefficient beta du frottement basal
-  ! betamx=cf*neffmx et tobmx=betamx*uxbar
-  real,dimension(nx,ny) :: betamy       !< coefficient beta du frottement basal
-  ! betamy=cf*neffmy et tobmy=betamy*uybar
-  !hassine
+  real,dimension(nx,ny) :: betamx       !< coefficient beta du frottement basal, betamx=cf*neffmx et tobmx=betamx*uxbar
+  real,dimension(nx,ny) :: betamy       !< coefficient beta du frottement basal, betamy=cf*neffmy et tobmy=betamy*uybar
   real, dimension(nx,ny) :: beta_centre !< beta on major node (average)
   real, dimension(nx,ny) :: betamax_2d  !< (Pa) frottement maxi sous les streams
@@ -116,6 +107,4 @@
   real,dimension(nx,ny) :: CALV         !< calving
   
-  !real,dimension(:,:,:),allocatable :: DDX  ! pour le calcul de UX 
-  !real,dimension(:,:,:),allocatable :: DDY  ! pour le calcul de UY
   real,dimension(nx,ny) :: DDBX         !< UBX/(-SDX) '>'
   real,dimension(nx,ny) :: DDBY         !< UBY/(-SDY) '^'
@@ -123,13 +112,5 @@
   real,dimension(nx,ny) :: Diffmy       !< partie diffusive en y
 
-  !REAL,dimension(nx,ny) :: DPHIX        !< ROG*HMX**2  '>'
-  !REAL,dimension(nx,ny) :: DPHIY        !< ROG*HMY**2  '^'
   real,dimension(nx,ny) :: DHDT         !< derive lagrangienne de l'epaisseur 'o'
-  real,dimension(nx,ny) :: DIF1         !< UZR(I,J,1)-BM(I,J)
-  real,dimension(nx,ny) :: DIF2         !< UZR(I,J,1)-BMELT(I,J)
-  real,dimension(nx,ny) :: DIFT         !< pour sorties eismint
-  real,dimension(nx,ny) :: DIVU         !< divergence de UH
-  ! real,dimension(nx,ny) :: DISTCENT     ! 
-  real,dimension(nx,ny) :: ELA          !< Equilibrium line altitude, for mass balance
   real,dimension(nx,ny) :: epsxx     !< vitesse de deformation selon x
   real,dimension(nx,ny) :: epsyy     !< vitesse de deformation selon y
@@ -141,5 +122,4 @@
   real,dimension(nx,ny) :: FROTMY       !< 
   real,dimension(nx,ny) :: GHF          !< geothermal heat flux J/m2/a 'o'
-  real,dimension(nx,ny) :: GHF0         !< geothermal heat flux J/m2/a 'o'
   double precision,dimension(nx,ny) :: H            !< ice thickness  'o'
   double precision,dimension(nx,ny) :: H0           !< initial ice thickness, must be initialized before init_iso
@@ -148,35 +128,19 @@
   real,dimension(nx,ny)    :: Hp           !< H value if prescribed
   real,dimension(nx,ny)    :: Hp0          !< H value if prescribed (reference value)
-!  real,dimension(nx,ny)    :: Delta_H      !< Delta_H value if prescribed
-!  integer,dimension(nx,ny) :: i_delta_H    !< 1 if Delta_H is prescribed on this node, else 0
   integer,dimension(nx,ny) :: i_Hp         !< 1 if H is prescribed on this node, else 0
   integer,dimension(nx,ny) :: i_Hp0        !< i_hp mask reference value does not change with time
   integer, dimension(nx,ny) :: imx_diag    !< masque pour eq elliptique  
   integer, dimension(nx,ny) :: imy_diag    !< masque pour eq elliptique 
-  integer,dimension(nx,ny) :: MK_gl0       !< mask grounding line initial
-  integer,dimension(nx,ny) :: MK_flot0     !< mask float initial
-
 
   double precision,dimension(nx,ny) :: HDOT         !< ice thickness derivee / t
   real,dimension(nx,ny) :: HDOTWATER
-  real,dimension(nx,ny) :: H1           !< 
   real,dimension(nx,ny) :: HMX          !< ice thickness moy selon x '>'
   real,dimension(nx,ny) :: HMY          !< ice thickness moy selon y '^'
-  real,dimension(nx,ny) :: HDOTPREC     !< 
-  real,dimension(nx,ny) :: HDOTRAP      !< 
-
-  real,dimension(nx,ny) :: KOND         !< basale hydro. conductivity 'o'
-  real,dimension(nx,ny) :: PDD          !< Positive degree day
-  real,dimension(nx,ny) :: PRECIP       !< precipitation
-  real,dimension(nx,ny) :: PRECIP0      !< initial precipitation  (used in 'heminord')
+
   real,dimension(nx,ny) :: PHID         !< flux de chaleur lie a la deformation et glissement basal
-  real,dimension(nx,ny) :: PWATER       !< basal water pressure
-  real,dimension(nx,ny) :: pgx,pgy      !< hydro. potential gradient  / x '>' and / y '^'
-  real,dimension(nx,ny) :: phiWx,phiWy  !< flux d'eau sous glaciaire  / x '>' and / y '^'
   real,dimension(nx,ny) :: NEFFMX       !< pression effective '>'
   real,dimension(nx,ny) :: NEFFMY       !< pression effective '^'
   real,dimension(nx,ny) :: TOBMX        !< cisaillement basal '>'
   real,dimension(nx,ny) :: TOBMY        !< cisaillement basal '^'
-  real,dimension(nx,ny) :: SW           !< for bedrock isostasy
   double precision,dimension(nx,ny) :: S            !< altitude of ice sheet surface
   real,dimension(nx,ny) :: sealevel_2d  !< local sea surface elevation
@@ -184,6 +148,4 @@
   real,dimension(nx,ny) :: SDX          !< slope derivee / x '>'
   real,dimension(nx,ny) :: SDY          !< slope derivee / y '^'
-  real,dimension(nx,ny) :: SDXMY        !< slope selon x moy selon y '^' remplace SDMX
-  real,dimension(nx,ny) :: SDYMX        !< slope selon y moy selon x '>' remplace SDMY
   real,dimension(nx,ny) :: SLOPE2mx     !< = Sdx**2 + Sdymx**2 '>'
   real,dimension(nx,ny) :: SLOPE2my     !< = Sdy**2 + Sdxmy**2 '^'
@@ -193,12 +155,8 @@
   real,dimension(nx,ny) :: TANN         !< Ground air temperature annual
   real,dimension(nx,ny,12) :: Tmois        !< Ground air temperature monthly
-  real,dimension(nx,ny) :: TSHELF       !< temperature des shelfs pour viscosite
-  real,dimension(nx,ny) :: TJ0          !< initial air temperature at sea level July
-!  real,dimension(nx,ny) :: TA0          !< initial air temperature at sea level annual
   real,dimension(nx,ny) :: TAUB         !< basal shear stress (for output)
   real,dimension(nx,ny) :: TAUSHELF     !< effective stress in ice shelves (vertical av.) 
   real,dimension(nx,ny) :: TS           !< surface ice temperature  'o'
   real,dimension(nx,ny) :: TB           !< basal ice temperature  'o'
-  real,dimension(nx,ny) :: TG           !< degrees above melting point at the base
   real,dimension(nx,ny) :: TBDOT        !< variation in time of basal temperature 
   real,dimension(nx,ny) :: UZK          !< vertical velocity at the ice surface (kinematic boundary)
@@ -212,23 +170,12 @@
   real,dimension(nx,ny) :: UBY          !< basal sliding '^'
   real,dimension(nx,ny) :: UZSDOT       !< variation in time of surface vertical velocity 
-  real,dimension(nx,ny) :: UX1          !< 
-  real,dimension(nx,ny) :: UY1          !< 
-  real,dimension(nx,ny) :: VBAR         !< depth averaged velocity magnitude (for output)
-!!!!!!!real,dimension(nx,ny) :: VSDOT        ! variation in time of surface velocity magnitude 
-  real,dimension(nx,ny) :: W0           !< enfoncement du socle a l'equilibre isostatique
-  real,dimension(nx,ny) :: W1           !< enfoncement du socle courant
   real,dimension(nx,ny) :: XLONG        !< longitude
   real,dimension(nx,ny) :: YLAT         !< latitude
   real,dimension(nx,ny) :: xcc          ! grille Xkm
   real,dimension(nx,ny) :: ycc          ! grille Ykm
-!  real,dimension(nx,ny) :: ZS           !< surface topography above sea level
-  real,dimension(nx,ny) :: ETABAR       !< 
   real,dimension(nx,ny) :: BSOC         !< altitude (ou bathymetrie) du socle 'o'
   real,dimension(nx,ny) :: PVI          !< viscosite ice shelf
   real,dimension(nx,ny) :: PVM          !< viscosite ice shelf
-  real,dimension(nx,ny) :: ramollo      !< pour ramollir les ice shelves
   real,dimension(nx,ny) :: Abar         !< coefficient de Glen integre
-  real,dimension(nx,ny) :: Uiter_centre !< pour iterations equation diagnostique
-  real,dimension(nx,ny) :: tabtest      !< tableau de travail
 
 
@@ -241,29 +188,9 @@
   integer,dimension(nx,ny) :: mstream_my  !< masque stream selon y
   integer,dimension(nx,ny) :: mstream     !< masque stream sur les noeuds majeurs
-  real,dimension(nx,ny) :: socle_cry      !< courbure du socle, negatif -> vallees
   real,dimension(nx,ny) :: drag_mx        !< coefficient pour tenir compte de variations        
   real,dimension(nx,ny) :: drag_my        !< geographiques eventuelles dans le basal drag
 
-  integer,dimension(nx,ny) :: mslid_mx    !< masque glissement selon x
-  integer,dimension(nx,ny) :: mslid_my    !< masque glissement selon y
-  real,dimension(nx,ny) :: slid_mx        !< coefficient pour tenir compte de variations        
-  real,dimension(nx,ny) :: slid_my        !< geographiques eventuelles dans le glissement
-
-
-
-  !REAL,dimension(-NL:NL,-NL:NL) :: WE   ! enfoncement du socle autour d'une charge unitaire
-
-
-
   !      ********** common des tableaux a 3 dimensions *****
-  !real,dimension(:,:,:,:),allocatable :: BTT !< deformation parameter for flow law 
-
-  !real,dimension(:,:,:,:),allocatable :: SA ! effet temperature sur la deformation (loi n=3)
-  !real,dimension(:,:,:,:),allocatable :: S2A ! effet integre sur l'epaisseur temperature sur la deformation (loi n=3)
-  !real,dimension(:,:,:,:),allocatable :: SA_mx ! effet temperature sur la deformation (loi n=3)
-  !real,dimension(:,:,:,:),allocatable :: S2A_mx ! effet integre sur l'epaisseur temperature sur la deformation (loi n=3)
-
-  !real,dimension(:,:,:,:),allocatable :: SA_my ! effet temperature sur la deformation (loi n=3)
-  !real,dimension(:,:,:,:),allocatable :: S2A_my ! effet integre sur l'epaisseur temperature sur la deformation (loi n=3)
+
   real,dimension(nx,ny,nz) :: SUX       !< 
   real,dimension(nx,ny,nz) :: SUY       !< 
@@ -299,7 +226,5 @@
   ! ===================== Booleens ========================================
   logical :: SHELFY                     !< 
-  logical :: BOOST                      !< 
   logical :: MARINE                     !< 
-  logical :: appel_new_flot             !< pour appeler la routine new_flot
 
   logical,dimension(nx,ny) :: FLOT      !< vrai si flottant (test d'archimede) 'o'
@@ -308,7 +233,4 @@
   logical,dimension(nx,ny) :: FRONT_RESC!< comme front (mais boolean) pour remplimat rescue
   logical,dimension(nx,ny) :: FBM       !< 
-  logical,dimension(nx,ny) :: OKUMAT    !< 
-  logical,dimension(nx,ny) :: OKVMAT    !< 
-  logical,dimension(nx,ny) :: GRZONE    !< 
   logical,dimension(nx,ny) :: GZMX      !< point stream ">"
   logical,dimension(nx,ny) :: GZMY      !< point stream "^"
@@ -317,13 +239,6 @@
   logical,dimension(nx,ny) :: ILEMX     !< points ile ">"
   logical,dimension(nx,ny) :: ILEMY     !< points ile "^"
-  logical,dimension(nx,ny) :: cotemx    !< points cotiers ">"
-  logical,dimension(nx,ny) :: cotemy    !< points cotiers "^"
   logical,dimension(nx,ny) :: fleuvemx  !< actual grounded stream
   logical,dimension(nx,ny) :: fleuvemy  !< actual grounded stream
-  logical,dimension(nx,ny) :: isolx     !< designe un point sans voisin sur ses faces x
-  logical,dimension(nx,ny) :: isoly     !< designe un point sans voisin sur ses faces y
-  logical,dimension(nx,ny) :: new_flot_point  !< pour signaler les points qui se mettent a flotter entre 2 pas de temps dtt
-  logical,dimension(nx,ny) :: new_flotmx !< pour signaler les points qui deviennent flottantmx entre 2 dtt 
-  logical,dimension(nx,ny) :: new_flotmy !< pour signaler les points qui deviennent flottantmy entre 2 dtt 
   logical,dimension(nx,ny) :: flot_marais  !< afq -- vrai si flottant et coince entre points poses 'o'
   logical,dimension(nx,ny) :: iceberg   !< point iceberg
@@ -359,9 +274,6 @@
   integer :: num_file4     = 995  !< Id of 
   integer :: num_coor      = 2004 !< Id of coord-Ant-40km
-!  integer :: num_tracebug        !< numero de l'unite itracebug
 
   ! Variables communes au main et aux subroutines du step
-  logical :: base_froide
-  logical :: base_temp
   real  :: timemax
 
Index: /branches/GRISLIv3/SOURCES/Climate_modules/climat-perturb_mod-0.4.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/Climate_modules/climat-perturb_mod-0.4.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/Climate_modules/climat-perturb_mod-0.4.f90	(revision 467)
@@ -13,5 +13,5 @@
 
 
-use module3d_phy,only: S,S0,Tann,Tjuly,precip,acc,Ylat,num_forc,num_param,num_rep_42,tafor,time,sealevel,sealevel_2d,coefbmshelf
+use module3d_phy,only: S,S0,Tann,Tjuly,acc,Ylat,num_forc,num_param,num_rep_42,tafor,time,sealevel,sealevel_2d,coefbmshelf
 use geography, only: nx,ny,dirforcage,dirnameinp
 use io_netcdf_grisli, only: read_ncdf_var
@@ -25,4 +25,5 @@
 
 real,dimension(nx,ny) :: ta0          !< initial air temperature at sea level annual
+real,dimension(nx,ny) :: precip
 
 real :: coefT                   !< pour modifier l'amplitude de la perturb. T
Index: /branches/GRISLIv3/SOURCES/Netcdf-routines/sortie_netcdf_GRISLI_mod.0.2-hassine.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/Netcdf-routines/sortie_netcdf_GRISLI_mod.0.2-hassine.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/Netcdf-routines/sortie_netcdf_GRISLI_mod.0.2-hassine.f90	(revision 467)
@@ -665,11 +665,12 @@
     use geography,    only: nx,ny,nz,nzm,dx,dy
     use runparam,     only: itracebug,xmin,ymin,dirsource
-    use module3d_phy, only: s,s0,h,h0,b,socle_cry,sealevel,mk,mk_init,hdot,dhdt,tann,tjuly,t,tpmp, &
+    use module3d_phy, only: s,s0,h,h0,b,sealevel,mk,mk_init,hdot,dhdt,tann,tjuly,t,tpmp, &
                             acc,abl,ghf,phid,bmelt,uxbar,uybar,ux,uy,frotmx,frotmy,tobmx,tobmy,    &
                             taushelf,epsxx,epsyy,epsxy,eps,abar,pvi,pvm,betamx,betamy,beta_centre, &
-                            ablbord_dtt,ice,front,gr_line_schoof,hwater,hdotwater,kond,bsoc,bm,    &
-                            pgx,pgy,phiwx,phiwy,flot_marais,neffmx,neffmy,gzmx,gzmy,ilemx,ilemy,   &
-                            fleuvemx,fleuvemy,flotmx,flotmy,hmx,hmy,frontfacex,frontfacey,sux,suy, &
+                            ablbord_dtt,ice,front,gr_line_schoof,hwater,hdotwater,bsoc,bm,    &
+                            flot_marais,neffmx,neffmy,gzmx,gzmy,ilemx,ilemy,   &
+                            fleuvemx,fleuvemy,flotmx,flotmy,hmx,hmy,sux,suy, &
                             tpmp,ux,uy,uzr,debug_3d,xlong,ylat,time,dtt
+    use eau_basale, only: kond,pgx,pgy,phiwx,phiwy
     use tracer_vars,  only: xdep_out,ydep_out,tdep_out
     use bilan_eau_mod,only: tot_water,calv_dtt
@@ -788,7 +789,7 @@
                    tab(:,:) = b(:,:)
                 end if
-                if (itab.eq.8) then
-                   tab(:,:) = socle_cry(:,:)
-                end if
+                !if (itab.eq.8) then
+                !   tab(:,:) = socle_cry(:,:)  !no longer computed
+                !end if
                 if (itab.eq.9) then
                    tab(:,:) = mk_init(:,:)
@@ -1015,10 +1016,10 @@
                    end do
                 end if
-                if (itab.eq.72) then 
-                   tab(:,:) = frontfacex(:,:)
-                end if
-                if (itab.eq.73) then 
-                   tab(:,:) = frontfacey(:,:)
-                end if
+                !if (itab.eq.72) then !no longer computed
+                !   tab(:,:) = frontfacex(:,:)
+                !end if
+                !if (itab.eq.73) then 
+                !   tab(:,:) = frontfacey(:,:)
+                !end if
 
                 !SORTIE 3D 
Index: /branches/GRISLIv3/SOURCES/New-remplimat/diagno-L2_mod.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/New-remplimat/diagno-L2_mod.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/New-remplimat/diagno-L2_mod.f90	(revision 467)
@@ -6,5 +6,5 @@
 implicit none
 
-
+real ::  pvimin                       !< valeur de pvi pour les noeuds fictifs
 
 real                   :: somint,test,delp,prec
@@ -12,4 +12,5 @@
 real, dimension(nx,ny) :: uyb1
 
+real, dimension(nx,ny) :: ramollo      !< pour ramollir les ice shelves
 real, dimension(nx,ny) :: uxb1ramollo
 real, dimension(nx,ny) :: uyb1ramollo
@@ -49,5 +50,5 @@
 subroutine init_diagno
 
-use module3D_phy,only: num_rep_42,num_param,pvimin
+use module3D_phy,only: num_rep_42,num_param
 use geography, only: geoplace
 use deformation_mod_2lois,only: sf ! afq - 07/2023 - not clean since it is in module_choix
@@ -111,5 +112,5 @@
 
   use module3D_phy, only: niter_nolin,taushelf,flot,h,pvi, & 
-                          gr_line_schoof,uxbar,uybar,pvimin,debug_3D,V_limit, & 
+                          gr_line_schoof,uxbar,uybar,debug_3D,V_limit, & 
                           imx_diag,imy_diag,tobmx,tobmy,neffmx,neffmy,gr_line, schoof, & 
                           betamx,betamy,flgzmx,flgzmy,uxflgz,uyflgz,taub 
@@ -378,6 +379,6 @@
 subroutine calc_pvi
 
-use module3d_phy, only: pvi,pvimin,abar,flot,gzmx,gzmy, & 
-                        ilemx,ilemy,eps,taushelf,h,ramollo,debug_3d,pvm
+use module3d_phy, only: pvi,abar,flot,gzmx,gzmy, & 
+                        ilemx,ilemy,eps,taushelf,h,debug_3d,pvm
 use runparam, only: itracebug
 use geography, only: nx,ny,nz
Index: /branches/GRISLIv3/SOURCES/New-remplimat/remplimat-shelves-tabTu.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/New-remplimat/remplimat-shelves-tabTu.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/New-remplimat/remplimat-shelves-tabTu.f90	(revision 467)
@@ -44,5 +44,5 @@
 use module3d_phy, only: debug_3D,hmx,hmy,pvi,  &
                         flgzmx,flgzmy,frotmx,frotmy,betamx,betamy,betamax,betamax_2d, &
-                        sdx,sdy,h,b,pvimin,sealevel_2d,flotmx,flotmy,drag_mx,drag_my, &
+                        sdx,sdy,h,b,sealevel_2d,flotmx,flotmy,drag_mx,drag_my, &
                         beta_centre,pvm
 use runparam, only: itracebug,num_tracebug
@@ -1110,4 +1110,6 @@
 ! et qui ne dépendent d'aucun noeud (ligne Tu,Tv nulle sauf diagonale)
 
+use diagno_mod, only: pvimin
+
 implicit none
 
Index: /branches/GRISLIv3/SOURCES/Temperature-routines/icetemp_mod.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/Temperature-routines/icetemp_mod.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/Temperature-routines/icetemp_mod.f90	(revision 467)
@@ -70,5 +70,5 @@
   subroutine icetemp
     !$ use OMP_LIB 
-    use module3D_phy, only: uxbar,uybar,dtt,nzm,time,ux,uy,uzr,ts,h,h1,b,b1, &
+    use module3D_phy, only: uxbar,uybar,dtt,nzm,time,ux,uy,uzr,ts,h,b, &
                              ghf,t,tpmp,phid,flot,ibase,debug_3D,tbdot
     use runparam, only: itracebug,num_tracebug
@@ -223,6 +223,4 @@
     Do J=1,Ny
        Do I=1,Nx
-          H1(I,J)=H(I,J)
-          B1(I,J)=B(I,J)
           Tpmp(I,J,1)=0.
        End Do
Index: /branches/GRISLIv3/SOURCES/ablation_mod.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/ablation_mod.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/ablation_mod.f90	(revision 467)
@@ -84,10 +84,10 @@
 
 
-  use module3d_phy,only:Tjuly,Tann,Tmois,acc,pdd,TS,Tshelf,precip,BM,Abl,S
+  use module3d_phy,only:Tjuly,Tann,Tmois,acc,TS,BM,Abl,S
   use param_phy_mod,only:dice,cl
 
   IMPLICIT NONE
 
-  real, dimension(nx,ny) :: pds, simax, pdsi, sif
+  real, dimension(nx,ny) :: precip, pdd, pds, simax, pdsi, sif
   integer :: i,j,k,mo,nday
   real :: summ
@@ -285,5 +285,4 @@
   TS(:,:)=(TANN(:,:)+26.6*SIF(:,:))
   TS(:,:)=min(0.0,TS(:,:))
-  tshelf(:,:)=TS(:,:)
   Abl(:,:)=BM(:,:)-Acc(:,:)
 
Index: /branches/GRISLIv3/SOURCES/conserv-mass-adv-diff_sept2009_mod.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/conserv-mass-adv-diff_sept2009_mod.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/conserv-mass-adv-diff_sept2009_mod.f90	(revision 467)
@@ -19,5 +19,5 @@
 module equat_adv_diff_2D_vect                          ! Cat nouvelle mouture juin 2009
   use module3D_phy, only: V_limit,num_param,num_rep_42,dx1,mk0,i_Hp,Hp,H,uxbar,uybar,testdiag,&
-       dtmax,dt,dtmin,time,dtt,isynchro,diffmx,diffmy,sdx,sdy,hmx,hmy,flgzmx,flgzmy,flot,tabtest,timemax,&
+       dtmax,dt,dtmin,time,dtt,isynchro,diffmx,diffmy,sdx,sdy,hmx,hmy,flgzmx,flgzmy,flot,timemax,&
        marine,dtdx2,dtdx,bm,bmelt,igrdline,ibmelt_inv,ice,ablbord,hdot
   use geography, only: nx,ny,dx,geoplace
@@ -89,4 +89,5 @@
     implicit none
     integer :: i,j
+    real,dimension(nx,ny) :: tabtest      ! tableau de travail
     real,dimension(nx,ny) :: Dminx,Dminy
     real,dimension(nx,ny) :: Uxdiff,Uydiff           ! vitesse due a la diffusion
Index: /branches/GRISLIv3/SOURCES/eaubasale-0.5_mod.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/eaubasale-0.5_mod.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/eaubasale-0.5_mod.f90	(revision 467)
@@ -34,4 +34,8 @@
   real :: nefflocal
 
+  real,dimension(nx,ny) :: KOND         !< basale hydro. conductivity 'o'
+  real,dimension(nx,ny) :: pgx,pgy      !< hydro. potential gradient  / x '>' and / y '^'
+  real,dimension(nx,ny) :: phiWx,phiWy  !< flux d'eau sous glaciaire  / x '>' and / y '^'
+  real,dimension(nx,ny) :: PWATER       !< basal water pressure
   REAL,dimension(NX,NY) :: limit_hw    !< conditions aux limites
   integer,dimension(NX,NY) :: klimit    !< ou appliquer les conditions 
@@ -57,5 +61,5 @@
   subroutine init_eaubasale
 
-    use module3d_phy, only:num_param,num_rep_42,kond,secyear,hdotwater,pgx,pgy
+    use module3d_phy, only:num_param,num_rep_42,secyear,hdotwater
     
     namelist/eaubasale1/ecoulement_eau,hwatermax,infiltr 
@@ -126,6 +130,6 @@
   subroutine eaubasale !(pwater)   version correspondant à la thèse de Vincent
 
-    use module3d_phy, only:hwater,kond,secyear,flot,sealevel_2D,Bsoc,ibase,S,H,B,bmelt,&
-         debug_3d,flotmx,flotmy,pgx,pgy,phiwx,phiwy,isynchro,dtt,dt,hdotwater,pwater
+    use module3d_phy, only:hwater,secyear,flot,sealevel_2D,Bsoc,ibase,S,H,B,bmelt,&
+         debug_3d,flotmx,flotmy,isynchro,dtt,dt,hdotwater
     use geography, only:dx,dy
     use param_phy_mod, only:rowg,ro,rog,rofreshg,rofresh
Index: /branches/GRISLIv3/SOURCES/flottab2-0.7.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/flottab2-0.7.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/flottab2-0.7.f90	(revision 467)
@@ -104,7 +104,7 @@
     
     use runparam, only:itracebug,nt
-    use module3D_phy, only:shelfy,igrdline,mk_init,flot,H,sealevel_2d,Bsoc,S,H,B,appel_new_flot,&
-         new_flot_point,new_flotmx,new_flotmy,ice,front,frontfacex,frontfacey,isolx,isoly,cotemx,&
-         cotemy,boost,iceberg,uxbar,uybar,mk,gzmx,gzmy,flotmx,flotmy,hmx,hmy,isynchro,ilemx,ilemy,&
+    use module3D_phy, only:shelfy,igrdline,mk_init,flot,H,sealevel_2d,Bsoc,S,H,B,&
+         ice,front,&
+         iceberg,uxbar,uybar,mk,gzmx,gzmy,flotmx,flotmy,hmx,hmy,isynchro,ilemx,ilemy,&
          flgzmx,flgzmy,marine,fbm,bm,bmelt,debug_3D,dt
     use param_phy_mod, only:row,ro
@@ -113,4 +113,8 @@
     implicit none
     
+    logical,dimension(nx,ny) :: new_flot_point  !< pour signaler les points qui se mettent a flotter entre 2 pas de temps dtt
+    logical,dimension(nx,ny) :: new_flotmx !< pour signaler les points qui deviennent flottantmx entre 2 dtt
+    logical,dimension(nx,ny) :: new_flotmy !< pour signaler les points qui deviennent flottantmy entre 2 dtt 
+
     integer :: i,j
     
@@ -141,5 +145,4 @@
     ! a flotter entre 2 dtt
 
-    appel_new_flot=.false.
     !$OMP DO
     do j=1,ny
@@ -157,11 +160,4 @@
     ICE(:,:)=0
     front(:,:)=0
-    frontfacex(:,:)=0
-    frontfacey(:,:)=0
-    isolx(:,:)=.false.
-    isoly(:,:)=.false.
-    cotemx(:,:)=.false.
-    cotemy(:,:)=.false.
-    boost=.false.
     iceberg(:,:)=.false.
     !$OMP END WORKSHARE
@@ -190,5 +186,4 @@
              ex_pose: if ((.not.FLOT(I,J)).and.(isynchro.eq.1)) then  !  il ne flottait pas avant
                 FLOT(I,J)=.true.
-                BOOST=.false.
 
                 if (igrdline.eq.1) then   ! en cas de grounding line prescrite
@@ -220,5 +215,4 @@
              if(FLOT(I,J)) then         !  mais il flottait avant
                 FLOT(I,J)=.false.
-                BOOST=.false.
              endif
              !cdc correction topo pour suivre  variations sealevel
@@ -283,6 +277,5 @@
           if (flotmx(i,j).and.(new_flot_point(i,j).or. &
                new_flot_point(i-1,j))) then
-             appel_new_flot=.true.
-             new_flotmx(i,j)=.true.
+               new_flotmx(i,j)=.true.
           endif
 
@@ -301,5 +294,4 @@
           archim=(Bsoc(i,j)+Bsoc(i-1,j))*0.5-(sealevel_2d(i,j)+sealevel_2d(i-1,j))*0.5+ro/row*Hmx(i,j)
           gzmx(i,j)=gzmx(i,j).and.(archim.le.100.) 
-          cotemx(i,j)=gzmx(i,j)
 
        end do
@@ -323,6 +315,5 @@
           if (flotmy(i,j).and.(new_flot_point(i,j).or.     &
                new_flot_point(i,j-1))) then
-             appel_new_flot=.true.
-             new_flotmy(i,j)=.true.
+               new_flotmy(i,j)=.true.
           endif
 
@@ -339,6 +330,5 @@
 
           archim=(Bsoc(i,j)+Bsoc(i,j-1))*0.5-(sealevel_2d(i,j)+sealevel_2d(i,j-1))*0.5+ro/row*Hmy(i,j)
-          gzmy(i,j)=gzmy(i,j).and.(archim.le.100.) 
-          cotemy(i,j)=gzmy(i,j)
+          gzmy(i,j)=gzmy(i,j).and.(archim.le.100.)
 
        end do
@@ -820,5 +810,5 @@
   subroutine determin_front
 
-    use module3D_phy, only:ice,H,front,frontfacex,frontfacey,isolx,isoly
+    use module3D_phy, only:ice,H,front
 
     
@@ -1039,100 +1029,4 @@
     !$OMP END DO
 
-    !isolx signifie pas de voisins en x
-    !isoly signifie pas de voisins en y
-    !remarque : 
-    !si isolx/y=.true. alors frontfacex/y=0 (a la fois +1 & -1 or +1-1=0)
-
-    ! calcul de frontfacex et isolx
-    !$OMP DO
-    do j=1,ny
-       do i=2,nx-1
-
-          if (front(i,j).ge.1.and.front(i,j).le.3) then    !front(entre 1 et 3)
-
-             if ((ice(i-1,j)+ice(i+1,j)).lt.2) then        ! il y a un front // a x
-
-                if ((ice(i-1,j)+ice(i+1,j)).eq.0) then
-                   isolx(i,j)=.true.
-                elseif (ice(i-1,j).eq.0) then
-                   frontfacex(i,j)=-1                      ! front  i-1 |i  i+1
-                else
-                   frontfacex(i,j)=+1                      ! front  i-1  i| i+1
-                endif
-             endif
-          end if !fin du test il y a un front
-
-       end do
-    end do
-    !$OMP END DO
-
-    ! calcul de frontfacey et isoly
-    !$OMP DO
-    do j=2,ny-1
-       do i=1,nx
-
-          if (front(i,j).ge.1.and.front(i,j).le.3) then   !front(entre 1 et 3)
-
-             if ((ice(i,j-1)+ice(i,j+1)).lt.2) then       ! il y a un front // a y
-
-                if ((ice(i,j-1)+ice(i,j+1)).eq.0) then
-                   isoly(i,j)=.true.                      !front   j-1 |j| j+1
-                elseif (ice(i,j-1).eq.0) then
-                   frontfacey(i,j)=-1                     !front   j-1 |j j+1
-                else
-                   frontfacey(i,j)=+1                     !front   j-1  j| j+1
-                endif
-             endif
-          end if !fin du test il y a un front
-
-       end do
-    end do
-    !$OMP END DO
-
-
-    ! traitement des bords. On considere que l'exterieur n'a pas de glace
-    ! attention ce n'est vrai que sur la grande grille
-
-    !$OMP DO PRIVATE(i)
-    do j=2,ny-1
-       i=1
-       if (front(i,j).ge.1)  then 
-          if (ice(i+1,j).eq.0) then 
-             isolx(i,j)=.true.
-          else
-             frontfacex(i,j)=-1
-          endif
-       end if
-       i=nx
-       if (front(i,j).ge.1)  then 
-          if (ice(i-1,j).eq.0) then 
-             isolx(i,j)=.true.
-          else
-             frontfacex(i,j)=1
-          endif
-       end if
-    end do
-    !$OMP END DO
-
-    !$OMP DO PRIVATE(j)
-    do i=2,nx-1
-       j=1 
-       if (front(i,j).ge.1)  then 
-          if (ice(i,j+1).eq.0) then 
-             isoly(i,j)=.true.
-          else
-             frontfacey(i,j)=-1
-          endif
-       end if
-       j=ny
-       if (front(i,j).ge.1)  then 
-          if (ice(i,j-1).eq.0) then 
-             isoly(i,j)=.true.
-          else
-             frontfacey(i,j)=1
-          endif
-       end if
-    end do
-    !$OMP END DO
     !$OMP END PARALLEL
 
Index: /branches/GRISLIv3/SOURCES/initial2-0.4.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/initial2-0.4.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/initial2-0.4.f90	(revision 467)
@@ -17,10 +17,10 @@
 subroutine INITIAL2()
 
-  USE module3D_phy, only: cde,abl,bdot,bmelt,epsxx,epsyy,epsxy,hwater,calv,hdot,hdotwater,pdd, &
-                    slope,taub,tg,ubx,uby,uzk,uxbar,uybar,vbar,ibase,tpmp,mk,mk0,h,grzone,front, &
-                    frontfacex,frontfacey,bsoc,flot,sealevel_2D,flotmx,flotmy,okumat,okvmat,gzmx, &
+  USE module3D_phy, only: cde,abl,bdot,bmelt,epsxx,epsyy,epsxy,hwater,calv,hdot,hdotwater, &
+                    slope,taub,ubx,uby,uzk,uxbar,uybar,ibase,tpmp,mk,mk0,h,front, &
+                    bsoc,flot,sealevel_2D,flotmx,flotmy,gzmx, &
                     gzmy,flgzmx,flgzmy,ilemx,ilemy,sdx,sdy,ux,uy,uzr,t,s
   use geography, only: nx,ny,nz,nzm,dx
-  USE param_phy_mod
+  USE param_phy_mod,only: ro,row
   implicit none
   real,dimension(NZ) ::cord_vert 
@@ -47,8 +47,6 @@
   HDOT(:,:)=0. 
   HDOTWATER(:,:)=0.
-  PDD(:,:)=0.
   SLOPE(:,:)=0. 
   TAUB(:,:)=0. 
-  TG(:,:)=0.
   UBX(:,:)=0.
   UBY(:,:)=0.
@@ -56,5 +54,4 @@
   UXBAR(:,:)=0.
   UYBAR(:,:)=0.
-  VBAR(:,:)=0.
   IBASE(:,:)=1
   TPMP(:,:,1)=0
@@ -62,11 +59,8 @@
 
   MK(:,:)=MK0(:,:)
-  grzone(:,:)=.false.
 
 
   !-------Initialisation des fronts.
   FRONT(:,:)  =0 
-  FRONTFACEX(:,:)=0
-  FRONTFACEY(:,:)=0
 
   !-------Initialisation des flot... devrait etre mise dans input_topo.
@@ -86,6 +80,4 @@
   FLOTMX(:,:)=FLOT(:,:)
   FLOTMY(:,:)=FLOT(:,:)
-  OKUMAT(:,:)=.FALSE.
-  OKVMAT(:,:)=.FALSE.
   GZMX(:,:)=.FALSE.
   GZMY(:,:)=.FALSE.
Index: /branches/GRISLIv3/SOURCES/isostasie_mod-0.3.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/isostasie_mod-0.3.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/isostasie_mod-0.3.f90	(revision 467)
@@ -16,6 +16,12 @@
 
 module isostasie_mod
- 
+
+   use geography, only: nx,ny 
    use iso_declar,only: nbed 
+
+   implicit none
+
+   real,dimension(nx,ny) :: W0           !< enfoncement du socle a l'equilibre isostatique
+   real,dimension(nx,ny) :: W1           !< enfoncement du socle courant
 
 contains
@@ -27,5 +33,5 @@
 
     use iso_declar,only: nlith,dt_iso,tausoc,dl,rl,lbloc,we,charge
-    use module3D_phy, only: icouple,marine,err,h0,bsoc0,sealevel_2d,w0,w1
+    use module3D_phy, only: icouple,marine,err,h0,bsoc0,sealevel_2d
     use geography, only: geoplace,nx,ny,dx,dy
     use param_phy_mod, only: ro,row,rog,rowg,romg
Index: /branches/GRISLIv3/SOURCES/lineartemp-0.2.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/lineartemp-0.2.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/lineartemp-0.2.f90	(revision 467)
@@ -13,5 +13,5 @@
 
   use geography, only:nx,ny,nz,nzm
-  use module3D_phy, only:H1,H,B1,B,TG,bmelt,T,TS,tpmp,ghf
+  use module3D_phy, only:H,B,bmelt,T,TS,tpmp,ghf
   use icetemp_declar, only:Dzm,Cm
   
@@ -19,9 +19,8 @@
 
   integer :: i,j,k
+  real,dimension(nx,ny) :: TG           !< degrees above melting point at the base
 
   do I=1,NX
      do J=1,NY 
-        H1(I,J)=H(I,J) 
-        B1(I,J)=B(I,J) 
         TG(I,J)=0.
         BMELT(I,J)=0.
Index: /branches/GRISLIv3/SOURCES/litho-0.4.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/litho-0.4.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/litho-0.4.f90	(revision 467)
@@ -39,6 +39,7 @@
 subroutine litho
   !$ USE OMP_LIB
-  USE module3D_phy, only: err,w1
+  USE module3D_phy, only: err
   use geography, only: nx,ny
+  use isostasie_mod, only: w1
   USE iso_declar, only: lbloc,we,charge ! module de declaration des variables specifiques a l'isostasie
 
Index: /branches/GRISLIv3/SOURCES/main3D-0.4-40km.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/main3D-0.4-40km.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/main3D-0.4-40km.f90	(revision 467)
@@ -133,5 +133,5 @@
   USE module3D_phy, only: itemp,icouple,isynchro,icompteur,iglen,marine,num_sealevel, &
                           num_ts_ritz,num_ic_vo,num_ic_by,num_ic_dm,num_ic_dc,num_ic_df, &
-                          s,h,b,bsoc,flot,mk,mk0,uxbar,uybar,hwater,time,timemax,boost,ndebug,ndebug_max
+                          s,h,b,bsoc,flot,mk,mk0,uxbar,uybar,hwater,time,timemax,ndebug,ndebug_max
   use runparam, only: nt,tbegin,tgrounded,dtprofile,dtcpt,dirnameout,runname,itracebug
   use geography, only: nx,ny,geoplace
@@ -349,6 +349,4 @@
   endif
 
-  boost = .false.
-
   do i=2,nx-1
      do j=2,ny-1
Index: /branches/GRISLIv3/SOURCES/prescribe-H-i2s_mod.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/prescribe-H-i2s_mod.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/prescribe-H-i2s_mod.f90	(revision 467)
@@ -21,4 +21,6 @@
   !  use toy_retreat_mod
 
+  use geography,only: nx,ny
+
   implicit none
   !  real,dimension(nx,ny)    :: Hp           !< H value if prescribed
@@ -28,6 +30,6 @@
   !  integer,dimension(nx,ny) :: i_Hp         !< 1 if H is prescribed on this node, else 0
   !  integer,dimension(nx,ny) :: i_Hp0        !< i_hp mask reference value does not change with time
-  !  integer,dimension(nx,ny) :: MK_gl0       !< mask grounding line initial
-  !  integer,dimension(nx,ny) :: MK_flot0     !< mask float initial
+    integer,dimension(nx,ny) :: MK_gl0       !< mask grounding line initial
+    integer,dimension(nx,ny) :: MK_flot0     !< mask float initial
 
   ! pour grounding line retreat, ice2sea
@@ -48,5 +50,5 @@
   subroutine init_prescribe_H
 
-    use module3D_phy, only: flot,MK_flot0,MK_gl0,Hp0,H0,Mk_init,i_hp0
+    use module3D_phy, only: flot,Hp0,H0,Mk_init,i_hp0
     use geography, only: nx,ny,geoplace
     use runparam, only: itracebug 
@@ -103,5 +105,5 @@
   subroutine prescribe_present_H_gl
     
-    use module3D_phy, only: MK_flot0,MK_gl0,i_hp,Hp,Hp0
+    use module3D_phy, only: i_hp,Hp,Hp0
     use runparam, only: itracebug 
     
@@ -170,5 +172,5 @@
   subroutine prescribe_paleo_gl_shelf
     
-    use module3D_phy, only: MK_flot0,MK_gl0,i_hp,hp,sealevel_2d,Bsoc
+    use module3D_phy, only: i_hp,hp,sealevel_2d,Bsoc
     use runparam, only: itracebug 
     use param_phy_mod, only: row,ro
@@ -245,5 +247,5 @@
   subroutine break_all_ice_shelves
     
-    use module3D_phy, only: i_hp,hp,hp0,H,debug_3D,flot,MK_flot0
+    use module3D_phy, only: i_hp,hp,hp0,H,debug_3D,flot
     use runparam, only: itracebug
     
@@ -336,5 +338,5 @@
   subroutine prescribe_present_H_gl_copy
     
-    use module3D_phy, only: i_hp,hp,Hp0,MK_flot0,MK_gl0
+    use module3D_phy, only: i_hp,hp,Hp0
     use runparam, only: itracebug 
     use geography, only: nx,ny
Index: /branches/GRISLIv3/SOURCES/slope_surf.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/slope_surf.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/slope_surf.f90	(revision 467)
@@ -13,9 +13,12 @@
 subroutine slope_surf
 
-  use module3D_phy, only: sdx,sdy,S,slope,sdxmy,sdymx,slope2mx,slope2my,debug_3D
+  use module3D_phy, only: sdx,sdy,S,slope,slope2mx,slope2my,debug_3D
   use geography, only: nx,ny,dx,dy
   !$ USE OMP_LIB
   
   implicit none
+
+  real,dimension(nx,ny) :: SDXMY        !< slope selon x moy selon y '^' remplace SDMX
+  real,dimension(nx,ny) :: SDYMX        !< slope selon y moy selon x '>' remplace SDMY
 
   real :: inv_4dx        ! inverse de dx pour eviter les divisions =1/(4*dx)
Index: /branches/GRISLIv3/SOURCES/spinup_mod.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/spinup_mod.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/spinup_mod.f90	(revision 467)
@@ -682,8 +682,9 @@
   subroutine ajust_ghf
     
-    use module3D_phy, only: debug_3D,flot,ibase,ghf,ghf0,secyear
+    use module3D_phy, only: debug_3D,flot,ibase,ghf,secyear
     
     implicit none
-    
+
+    real,dimension(nx,ny) :: ghf0         !< geothermal heat flux J/m2/a 'o'
     real,dimension(nx,ny) :: coefdef_maj     !< coefficient deformation
     real :: increment_ghf
Index: /branches/GRISLIv3/SOURCES/taubed-0.3.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/taubed-0.3.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/taubed-0.3.f90	(revision 467)
@@ -44,7 +44,8 @@
 
   !$USE OMP_LIB 
-  USE module3D_phy, only: h,bsoc0,bsoc,bdot,sealevel_2d,w0,w1
+  USE module3D_phy, only: h,bsoc0,bsoc,bdot,sealevel_2d
   USE geography, only: nx,ny
   USE param_phy_mod, only: ro,row,rog,rowg,rom
+  USE isostasie_mod, only:w0,w1
   USE iso_declar,only: nlith,lbloc,charge,dt_iso,tausoc ! module de declaration des variables de l'isostasie
 
Index: /branches/GRISLIv3/SOURCES/velocities-polyn-0.3.f90
===================================================================
--- /branches/GRISLIv3/SOURCES/velocities-polyn-0.3.f90	(revision 466)
+++ /branches/GRISLIv3/SOURCES/velocities-polyn-0.3.f90	(revision 467)
@@ -26,5 +26,5 @@
   use runparam, only : itracebug      
   use module3d_phy, only: flotmx,flotmy,ux,uy,sux,suy,sdx,sdy,ubx,uby,&
-       hmx,hmy,uxbar,uybar,iglen,cde,divu,uzr,bmelt,bm,flot,front,flgzmx,flgzmy,&
+       hmx,hmy,uxbar,uybar,iglen,cde,uzr,bmelt,bm,flot,front,flgzmx,flgzmy,&
        uzsdot,dtt,uzk,bdot,hdot
   use geography, only : nx,ny,nz,dx
@@ -36,4 +36,6 @@
   real, dimension(nx,ny) :: hdd
   real, dimension(nx,ny) :: xx
+  real, dimension(nx,ny) :: divu         !< divergence de UH
+
 
   if (itracebug.eq.1)  call tracebug(' Entree dans routine SIA_velocity')
