# Can someone help me figure out what is wrong with my function c1 = ult_trid(m,b) ??

1 view (last 30 days)
daniel choudhry on 8 Oct 2020
Answered: James Tursa on 9 Oct 2020
I cannot get my code to run because it tells me:Index in position 2 exceeds array bounds (must not exceed 1).
which i cannot figure out what part of my function c1 = ult_trid(m,b) is wrong. can someone help me?
function Math_5532_6_6_P1
rand('seed',1);
b = rand(4,1);
c = rand(4,1);
d = rand(4,1);
e = rand(4,1);
c(1) = 0;
e(4) = 0;
d1 = d;
[m,d,cri] = lu_fac_trid(c,d,e);
if cri == 0
disp('A is not invertible');
return
end
disp('--------------------------');
disp('max_multiplier');
mm =max(abs(m))
disp('--------------------------');
c1 = ult_trid(m,b);
x = ut_trid(d,e,c1);
%check solution:
r = zeros(4,1);
r(1) = b(1)- (d1(1) * x(1) + e(1) * x(2));
r(2) = b(2)- (c1(2) * x(1) + d1(2) * x(2) + e(2) * x(3));
r(3) = b(3)- (c1(3) * x(2) + d1(3) * x(3) + e(3) * x(4));
r(4) = b(4)- (c1(4) * x(3) + d1(4) * x(4));
disp('The residual vector');
r
disp('--------------------------');
%%%%%%%%%%%%%%%
function [m,d,cri] = lu_fac_trid(c,d,e)
n = length(d);
cri = 1;
m = zeros(n,1);
%Check for invertibility
for k = 1:n-1
if max(abs(d(k)),abs(c(k+1))) < n*10^(-15)
cri = 0;
return
end
m(k+1) = d(k);
d(k+1) = c(k+1);
end
if abs(d(n)) < n*10^(-15)
cri = 0;
return
end
for i=n:-1:k+1
factor= c(i)/c(k);
factor2= e(i)*d(k);
m(i)=c(i)-c(k)*(factor);
d(i)=d(i)-d(k)*(factor2);
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%
function c1 = ult_trid(m,b)
[brows, bcol] = size(b);
[mrows, mcol] = size(m);
c1 = b;
for k = 1:mrows
for j = 1:k-1
c1(k) = c1(j) - m(j)*c1(j);
end
c1(k) = c1(k)/m(k);
end
end
%%
function x=ut_trid(d,e,c1)
% upper triangle
n=length(d);
x=zeros(n,1);
x(n)=c1(n)/d(n,n);
for i=n-1:-1:1
x(i) = (c1(i)-e(i)*x(i+1)) /d(i);
end
end
end

James Tursa on 9 Oct 2020
d is a 4x1 vector and you are trying to access the d(4,4) element, hence the error.