  SUBROUTINE dyn2dx_cyc_solve(xeuv)

!&E---------------------------------------------------------------------
!&E                 ***  ROUTINE dyn2dx_cyc_solve  ***
!&E
!&E ** Purpose :
!&E
!&E ** Description :
!&E
!&E ** Called by :
!&E
!&E ** External calls :
!&E
!&E ** Used ij-arrays :
!&E
!&E ** Modified variables : xeuv (1D)
!&E
!&E ** Reference :
!&E
!&E ** History :
!&E       !  2010-07  (L. Debreu) rewritting of ADI - substitution of Lapack routine
!&E       !  2011-11  (V. Garnier) Spatial extension of ssh(liminm2:limaxp2,ljminm2:ljmaxp2)
!&E
!&E---------------------------------------------------------------------
!! * Modules used
   USE parameters
   USE comvars2d,   ONLY : id,ig

   USE toolmpi,     ONLY : rjminp2,rjmaxm1
   USE comcolig_mpi,ONLY : y2dxe=>ry2dxe,axe=>raxe,bxe=>rbxe,cxe=>rcxe


   IMPLICIT NONE

!! * Arguments

   REAL(kind=rsh),DIMENSION(imin:imax,rjminp2:rjmaxm1),INTENT(inout)       :: xeuv


!! * Local declarations
   INTEGER                                                   :: i,j,k,n,ktemp
   REAL(kind=rlg)                                            :: xk,p0
   REAL(kind=rlg),DIMENSION(2*(imax-imin)-3)                 :: q0,s0,u0

!!----------------------------------------------------------------------
!! * Executable part


!$OMP DO SCHEDULE(RUNTIME)
   DO j=rjminp2,rjmaxm1


     
       n=id(j)-ig(j)-1

! introduction des conditions aux limites
!----------------------------------------
! factorisation du systeme lineaire tridiagonal en LU
!----------------------------------------------------

! resolution du systeme triangulaire inferieur
!---------------------------------------------


       k = ig(j)+1
       
       q0(1)=-cxe(k,j)/bxe(k,j)
       s0(1)=-axe(k,j)/bxe(k,j)
       u0(1)=y2dxe(k,j)/bxe(k,j)
       
       ktemp = 1
       
       DO k=ig(j)+2,id(j)-1
         ktemp = ktemp + 1
         p0= 1./(axe(k,j)*q0(ktemp-1)+bxe(k,j))
         q0(ktemp)=-cxe(k,j)*p0
         s0(ktemp)=-axe(k,j)*s0(ktemp - 1 )*p0
         u0(ktemp )=(y2dxe(k,j)-axe(k,j)*u0(ktemp-1))*p0
       END DO



! resolution du systeme triangulaire superieur
!---------------------------------------------

       k=n
       q0(k)=0.0_rlg
       s0(k)=1.0_rlg
       
       DO i=id(j)-2,ig(j)+1,-1
         k=k-1
         s0(k) = s0(k)  + q0(k)*s0(k+1)
         q0(k) = u0(k)  + q0(k)*q0(k+1)
       END DO

       k=id(j)-1
       
       xk=(y2dxe(k,j)-cxe(k,j)*q0(1)-axe(k,j)*q0(n-1))/(cxe(k,j)*s0(1)+axe(k,j)*s0(n-1)+bxe(k,j))
       
       xeuv(k,j)=xk
       
       k = 0
       
       DO i=ig(j)+1,id(j)-2
         k=k+1
         xeuv(i,j)=xk*s0(k)+q0(k)
       END DO

   END DO
!$OMP END DO

   RETURN

  END SUBROUTINE dyn2dx_cyc_solve
