Problem with finite difference
3 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
Paul Eluard
il 12 Nov 2017
Risposto: David Goodmanson
il 13 Nov 2017
Hello !
In order to solve the following ODE : $U_{xx} + U_x - 2U + 2 = 0$ with boundary condition $U(-1) = 0 = U(1)$, I use centered differences of order 2 :
$U_{xx} = \frac{U_{i+1} - 2 U_i + U_{i-1}}{h^2} + O(h^2)$
$U_x = \frac{U_{i+1} - U_{i-1}}{2h} + O(h^2)$
Putting this altogether gives me the following code :
Mxx = gallery('tridiag',n+1,1,-2,1)*1/h^2; Mx = gallery('tridiag',n+1,1,0,-1)*1/(2*h); M0 = speye(n+1); M = Mxx + Mx - 2*M0; M(1,1) = 1; M(1,2) = 0; M(end,end) = 1; M(end,end-1) = 0; M = sparse(M); B = [0;-2*ones(n-1,1);0]; U = M\B;
But the solution doesn't fit the analytical one at all... $u(x) = 1 - \frac{\sinh(2)e^x + \sinh(1)e^{-2x}}{\sinh(3)}$
![](/matlabcentral/answers/uploaded_files/94340/untitled.jpg)
Where is the mistake ?
0 Commenti
Risposta accettata
David Goodmanson
il 13 Nov 2017
Hi Paul,
You did not say what h is, but I took it to be 2/n. After replacing Mx with its transpose, which changes the sign of the derivative (the results say it's now the correct sign but I didn't sit down and prove it) you get really good agreement. The plot is basically an overlay.
n = 1000;
h=2/n;
Mxx = gallery('tridiag',n+1,1,-2,1)*1/h^2;
Mx = gallery('tridiag',n+1,1,0,-1)'*1/(2*h); % transposed
M0 = speye(n+1);
M = Mxx + Mx - 2*M0;
M(1,1) = 1; M(1,2) = 0;
M(end,end) = 1; M(end,end-1) = 0;
M = sparse(M);
B = [0;-2*ones(n-1,1);0];
U = M\B;
x = linspace(-1,1,1001);
u = 1 - (sinh(2)*exp(x) + sinh(1)*exp(-2*x))/sinh(3);
figure(1)
plot(1:n+1,U',1:n+1,u)
max(abs(U'-u))
0 Commenti
Più risposte (0)
Vedere anche
Categorie
Scopri di più su Ordinary Differential Equations in Help Center e File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!