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

1 Ansicht (letzte 30 Tage)
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

Antworten (1)

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

Kategorien

Mehr zu Matrix Indexing finden Sie in Help Center und File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by