Question

MatLab

Is this algorithm correct for solving the heat equation with a Crank-Nicolson scheme? Thanks!

Answer 1

clc;
clear;

n = 100;
m = 100;
dt = 0.0000049;
dx = 0.01;
k = 1;
s = k * dt / ( dx ^2 );
ui = 10;
ud = 5;

A = zeros( n, m );

A( 1, 1 ) = ( 1 + s );
A( 1, 2 ) = -s / 2;
A( n, m ) = ( 1 + s );
A( n, m - 1 ) = -s / 2;

for i = 2 : n - 1

A( i, i - 1 ) = -s / 2;
A( i, i ) = ( 1 + s );
A( i, i + 1 ) = -s / 2;

end

u = zeros( n, 1 );
d = zeros( n, 1 );

d( 1 ) = ui;
d( n ) = ud;

for i = 2 : n - 1

d( i ) = s / 2 * u ( i + 1 ) + ( 1 - s ) * u( i ) + s / 2 * u ( i - 1 );

end

for i = 1 : m

u = A \ d;

d = u;

d( 1 ) = ui;
d( n ) = ud;

end