  SUBROUTINE dyn2dx_tls_solve(xeuv)

!&E---------------------------------------------------------------------
!&E                 ***  ROUTINE dyn2dx_tls_solve  ***
!&E
!&E ** Purpose :
!&E
!&E ** Description :
!&E
!&E ** Called by :
!&E
!&E ** External calls :
!&E
!&E ** Used ij-arrays : y2dxe,y2dv,axe,bxe,cxe,cv,bv,av
!&E
!&E ** Modified variables : xeuv,uv (1D)
!&E
!&E ** Reference :
!&E                T. Duhaut, M. Honnorat, L. Debreu, 2008 : Rapport
!&E                contrat PREVIMER - Ref : 06/2 210 290 (pages 100-125)
!&E
!&E ** History :
!&E       !  2010-07  (L. Debreu, F. Dumas) rewritting of ADI
!&E       !  2011-11  (V. Garnier) Spatial extension of ssh(liminm2:limaxp2,ljminm2:ljmaxp2)
!&E
!&E---------------------------------------------------------------------





























































!! * Modules used
   USE parameters
   USE comvars2d,   ONLY : ig,id

   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
   REAL(kind=rlg)                                    :: xk
   REAL(kind=rlg),DIMENSION(imax-imin-1)             :: p0,q0,u0

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


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


! nothing to do if no wet point all along the line
     IF (ig(j) < id(j)) THEN

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

! factorisation du systeme lineaire tridiagonal en LU
!----------------------------------------------------

       q0(1)=-cxe(ig(j)+1,j)/bxe(ig(j)+1,j)
       DO k=ig(j)+2,id(j)-1
!         write(111,*)'k : ',k,'  j : ',j,'  ig : ',ig(j)
!         write(222,*)'axe : ',axe(k,j),'  q0(k-ig(j)-1) : ',q0(k-ig(j)-1),'  bxe : ',bxe(k,j)
         p0(k-ig(j)  )= 1.0_rlg/(axe(k,j)*q0(k-ig(j)-1)+bxe(k,j))
         q0(k-ig(j)  )=-cxe(k,j)*p0(k-ig(j)  )
       END DO


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

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

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

       xk=0.0_rlg
       k=n+1

       DO i=id(j)-1,ig(j)+1,-1
         k=k-1
         xk=q0(k)*xk+u0(k)
         xeuv(i,j)=xk
       END DO

     END IF

   END DO
!$OMP END DO

   !PRINT*, 'exit DYN2DX_TLS_SOLVE routine'

   RETURN

  END SUBROUTINE dyn2dx_tls_solve
