Matlab code for the Formula

Following is the Analytical solution of a Heat Diffusion formula where x is length of Rod from 0 to 1 meter with an increment of 0.1. Also t is the time from 0 to 10 minute with an increment of 1 min. Value of r is 0.1. The result should almost match with the following table

. The result in the table was obtained using the Numerical Method (Explicit Method) of the problem.

I've tried the following code:

close all clear clc

x = 0:0.1:1; t = 0:1:10; r = 0.1; [X,T] = meshgrid(x,t) S = 0;

for n=1:1e6 S = S + sin(pi*n/2)*sin(n*pi*X).*exp(-n*n*r*r*t)/(n*n); end U = (8/pi*pi)*S

What is wrong with this code?

[John D'Errico]

What have you tried? If nothing, why not? It is by making an effort that you will learn, not by being given the solution to your homework on a platter. So show what you tried. If you do, you will have abetter chance of getting some help here.

[Rasel Munshi]

close all
clear
clc

x = 0:0.1:1;
t = 0:1:10;
r = 0.1;
[X,T] = meshgrid(x,t)
S = 0;

for n=1:1e6
S = S + sin(pi*n/2)*sin(n*pi*X).*exp(-n*n*r*r*t)/(n*n);
end
U = (8/pi*pi)*S

Answer 1

[Image Analyst]

One way is a more direct, though not vectorized, approach of using for loops over t and x:

x = 0:0.1:1;
t = 0:1:10;
r = 0.1;
U = zeros(length(t), length(x));
for kt = 1 : length(t)
  for kx = 1 : length(x)
    s = 0;
    for n=1:1e6
      s = s + sin(pi*n/2).*sin(n*pi*x(kx)).*exp(-n*n*r*r*t(kt))/(n*n);
    end
    U(kt, kx) = 8 * s / pi^2;
  end
end
U

Answer 2

[Abraham Boayue]

Hey Rasel, your code was almost right. You made just a few errors, instead of initializing the sum as you did, you should have just initialized U as a matrix. This is your code with a few changes.

clear variables
close all
x = 0:0.1:1;
t = 0:10;
r = 0.1;
N = length(x);
M = length(t);
U = zeros(N,M);   % This is what you should have done to start with.
[x,t] = meshgrid(x,t);
for n = 1: 100
    U = U + ((1/n^2)*sin(0.5*pi*n)*sin(pi*n*x).*exp(-n^2*r^2*t))';
end
U = (8/pi.^2)*U;
disp(U')
figure
surf(x,t,U')
grid