[n,a,b,c,d,u,z]=trisp(n,a,b,c,d,u,z);

function [n,a,b,c,d,u,z]=trisp(n,a,b,c,d,u,z);
persistent an bn den j k nm1 nm2 v ; 

if isempty(an), an=0; end;
if isempty(bn), bn=0; end;
if isempty(den), den=0; end;
if isempty(v), v=0; end;
if isempty(j), j=0; end;
if isempty(k), k=0; end;
if isempty(nm1), nm1=0; end;
if isempty(nm2), nm2=0; end;
%***BEGIN PROLOGUE  TRISP
%***SUBSIDIARY
%***PURPOSE  Subsidiary to SEPELI
%***LIBRARY   SLATEC
%***TYPE      SINGLE PRECISION (TRISP-S)
%***AUTHOR  (UNKNOWN)
%***DESCRIPTION
%
%     This subroutine solves for a non-zero eigenvector corresponding
%     to the zero eigenvalue of the transpose of the rank
%     deficient ONE matrix with subdiagonal A, diagonal B, and
%     superdiagonal C , with A(1) in the (1,N) position, with
%     C(N) in the (N,1) position, and all other elements zero.
%
%***SEE ALSO  SEPELI
%***ROUTINES CALLED  (NONE)
%***REVISION HISTORY  (YYMMDD)
%   801001  DATE WRITTEN
%   890831  Modified array declarations.  (WRB)
%   891214  Prologue converted to Version 4.0 format.  (BAB)
%   900402  Added TYPE section.  (WRB)
%***end PROLOGUE  TRISP
%
a_shape=size(a);a=reshape(a,1,[]);
b_shape=size(b);b=reshape(b,1,[]);
c_shape=size(c);c=reshape(c,1,[]);
d_shape=size(d);d=reshape(d,1,[]);
u_shape=size(u);u=reshape(u,1,[]);
z_shape=size(z);z=reshape(z,1,[]);
%***FIRST EXECUTABLE STATEMENT  TRISP
bn = b(n);
d(1) = a(2)./b(1);
v = a(1);
u(1) = c(n)./b(1);
nm2 = fix(n - 2);
for j = 2 : nm2;
den = b(j) - c(j-1).*d(j-1);
d(j) = a(j+1)./den;
u(j) = -c(j-1).*u(j-1)./den;
bn = bn - v.*u(j-1);
v = -v.*d(j-1);
end; j = fix(nm2+1);
den = b(n-1) - c(n-2).*d(n-2);
d(n-1) =(a(n)-c(n-2).*u(n-2))./den;
an = c(n-1) - v.*d(n-2);
bn = bn - v.*u(n-2);
den = bn - an.*d(n-1);
%
%     SET LAST COMPONENT EQUAL TO ONE
%
z(n) = 1.0;
z(n-1) = -d(n-1);
nm1 = fix(n - 1);
for j = 2 : nm1;
k = fix(n - j);
z(k) = -d(k).*z(k+1) - u(k).*z(n);
end; j = fix(nm1+1);
a_shape=zeros(a_shape);a_shape(:)=a(1:numel(a_shape));a=a_shape;
b_shape=zeros(b_shape);b_shape(:)=b(1:numel(b_shape));b=b_shape;
c_shape=zeros(c_shape);c_shape(:)=c(1:numel(c_shape));c=c_shape;
d_shape=zeros(d_shape);d_shape(:)=d(1:numel(d_shape));d=d_shape;
u_shape=zeros(u_shape);u_shape(:)=u(1:numel(u_shape));u=u_shape;
z_shape=zeros(z_shape);z_shape(:)=z(1:numel(z_shape));z=z_shape;
end
%DECK TRIX