Index exceeds the number of array elements (23).
1 次查看(过去 30 天)
显示 更早的评论
function main
clc;
clear all;
eps = 0.001;
omega = input(' enter the omega value: ');
beta = input (' enter the beta value: ');
n= 100000;
nx = 26;
ny = 26;
S(1:nx, 1:ny) = 0;
SN(1:nx, 1:ny) = 0;
S(7:nx, 1)= 100;
S(nx, 1:16) = 100;
SN(7:nx, 1)= 100;
SN(nx, 1:16) = 100;
% its number of iteration
coeff = ( 2*(1+beta^2));
for iterations = 1:n
for j = 2:ny-1
a(1:nx-2) = -coeff;
b(1:nx-3)= omega;
c(1:nx-3)= omega;
for i = 2:nx-1
r(i-1)= -coeff*(1-omega)*S(i,j)- omega*beta^2*S(i,j+1)-omega*beta^2*SN(i,j-1);
end
r(1)= r(1)-omega*SN(1,j);
r(nx-2)= r(nx-2)- omega*SN(nx,j);
y = Tridiagonal (c,a,b,r);
for k = 1:nx-2
SN(k+1,j) = y(k);
end
end
for i = 2:nx-1
a(1:ny-2) = -coeff;
b(1:ny-3) = beta*beta;
c(1:ny-3) = beta*beta;
for j = 2:ny-1
r(j-1) = -coeff*(1-omega)*S(i,j)- omega*S(i+1,j)-omega*SN(i-1,j);
end
r(1) = r(1)-SN(i,1);
r(ny-2)= r(ny-2)- SN(i,ny);
y = Tridiagonal (c,a,b,r);
for k = 1:ny-2
SN (i,k+1)= y(k);
end
end
error = abs(SN-S);
totalerror = sum(error,'all');
if totalerror <= eps
break
end
S = SN;
end
iterations
contour(SN');
end
function x = Tridiagonal(e,f,g,r)
% Tridiagonal: Tridiagonal equation solver banded system
% x = Tridiagonal(e,f,g,r): Tridiagonal system solver.
% input:
% e = subdiagonal vector
% f = diagonal vector
% g = superdiagonal vector
% r = right hand side vector
% output:
% x = solution vector
n=length(f);
% forward elimination
for k = 2:n
factor = e(k)/f(k-1);
f(k) = f(k) - factor*g(k-1);
r(k) = r(k) - factor*r(k-1);
end
% back substitution
x(n) = r(n)/f(n);
for k = n-1:-1:1
x(k) = (r(k)-g(k)*x(k+1))/f(k);
end
end
COMMAND WINDOW;
enter the omega value: 1.3
enter the beta value: 1
Index exceeds the number of array elements (23).
Error in adior1>Tridiagonal (line 71)
factor = e(k)/f(k-1);
Error in adior1 (line 28)
y = Tridiagonal (c,a,b,r);
0 个评论
回答(1 个)
Walter Roberson
2021-10-10
nx = 26;
a(1:nx-2) = -coeff;
You do not otherwise initialize a before you call the function, so a is going to be length nx-2 which would be 26-2 = 24
c(1:nx-3)= omega;
You do not otherwise initialize c before you call the function, so c is going to be length nx-3 which would be 26-3 = 23
y = Tridiagonal (c,a,b,r);
c (shorter vector) gets passed in first position, a (longer vector) gets passed in second position
function x = Tridiagonal(e,f,g,r)
first position (shorter vector) gets named e inside the function; second position (longer vector) gets named f inside the function.
n=length(f);
That is length() of the second position, the longer vector. So it would be n = 24
for k = 2:n
k can be a maximum of 24
factor = e(k)/f(k-1);
e is the shorter vector, and would be indexed up to n = 24. But e is shorter and only has 23 elements.
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Operating on Diagonal Matrices 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!