Error in Finite Difference Method for ODE Code

14 Ansichten (letzte 30 Tage)
Mel Hernandez
Mel Hernandez am 30 Apr. 2018
Kommentiert: Mel Hernandez am 2 Mai 2018
I'm working on an ODE Finite Difference Method. I think I have it down (so I deleted my previous question), but I am getting this error:
Array indices must be positive integers or logical values.
Error in FiniteDifference (line 18)
p(X) = (1 + 0*X); q(X) = (2 + 0*X); r(X) = cos(X);
Here is my code:
function [A,B,W] = FiniteDifference(N)
%equation
%y'' = y' + 2y + cosx, 0<x<pi/2
%y(0) = -0.3, y(pi/2) = -0.1
%breaking it down:
% a = 0; alpha = -0.3; b = pi/2; beta = -0.1;
% p(x) = 1; q(x) = 2; r(x) = cosx;
%to check: plot and residue
A = zeros(1,N+1); B = zeros(1,N+1); C = zeros(1,N+1); D = zeros(1,N+1);
L = zeros(1,N+1); U = zeros(1,N+1); Z = zeros(1,N+1); W = zeros(1,N+1);
aa = 0; alpha = -0.3;
bb = pi/2; beta = -0.1;
h = (bb-aa)/(N+1); %Step Size
X = aa + h;
p(X) = (1 + 0*X); q(X) = (2 + 0*X); r(X) = cos(X);
A(1) = 2+h^2*q(X);
B(1) = -1+0.5*h*p(X);
D(1) = -h^2*r(X)+(1+0.5*h*p(X))*alpha;
for i = 2:N-1
X = aa + i*h;
A(i) = 2+h^2*q(X);
B(i) = -1+0.5*h*p(X);
C(i) = -1-0.5*h*p(X);
D(i) = -h^2*r(X);
end
X = bb-h;
A(N) = 2 + h^2*q(X);
C(N) = -1 - 0.5*h*p(X);
D(N) = -h^2*r(X) + (1 - 0.5*h*p(X))*beta;
L(1) = A(1);
U(1) = B(1)/A(1);
Z(1) = D(1)/L(1);
for i = 2:N-1
L(i) = A(i) - C(i)*U(i-1);
U(i) = B(i)/L(i);
Z(i) = (D(i) - C(i)*Z(i-1))/L(i);
end
L(N) = A(N) - C(N)*U(N-1);
Z(N) = (D(N) - C(N)*Z(N-1))/L(N);
W(N) = Z(N);
for j = 1:N-1
i = N-j;
W(i) = Z(i) - U(i)*W(i+1);
end
i = 0;
for i = 1:N
X = aa + i*h;
end
plot(X,W);
end

Akzeptierte Antwort

Torsten
Torsten am 30 Apr. 2018
Bearbeitet: Torsten am 30 Apr. 2018
Maybe you mean
p(X) = @(X)(1 + 0*X); q(X) = @(X)(2 + 0*X); r(X) = @(X)cos(X);
?
Best wishes
Torsten.
  4 Kommentare
Jan
Jan am 2 Mai 2018
r = @(X)cos(X)? What about r=@cos? Anonymous functions are much slower than function handles. q = @(X)(2 + 0*X)? Looks like: q = 2. Do I oversee something?
Mel Hernandez
Mel Hernandez am 2 Mai 2018
Tortsen, thank you so much!!!

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (0)

Kategorien

Mehr zu Mathematics and Optimization 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