How can I solve this partial differentiation equation numerically?

Hi everyone

I'm trying to solve a PDE by Euler's method : the implicit one. I've written a code to solve it first by the explicit method and I got the results correctly, but the results given by the implicit method are not right. It must be better than the first method because it is unconditionally stable, but that's not happening...

I attached my code so that you can observe any errors and help me.

Answer 1

[Torsten]

A = zeros(M+1);

U2 = zeros(M+1,N+1);

A(1,1) = 1;

for i=2:M

A(i,i-1) = -sigma;

A(i,i) = 1+2*sigma;

A(i,i+1) = -sigma;

end

A(M+1,M+1) = 1;

U2(:,1) = linsolve(A,U(:,1));

for i=2,N+1

U2(:,i) = linsolve(A,U2(:,i-1));

end

Note that the code could be made faster by factorizing A once and solving the linear systems for only changing right-hand sides. This is possible since A is constant for all times.