Index: branches/GRISLIv3/SOURCES/deformation_mod_2lois_isotherme.f90
===================================================================
--- branches/GRISLIv3/SOURCES/deformation_mod_2lois_isotherme.f90	(revision 390)
+++ branches/GRISLIv3/SOURCES/deformation_mod_2lois_isotherme.f90	(revision 391)
@@ -27,42 +27,44 @@
 module deformation_mod_2lois_isoth
 
-use module3d_phy
-implicit none
+  use module3d_phy, only:nx,ny,nz,e,num_param,num_rep_42,secyear,iglen
+  use param_phy_mod, only: rog,ro,row
+  
+  implicit none
 
 
-! declarations ne dependant pas du nombre de lois
+  ! declarations ne dependant pas du nombre de lois
 
-real, dimension(nx,ny,nz)   :: visc           !< viscosite de la glace (pour les shelves)
-real, dimension(nz)         :: edecal         !< travail (decalage de E de 1 indice)
+  real, dimension(nx,ny,nz)   :: visc           !< viscosite de la glace (pour les shelves)
+  real, dimension(nz)         :: edecal         !< travail (decalage de E de 1 indice)
 
-! declaration des lois
-! la valeur de n1poly et de n2poly determine le nombre de lois additionnees
+  ! declaration des lois
+  ! la valeur de n1poly et de n2poly determine le nombre de lois additionnees
 
-integer, parameter :: n1poly=1     !< 2 lois numerotees de 1 a 2
-integer, parameter :: n2poly=2
+  integer, parameter :: n1poly=1     !< 2 lois numerotees de 1 a 2
+  integer, parameter :: n2poly=2
 
 
-! Tous les parametres de la loi sous forme de tableau (n1poly:n2poly)
+  ! Tous les parametres de la loi sous forme de tableau (n1poly:n2poly)
 
-real, dimension(n1poly:n2poly) :: glen             !< l'exposant
-real, dimension(n1poly:n2poly) :: sf               !< softening factor for flow law
-real, dimension(n1poly:n2poly) :: B_isoth          !< coefficient
+  real, dimension(n1poly:n2poly) :: glen             !< l'exposant
+  real, dimension(n1poly:n2poly) :: sf               !< softening factor for flow law
+  real, dimension(n1poly:n2poly) :: B_isoth          !< coefficient
 
 
 
 
-! declaration des tableaux utilises dans le calcul des vitesses, viscosites, ...
+  ! declaration des tableaux utilises dans le calcul des vitesses, viscosites, ...
 
-real, dimension(nx,ny,nz,n1poly:n2poly) :: Btt     !< coefficient au point considere
-real, dimension(nx,ny,nz,n1poly:n2poly) :: Sa      !< sert dans l'integration des vitesses
-real, dimension(nx,ny,nz,n1poly:n2poly) :: Sa_mx   !< Sa sur demi mailles >
-real, dimension(nx,ny,nz,n1poly:n2poly) :: Sa_my   !< Sa sur demi mailles ^
-real, dimension(nx,ny,nz,n1poly:n2poly) :: S2a     !< sert dans l'integration des vitesses
-real, dimension(nx,ny,nz,n1poly:n2poly) :: S2a_mx  !< S2a sur demi mailles >
-real, dimension(nx,ny,nz,n1poly:n2poly) :: S2a_my  !< S2a sur demi mailles ^
-real, dimension(nx,ny,n1poly:n2poly)    :: ddx     !< sert dans le calcul de conserv. masse
-real, dimension(nx,ny,n1poly:n2poly)    :: ddy     !< sert dans le calcul de conserv. masse
+  real, dimension(nx,ny,nz,n1poly:n2poly) :: Btt     !< coefficient au point considere
+  real, dimension(nx,ny,nz,n1poly:n2poly) :: Sa      !< sert dans l'integration des vitesses
+  real, dimension(nx,ny,nz,n1poly:n2poly) :: Sa_mx   !< Sa sur demi mailles >
+  real, dimension(nx,ny,nz,n1poly:n2poly) :: Sa_my   !< Sa sur demi mailles ^
+  real, dimension(nx,ny,nz,n1poly:n2poly) :: S2a     !< sert dans l'integration des vitesses
+  real, dimension(nx,ny,nz,n1poly:n2poly) :: S2a_mx  !< S2a sur demi mailles >
+  real, dimension(nx,ny,nz,n1poly:n2poly) :: S2a_my  !< S2a sur demi mailles ^
+  real, dimension(nx,ny,n1poly:n2poly)    :: ddx     !< sert dans le calcul de conserv. masse
+  real, dimension(nx,ny,n1poly:n2poly)    :: ddy     !< sert dans le calcul de conserv. masse
 
-real                                    :: de= 1./(nz-1)
+  real                                    :: de= 1./(nz-1)
 
 
@@ -71,114 +73,116 @@
 contains
 
-!-------------------------------------------------------------------------------------------
+  !-------------------------------------------------------------------------------------------
 
-!> SUBROUTINE: init_deformation
-!! Routine qui lit les variables de deformation 
-!>
-subroutine init_deformation 
+  !> SUBROUTINE: init_deformation
+  !! Routine qui lit les variables de deformation 
+  !>
+  subroutine init_deformation 
 
-real :: exposant_1
-real :: enhanc_fact_1
-real :: coef_iso_1
+    integer :: k
+    real :: exposant_1
+    real :: enhanc_fact_1
+    real :: coef_iso_1
 
-real :: exposant_2
-real :: enhanc_fact_2
-real :: coef_iso_2
+    real :: exposant_2
+    real :: enhanc_fact_2
+    real :: coef_iso_2
 
 
-namelist/loidef_1_iso/exposant_1,enhanc_fact_1,coef_iso_1
-namelist/loidef_2_iso/exposant_2,enhanc_fact_2,coef_iso_2
+    namelist/loidef_1_iso/exposant_1,enhanc_fact_1,coef_iso_1
+    namelist/loidef_2_iso/exposant_2,enhanc_fact_2,coef_iso_2
 
 
-! loi 1
-!--------------------------------------------------------------------------------
+    ! loi 1
+    !--------------------------------------------------------------------------------
 
-rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
-read(num_param,loidef_1_iso)
+    rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
+    read(num_param,loidef_1_iso)
 
-write(num_rep_42,*)'!___________________________________________________________' 
-write(num_rep_42,*)'! loi de deformation 1    module deformation_mod_2lois_isoth'
-write(num_rep_42,*)
-write(num_rep_42,loidef_1_iso)
-write(num_rep_42,*)'! exposant (glen), coefficient et enhancement factor (sf)   '
-write(num_rep_42,*)'!___________________________________________________________' 
+    write(num_rep_42,*)'!___________________________________________________________' 
+    write(num_rep_42,*)'! loi de deformation 1    module deformation_mod_2lois_isoth'
+    write(num_rep_42,*)
+    write(num_rep_42,loidef_1_iso)
+    write(num_rep_42,*)'! exposant (glen), coefficient et enhancement factor (sf)   '
+    write(num_rep_42,*)'!___________________________________________________________' 
 
-glen(1)    = exposant_1 
-sf(1)      = enhanc_fact_1  
-B_isoth(1) = coef_iso_1*2.*secyear
-write(6,*) 'loi de deformation B=',B_isoth(1)
-write(6,*) 'gamma=', rog*(1.-ro/row)
-write(6,*) 'epsilon = coef1 * H^3,   coef1 =',B_isoth(1)/2.*(rog*(1.-ro/row)/4.)**3
-write(6,*) 'pvi = coef2 /H,          coef2 =',(4./(rog*(1.-ro/row)))**2 / B_isoth(1)
+    glen(1)    = exposant_1 
+    sf(1)      = enhanc_fact_1  
+    B_isoth(1) = coef_iso_1*2.*secyear
+    write(6,*) 'loi de deformation B=',B_isoth(1)
+    write(6,*) 'gamma=', rog*(1.-ro/row)
+    write(6,*) 'epsilon = coef1 * H^3,   coef1 =',B_isoth(1)/2.*(rog*(1.-ro/row)/4.)**3
+    write(6,*) 'pvi = coef2 /H,          coef2 =',(4./(rog*(1.-ro/row)))**2 / B_isoth(1)
 
-! loi 2
-!--------------------------------------------------------------------------------
-  
-rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
-read(num_param,loidef_2_iso)
+    ! loi 2
+    !--------------------------------------------------------------------------------
 
-write(num_rep_42,*)'!___________________________________________________________' 
-write(num_rep_42,*)'! loi de deformation 2    module deformation_mod_2lois_isoth'
-write(num_rep_42,*)
-write(num_rep_42,loidef_2_iso)
-write(num_rep_42,*)'! exposant (glen), coefficient et enhancement factor (sf)   '
-write(num_rep_42,*)'!___________________________________________________________' 
+    rewind(num_param)        ! pour revenir au debut du fichier param_list.dat
+    read(num_param,loidef_2_iso)
 
-glen(2)    = exposant_2 
-sf(2)      = enhanc_fact_2  
-B_isoth(2) = coef_iso_2
+    write(num_rep_42,*)'!___________________________________________________________' 
+    write(num_rep_42,*)'! loi de deformation 2    module deformation_mod_2lois_isoth'
+    write(num_rep_42,*)
+    write(num_rep_42,loidef_2_iso)
+    write(num_rep_42,*)'! exposant (glen), coefficient et enhancement factor (sf)   '
+    write(num_rep_42,*)'!___________________________________________________________' 
 
-! calcul coordonnées reduites
+    glen(2)    = exposant_2 
+    sf(2)      = enhanc_fact_2  
+    B_isoth(2) = coef_iso_2
 
-e(1)=0.
-e(nz)=1.
-do k=1,nz
-   if ((k.ne.1).and.(k.ne.nz)) e(k)=(k-1.)/(nz-1.)
-end do
+    ! calcul coordonnées reduites
 
-do iglen= n1poly,n2poly
-   call flowlaw(iglen)
-end do
+    e(1)=0.
+    e(nz)=1.
+    do k=1,nz
+       if ((k.ne.1).and.(k.ne.nz)) e(k)=(k-1.)/(nz-1.)
+    end do
 
-ddx(:,:,:)=0.
-ddy(:,:,:)=0.
+    do iglen= n1poly,n2poly
+       call flowlaw(iglen)
+    end do
+
+    ddx(:,:,:)=0.
+    ddy(:,:,:)=0.
 
 
-end subroutine init_deformation
+  end subroutine init_deformation
 
 
-!--------------------------------------------------------------------------------
-subroutine flow_general  ! fake
-end subroutine flow_general
+  !--------------------------------------------------------------------------------
+  subroutine flow_general  ! fake
+  end subroutine flow_general
 
-!--------------------------------------------------------------------------------
-! application des sf et calcul des btt, Sa, ...
-! Tous independants de x et y
+  !--------------------------------------------------------------------------------
+  ! application des sf et calcul des btt, Sa, ...
+  ! Tous independants de x et y
 
-subroutine flowlaw (iiglen)
+  subroutine flowlaw (iiglen)
 
-  implicit none
-  integer ::  iiglen   !< compteur pour la boucle flowlaw
-  real    :: n1        !< glen+1
-  real    :: n2        !< glen+2
+    implicit none
+    integer :: k
+    integer ::  iiglen   !< compteur pour la boucle flowlaw
+    real    :: n1        !< glen+1
+    real    :: n2        !< glen+2
 
 
 
-  B_isoth(iiglen)    =  B_isoth(iiglen)*sf(iiglen)
-  Btt(:,:,:,iiglen)  =  B_isoth(iiglen)
-  n1                 =  glen(iiglen)+1.
-  n2                 =  glen(iiglen)+2.
- 
+    B_isoth(iiglen)    =  B_isoth(iiglen)*sf(iiglen)
+    Btt(:,:,:,iiglen)  =  B_isoth(iiglen)
+    n1                 =  glen(iiglen)+1.
+    n2                 =  glen(iiglen)+2.
 
-  do k = 1,nz
-     Sa (:,:,k,iiglen)    =  B_isoth(iiglen)/n1 * (1.-e(k)**n1 )
-     S2a(:,:,k,iiglen)    =  B_isoth(iiglen)/n2 * (1+ (e(k)**n2-n2*e(k))/n1 ) 
-  end do
-  Sa_mx(:,:,:,iiglen)  =  Sa (:,:,:,iiglen)
-  Sa_my(:,:,:,iiglen)  =  Sa (:,:,:,iiglen)
-  S2a_mx(:,:,:,iiglen) =  S2a (:,:,:,iiglen)
-  S2a_my(:,:,:,iiglen) =  S2a (:,:,:,iiglen)
 
-end subroutine flowlaw
+    do k = 1,nz
+       Sa (:,:,k,iiglen)    =  B_isoth(iiglen)/n1 * (1.-e(k)**n1 )
+       S2a(:,:,k,iiglen)    =  B_isoth(iiglen)/n2 * (1+ (e(k)**n2-n2*e(k))/n1 ) 
+    end do
+    Sa_mx(:,:,:,iiglen)  =  Sa (:,:,:,iiglen)
+    Sa_my(:,:,:,iiglen)  =  Sa (:,:,:,iiglen)
+    S2a_mx(:,:,:,iiglen) =  S2a (:,:,:,iiglen)
+    S2a_my(:,:,:,iiglen) =  S2a (:,:,:,iiglen)
+
+  end subroutine flowlaw
 
 end module deformation_mod_2lois_isoth
