anchor.diff:  
           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))-          &
  +                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
  +                    -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) - 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
  
positive.code:                                - Ro*(temp2*temp4 - temp1*temp5)
  
                      ! v (j=pasy)     
  
                      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         &
  
negative.code:                 do i = 2, pasx
                      do j = start+1, pasy-1, 2
                          temp1 = -2/dx**2 - 2/dy**2 - 1/Bu
                          temp1 = (q(i,j)/Bu - (eta(i+1,j) + eta(i-1,j))/dx**2   &
                                             - (eta(i,j+1) + eta(i,j-1))/dy**2  ) / temp1
                          temp1 = eta(i,j) + w*(temp1-eta(i,j))
                          if(abs(eta(i,j)) > 0.0D0) then
                              erreur = abs( (temp1-eta(i,j))/eta(i,j) )
  