[uplo,trans,diag,n,a,lda,x,incx]=ctrsv(uplo,trans,diag,n,a,lda,x,incx);

function [uplo,trans,diag,n,a,lda,x,incx]=ctrsv(uplo,trans,diag,n,a,lda,x,incx);
%***BEGIN PROLOGUE  CTRSV
%***PURPOSE  Solve a complex triangular system of equations.
%***LIBRARY   SLATEC (BLAS)
%***CATEGORY  D1B4
%***TYPE      COMPLEX (STRSV-S, DTRSV-D, CTRSV-C)
%***KEYWORDS  LEVEL 2 BLAS, LINEAR ALGEBRA
%***AUTHOR  Dongarra, J. J., (ANL)
%           Du Croz, J., (NAG)
%           Hammarling, S., (NAG)
%           Hanson, R. J., (SNLA)
%***DESCRIPTION
%
%  CTRSV  solves one of the systems of equations
%
%     A*x = b,   or   A'*x = b,   or   conjg( A')*x = b,
%
%  where b and x are n element vectors and A is an n by n unit, or
%  non-unit, upper or lower triangular matrix.
%
%  No test for singularity or near-singularity is included in this
%  routine. Such tests must be performed before calling this routine.
%
%  Parameters
%  ==========
%
%  UPLO   - CHARACTER*1.
%           On entry, UPLO specifies whether the matrix is an upper or
%           lower triangular matrix as follows:
%
%              UPLO = 'U' or 'u'   A is an upper triangular matrix.
%
%              UPLO = 'L' or 'l'   A is a lower triangular matrix.
%
%           Unchanged on exit.
%
%  TRANS  - CHARACTER*1.
%           On entry, TRANS specifies the equations to be solved as
%           follows:
%
%              TRANS = 'N' or 'n'   A*x = b.
%
%              TRANS = 'T' or 't'   A'*x = b.
%
%              TRANS = 'C' or 'c'   conjg( A' )*x = b.
%
%           Unchanged on exit.
%
%  DIAG   - CHARACTER*1.
%           On entry, DIAG specifies whether or not A is unit
%           triangular as follows:
%
%              DIAG = 'U' or 'u'   A is assumed to be unit triangular.
%
%              DIAG = 'N' or 'n'   A is not assumed to be unit
%                                  triangular.
%
%           Unchanged on exit.
%
%  N      - INTEGER.
%           On entry, N specifies the order of the matrix A.
%           N must be at least zero.
%           Unchanged on exit.
%
%  A      - COMPLEX          array of DIMENSION ( LDA, n ).
%           Before entry with  UPLO = 'U' or 'u', the leading n by n
%           upper triangular part of the array A must contain the upper
%           triangular matrix and the strictly lower triangular part of
%           A is not referenced.
%           Before entry with UPLO = 'L' or 'l', the leading n by n
%           lower triangular part of the array A must contain the lower
%           triangular matrix and the strictly upper triangular part of
%           A is not referenced.
%           Note that when  DIAG = 'U' or 'u', the diagonal elements of
%           A are not referenced either, but are assumed to be unity.
%           Unchanged on exit.
%
%  LDA    - INTEGER.
%           On entry, LDA specifies the first dimension of A as declared
%           in the calling (sub) program. LDA must be at least
%           max( 1, n ).
%           Unchanged on exit.
%
%  X      - COMPLEX          array of dimension at least
%           ( 1 + ( n - 1 )*abs( INCX ) ).
%           Before entry, the incremented array X must contain the n
%           element right-hand side vector b. On exit, X is overwritten
%           with the solution vector x.
%
%  INCX   - INTEGER.
%           On entry, INCX specifies the increment for the elements of
%           X. INCX must not be zero.
%           Unchanged on exit.
%
%***REFERENCES  Dongarra, J. J., Du Croz, J., Hammarling, S., and
%                 Hanson, R. J.  An extended set of Fortran basic linear
%                 algebra subprograms.  ACM TOMS, Vol. 14, No. 1,
%                 pp. 1-17, March 1988.
%***ROUTINES CALLED  LSAME, XERBLA
%***REVISION HISTORY  (YYMMDD)
%   861022  DATE WRITTEN
%   910605  Modified to meet SLATEC prologue standards.  Only comment
%           lines were modified.  (BKS)
%***end PROLOGUE  CTRSV
%     .. Scalar Arguments ..
%     .. Array Arguments ..
persistent i info ix j jx kx noconj nounit temp zero ; 

a_shape=size(a);a=reshape([a(:).',zeros(1,ceil(numel(a)./prod([lda])).*prod([lda])-numel(a))],lda,[]);
x_shape=size(x);x=reshape(x,1,[]);
%     .. Parameters ..
if isempty(zero), zero=complex(0.0e+0,0.0e+0) ; end;
%     .. Local Scalars ..
if isempty(temp), temp=0; end;
if isempty(i), i=0; end;
if isempty(info), info=0; end;
if isempty(ix), ix=0; end;
if isempty(j), j=0; end;
if isempty(jx), jx=0; end;
if isempty(kx), kx=0; end;
if isempty(noconj), noconj=false; end;
if isempty(nounit), nounit=false; end;
%     .. External Functions ..
%     .. External Subroutines ..
%     .. Intrinsic Functions ..
%***FIRST EXECUTABLE STATEMENT  CTRSV
%
%     Test the input parameters.
%
info = 0;
if( ~lsame(uplo,'U') && ~lsame(uplo,'L') )
info = 1;
elseif ( ~lsame(trans,'N') && ~lsame(trans,'T') &&~lsame(trans,'C') ) ;
info = 2;
elseif ( ~lsame(diag,'U') && ~lsame(diag,'N') ) ;
info = 3;
elseif( n<0 ) ;
info = 4;
elseif( lda<max(1,n) ) ;
info = 6;
elseif( incx==0 ) ;
info = 8;
end;
if( info~=0 )
[dumvar1,info]=xerbla('CTRSV ',info);
a_shape=zeros(a_shape);a_shape(:)=a(1:numel(a_shape));a=a_shape;
x_shape=zeros(x_shape);x_shape(:)=x(1:numel(x_shape));x=x_shape;
return;
end;
%
%     Quick return if possible.
%
if( n==0 )
a_shape=zeros(a_shape);a_shape(:)=a(1:numel(a_shape));a=a_shape;
x_shape=zeros(x_shape);x_shape(:)=x(1:numel(x_shape));x=x_shape;
return;
end;
%
[noconj ,trans]=lsame(trans,'T');
[nounit ,diag]=lsame(diag,'N');
%
%     Set up the start point in X if the increment is not unity. This
%     will be  ( N - 1 )*INCX  too small for descending loops.
%
if( incx<=0 )
kx = fix(1 -(n-1).*incx);
elseif( incx~=1 ) ;
kx = 1;
end;
%
%     Start the operations. In this version the elements of A are
%     accessed sequentially with one pass through A.
%
if( lsame(trans,'N') )
%
%        Form  x := inv( A )*x.
%
if( lsame(uplo,'U') )
if( incx==1 )
for j = n : -1: 1 ;
if( x(j)~=zero )
if( nounit )
x(j) = x(j)./a(j,j);
end;
temp = x(j);
for i = j - 1 : -1: 1 ;
x(i) = x(i) - temp.*a(i,j);
end; i = fix(1 -1);
end;
end; j = fix(1 -1);
else;
jx = fix(kx +(n-1).*incx);
for j = n : -1: 1 ;
if( x(jx)~=zero )
if( nounit )
x(jx) = x(jx)./a(j,j);
end;
temp = x(jx);
ix = fix(jx);
for i = j - 1 : -1: 1 ;
ix = fix(ix - incx);
x(ix) = x(ix) - temp.*a(i,j);
end; i = fix(1 -1);
end;
jx = fix(jx - incx);
end; j = fix(1 -1);
end;
elseif( incx==1 ) ;
for j = 1 : n;
if( x(j)~=zero )
if( nounit )
x(j) = x(j)./a(j,j);
end;
temp = x(j);
for i = j + 1 : n;
x(i) = x(i) - temp.*a(i,j);
end; i = fix(n+1);
end;
end; j = fix(n+1);
else;
jx = fix(kx);
for j = 1 : n;
if( x(jx)~=zero )
if( nounit )
x(jx) = x(jx)./a(j,j);
end;
temp = x(jx);
ix = fix(jx);
for i = j + 1 : n;
ix = fix(ix + incx);
x(ix) = x(ix) - temp.*a(i,j);
end; i = fix(n+1);
end;
jx = fix(jx + incx);
end; j = fix(n+1);
end;
%
%        Form  x := inv( A' )*x  or  x := inv( conjg( A' ) )*x.
%
elseif ( lsame(uplo,'U') ) ;
if( incx==1 )
for j = 1 : n;
temp = x(j);
if( noconj )
for i = 1 : j - 1;
temp = temp - a(i,j).*x(i);
end; i = fix(j - 1+1);
if( nounit )
temp = temp./a(j,j);
end;
else;
for i = 1 : j - 1;
temp = temp - conj(a(i,j)).*x(i);
end; i = fix(j - 1+1);
if( nounit )
temp = temp./conj(a(j,j));
end;
end;
x(j) = temp;
end; j = fix(n+1);
else;
jx = fix(kx);
for j = 1 : n;
ix = fix(kx);
temp = x(jx);
if( noconj )
for i = 1 : j - 1;
temp = temp - a(i,j).*x(ix);
ix = fix(ix + incx);
end; i = fix(j - 1+1);
if( nounit )
temp = temp./a(j,j);
end;
else;
for i = 1 : j - 1;
temp = temp - conj(a(i,j)).*x(ix);
ix = fix(ix + incx);
end; i = fix(j - 1+1);
if( nounit )
temp = temp./conj(a(j,j));
end;
end;
x(jx) = temp;
jx = fix(jx + incx);
end; j = fix(n+1);
end;
elseif( incx==1 ) ;
for j = n : -1: 1 ;
temp = x(j);
if( noconj )
for i = n : -1: j + 1 ;
temp = temp - a(i,j).*x(i);
end; i = fix(j + 1 -1);
if( nounit )
temp = temp./a(j,j);
end;
else;
for i = n : -1: j + 1 ;
temp = temp - conj(a(i,j)).*x(i);
end; i = fix(j + 1 -1);
if( nounit )
temp = temp./conj(a(j,j));
end;
end;
x(j) = temp;
end; j = fix(1 -1);
else;
kx = fix(kx +(n-1).*incx);
jx = fix(kx);
for j = n : -1: 1 ;
ix = fix(kx);
temp = x(jx);
if( noconj )
for i = n : -1: j + 1 ;
temp = temp - a(i,j).*x(ix);
ix = fix(ix - incx);
end; i = fix(j + 1 -1);
if( nounit )
temp = temp./a(j,j);
end;
else;
for i = n : -1: j + 1 ;
temp = temp - conj(a(i,j)).*x(ix);
ix = fix(ix - incx);
end; i = fix(j + 1 -1);
if( nounit )
temp = temp./conj(a(j,j));
end;
end;
x(jx) = temp;
jx = fix(jx - incx);
end; j = fix(1 -1);
end;
%
%
%     end of CTRSV .
%
a_shape=zeros(a_shape);a_shape(:)=a(1:numel(a_shape));a=a_shape;
x_shape=zeros(x_shape);x_shape(:)=x(1:numel(x_shape));x=x_shape;
end
%DECK CUCHK