Forward Time Centered Space Approach to Solve a Partial Differential Equation

Let us assume that there is a temperature distribution across a rectangular plate which is defined by the following unsteady one-dimensional convection-diffusion equation.

T_t + uT_x = αT_xx

Inside the plate, there is a flow of a fluid whose initial velocity is zero. The thickness of the plate is 1 cm, and the thermal diffusivity, α is 0.01 cm²/s. The initial temperature distribution is represented as,

T(x, 0) = 100 x / L;    0 ≤ x ≤ L

And, the boundary conditions are,

T(0, t) = 0;   T(L, t) = 100;

Assume, u = 0.1 cm/s. Develop a MATLAB program that solves this problem with the Forward Time Centered Space (FTCS) Approach. Consider, Δx = 0.1 and ΔT = 0.5. Plot the results for various times steps, such as, 0.5, 2.5, 5, 10, and 50 seconds.

Forward Time Centered Space (FTCS) Approach:

The following MATLAB program implements the FTCS approach to solve the unsteady one-dimensional convection-diffusion equation using the input parameters and boundary conditions.

% Forward Time Centered Space (FTCS) Approach

close all;

clc;

% Input Properties from Problem

L=1.0; % (dimension in cm)

alpha=0.01;

u=0.1;

Tt=50.0; % Maximum simulation time

% Mesh/Grid size calculation

% Gird size definition along x and y axes

dx=0.1;

dt=0.5;

c=u*dt/dx;

A=alpha*(dt/dx^2);

% No of grid points

imax=(L/dx)+1; % Array dimension along the width

tmax=Tt/dt; % Array dimension along the height

ndim=imax-2;

% Convection-Diffusion Equation by FTCS Approach

T(:,1)=0.0;

for i=1:imax

    x(i)=(i-1)*dx;

    T(i,1)=100*x(i)/L; % Initial temperature distribution

end

for t=1:tmax-1

    time(t+1)=t*dt;

    T(1,t+1)=T(1,t);

    for i=2:imax-1

        T(i,t+1)=T(i,t)-(c/2.0)*(T(i+1,t)-T(i-1,t))+A*(T(i+1,t)-2*T(i,t)+T(i-1,t));

    end

    T(imax,t+1)=T(imax,1);

    

end

hold on

plot(x,T(:,1),'LineWidth',2)

plot(x,T(:,2),'-p','LineWidth',2)

plot(x,T(:,6),'-*','LineWidth',2)

plot(x,T(:,11),'-s','LineWidth',2)

plot(x,T(:,21),'--','LineWidth',2)

plot(x,T(:,100),'^-','LineWidth',2)

xlabel('x, cm')

ylabel('Temperature, \circ C')

legend('t=0s','t=0.5s','t=2.5s','t=5s', 't=10s', 't=50s')

Program Output: