ODE 2-го порядка приблизительно с использованием конечной разности - PullRequest
0 голосов
/ 22 февраля 2020

Я пытаюсь приблизить и построить решение для u "(x) = exp (x) в интервале 0-3, с граничными условиями x (0) = 1, x (3) = 3. Я могу построить приблизительное решение против точного, но график выглядит немного не так:

% Interval
a=0;
b=3;
n=10;

% Boundary vals
alpha=1;
beta=3;
%grid size
h=(b-a)/(n+1); 


%Matrix generation
m = -2;
u = 1;
l = 1;

% Obtained from (y(i-1) -2y(i) + y(i+1)) = h^2 exp(x(i)


M = (1/h^2).*(diag(m*ones(1,n)) + diag(u*ones(1,n-1),1) + diag(l*ones(1,n-1),-1));

B=[];
xjj=[];


for j=1:n
    xjj=[xjj,j*h];
    if j==1
        B=[B,f(j*h)-(alpha/h^2)];

        continue 
    end
    if j==n
        B=[B,f(j*h)-(beta/h^2)];
        continue 
    else
        B=[B,f(j*h)];
    end
end

X=M\B';



x=linspace(0,3,101);
plot(xjj',X,'r*')
hold on 
plot(x,exp(x),'b-')

Я ценю все советы и объяснения. Вот схема, которой я следую: http://web.mit.edu/10.001/Web/Course_Notes/Differential_Equations_Notes/node9.html

1 Ответ

1 голос
/ 22 февраля 2020

Вы можете сократить большое l oop до простого

x=linspace(a,b,n+2); 
B = f(x(2:end-1)); 
B(1)-=alpha/h^2; 
B(n)-=beta/h^2;

Точное решение - u(x)=C*x+D+exp(x), граничные условия дают D=0 и 3*C+exp(3)=3 <=> C=1-exp(3)/3.

Plotting это точное решение по сравнению с численным решением дает довольно хорошее соответствие для этого большого размера шага:

enter image description here

f=@(x)exp(x)
a=0; b=3;
n=10;

% Boundary vals
alpha=1; beta=3;
%grid 
x=linspace(a,b,n+2); 
h=x(2)-x(1); 

% M*u=B obtained from (u(i-1) -2u(i) + u(i+1)) = h^2 exp(x(i))
M = (1/h^2).*(diag(-2*ones(1,n)) + diag(1*ones(1,n-1),1) + diag(1*ones(1,n-1),-1));
B = f(x(2:end-1)); 
B(1)-=alpha/h^2; B(n)-=beta/h^2;

U=M\B';

U = [ alpha; U; beta ];
clf;
plot(x',U,'r*')
hold on 
x=linspace(0,3,101);
C = 1-exp(3)/3
plot(x,exp(x)+C*x,'b-')
...