anchor.diff:                                -Ro*(temp2*temp4 - temp1*temp5)
   
                       ! v (j=1)      
  -
                       temp3 = (4*leta(i,j+1) - leta(i,j+2) - 3*leta(i,j))/2/dy 
                       temp5 = (eta(i,j+2) + eta(i,j) - 2*eta(i,j+1))/(dy*dy)
                       v2(i,j) = temp2-Et*(temp1 - temp3)/dt             &
                                - Ro*(temp2*temp5 - temp1*temp4)
                       
  -           
                       ! u (j=pasy)            
                       j = pasy
  -
                       temp1 = (-4*eta(i,j-1) + eta(i,j-2) + 3*eta(i,j))/2/dy
                       temp2 = (eta(i+1,j) - eta(i-1,j))/2/dx
                       temp3 = (leta(i+1,j) - leta(i-1,j))/2/dx
  
positive.code: 
          do i = 2, pasx
              do j = 2, pasy-1
                  temp1 = ((eta(i+1,j)-eta(i-1,j))*(b(i,j+1)-b(i,j-1))-          &
                  (eta(i,j+1) - eta(i,j-1))*(b(i+1,j) - b(i-1,j)))/4/dx/dy
                  temp2 = (eta(i+1,j)*(b(i+1,j+1) - b(i+1,j-1)) - eta(i-1,j)*        &
                      (b(i-1,j+1) - b(i-1,j-1)) - eta(i,j+1)*(b(i+1,j+1)             &
                              -b(i-1,j+1)) + eta(i,j-1)*(b(i+1,j-1) - b(i-1,j-1)))   &
                                      /4/dx/dy
                  temp3 = (b(i,j+1)*(eta(i+1,j+1) - eta(i-1,j+1)) - b(i,j-1)*        &
                      (eta(i+1,j-1) - eta(i-1,j-1)) - b(i+1,j)*(eta(i+1,j+1)-        &
                      eta(i+1,j-1)) + b(i-1,j)*(eta(i-1,j+1) - eta(i-1,j-1)))        &
                      /4/dx/dy
                  jacob(i,j) = (temp1 + temp2 + temp3)/3.0D0
  
negative.code: end module GEOPHY    
  
  
  ! ------------------------------------------------------------------
  subroutine init_xy(x, y, dx, dy, Ly)
      use GEOPHY
      implicit none
  
      real(DP) :: x(pasx+1), y(pasx+1)   !! [RB] taille exacte necessaire sinon ca plante au runtime
      real(DP) :: dx, dy, Ly
  
      integer :: i, j
  
      do i = 2, pasx+1
          print *,"i=",i
          x(i) = (i-2)*dx
      end do
  
      x(1) = x(pasx)
      x(2) = x(pasx+1)
  
      do j = 1, pasy      
          y(j) = -Ly + (j-1)*dy
      end do
  end subroutine init_xy
  
  ! ------------------------------------------------------------------
  
  program NIHOUL
  
  