# One Dimensional Heat Equation

**Disclaimer:** This dissertation has been submitted by a student. This is not an example of the work written by our professional dissertation writers. You can view samples of our professional work here.

Any opinions, findings, conclusions or recommendations expressed in this material are those of the authors and do not necessarily reflect the views of UK Essays.

**ABSTRACT**

**CHAPTER ONE**

**INTRODUCTION**

It's deeply truth that the search for the exact solution in our world-problems is needed for each of us, but unfortunately, not all problems can be solved exactly; because of nonlinearity and complicated geometry. In many interesting applications, exact solutions may be impossible and unattainable, so we need to approximate. One of the most beautiful branch in mathematics science can deal with these kind of problems, which is (numerical methods). From the earliest times, the development of the high speed-digital computers and modern electronic devices is widely enhances the use of numerical methods in many various of branch in science and engineering (Fausett, 2002; Ozisik, 1994).

In this dissertation we study the application (in particular accuracy and convergence) of two kinds of numerical methods based on Adomian Decomposition Method (ADM) and Forward Time Central Space (FTCS) method when applied to the heat equation.

Heat equation is an important partial differential equation (pde) used to describe various phenomena in many applications of our daily life.

**1.1** **Numerical methods**

One of the earliest mathematical writers in this field was by the Babylonians (3,700 years ago). There is evidence that they knew how to find the numerical solutions for quadratic equations, also they approximated the root for an integer number. **Numerical methods are:** the methods often used to solve ordinary differential equations (ode's) and partial differential equations (ode's) (Fauset, 2002).

The finite difference method and finite element method are an example of the numerical methods.

The finite difference method will be considered in more details in chapter 2.

**PDE's (background)**

Many physical phenomena can be described by partial differential equations (pde's), to better understanding the phenomena, we need to formulate and solve it using the pde's. so , what is the pde's?

**Partial differential equations (PDE's****) :** are the class of differential equations involving unknown functions of several variables (Wong Y.Y,W,.T.C,J.M,2005).

A pde for a function u(x,y) can be written as:

Fx,y,u,âˆ‚uâˆ‚x,âˆ‚uâˆ‚y,âˆ‚2uâˆ‚x2,âˆ‚2uâˆ‚y2,âˆ‚uâˆ‚xâˆ‚y,â€¦=0

The order of the pde is the highest order of its derivatives.

There are three types of the pde's, elliptic, parabolic and hyperbolic, where they are classified in order of their discriminates.

For example, consider the second order pde in 2 variables:

Auxx+Buxy+Cuyy=f(x,y,u,ux,uy)

Where A, B, C and f are functions of u. So the pde is said to be:

Hyperbolic pde if: B2-4AC>0 .For example:Â uxx-utt=0 .(wave equation) .

Elliptic pde if : B2-4AC<0 .For example uxx+utt=0.(laplace equation)

Parabolic pde if : B2-4AC=0.For example uxx-ut=0.(Heat equation).Â (Wong Y.Y,W,.T.C,J.M,2005).

**Heat Equation**

One dimensional heat equation that will be used in this dissertation is:

âˆ‚uâˆ‚t=Î±âˆ‚2uâˆ‚x2Â , Â for 0â‰¤xâ‰¤xfÂ ,Â 0â‰¤tâ‰¤T

Â Where Î± is positive constant. In order for this equation to be solved, the initial conditions (IC) and the boundary conditions (BC) should be found. (Wong Y.Y,W,.T.C,J.M,2005).

Â The one dimensional heat equation describes the distribution of heat, heat equation almost known as *diffusion equation*; it can arise in many fields and situations such as: physical phenomena, chemical phenomena, biological phenomena.

Heat equation will be considered in our study under specific conditions.

**Objectives**

The aims of this dissertation are:

1. To study the recent literature survey on FTCS method and Adomian Decomposition Method for heat equation.

2. To investigate the performance (in particular accuracy and stability) of these methods when applied to a real problem rather than a simple illustrative problem.

3. To compare between FTCS method and Adomian Decomposition Method under different values.

**Overview**

Â This study consists of six chapters. Chapter one provides the topic to be studied and the aims to be achieved. Chapter two states basis of the finite difference method and some illustrative examples. Chapter three and chapter four we will discuss in greater detail the FTCS method and Adomian decomposition method respectively. Recent literature review on FTCS method and ADM method for heat equation are presented.

Â Chapter five presents numerical experiments for one dimensional heat equation using MATLAB software for these methods. Chapter six provides the conclusion and suggestions for further future work.

**CHAPTER TWO**

**Finite Difference Method**

**Introduction**

Finite difference method (FDM) is widely used for the solution of partial differential equations of heat, mass and momentum transfer. Furthermore, it is very easy to learn and apply to many problems in physics and engineering. Extensive amount of literature exist on the applications of this method for solving such problems. The accuracy of the (FDM) can be easily tested by the order of the truncation error in the Taylor's series expansion (Ozisik, 1994).

Â The idea of the (FDM) proceeds by replacing the derivatives in the partial differential equations with finite difference approximations. This gives large systems of algebraic equations which usually can be solved in the area of partial differential equations. Algebraic equations can readily solved using a computer (LeVeqe, 2007).

In this chapter we will treat to the Finite difference method and other related topics in more details. Some of what this chapter includes is adapted from the Lecture notes book MAT518 (2008).

**Finite Difference Approximations**

Â Suppose a function u=u(x,y). Divide the first quadrant of the xy plan into uniform rectangles by grid lines parallel to x axis (uniform length âˆ†y) and grid lines parallel to y axis (uniform length âˆ†x). See Figure (2.1).

Â Assume we have a continues function u, so the Tylor's series expansion for uxi+âˆ†x,yj about (xi,yj) is :

uxi+âˆ†x,yj=uxi,yj+âˆ†x1!uxxi,yj+âˆ†x22!uxxxi,yj+âˆ†x33!uxxxxi,yj+âˆ†x44!uxxxxxi,yj+â€¦ (2.1)

Â Suppose we decide to truncate the series on right hand side (RHD) of (2.1) beginning with the 3rd term. If âˆ†x is sufficiently small, he 4thÂ and higher terms are much smaller than the 3rd term. Then we write:

uxi+âˆ†x,yj=uxi,yj+âˆ†x1!uxxi,yj+0âˆ†x2 (2.2)

0(âˆ†x2) Means the sum of the truncated terms (i.e. the truncation error) is in absolute term most a constant multiple ofÂ âˆ†x2 .Divide (2.2) by (âˆ†x) and after rearranging the terms once we get:

uxxi+âˆ†x,yj=uxi+âˆ†x,yj-uxi,yjâˆ†x+0âˆ†x Â (2.3)

uxi+âˆ†x,yj-uxi,yjâˆ†x Which is called a forward difference approximation to the partial derivative ofÂ âˆ‚uâˆ‚x at the points xi,yj and it is known to be the first order accurate or {0(âˆ†x)} accurate. IfÂ âˆ†x is reduced by fifty-present (50%), then the truncation error is likewise reduced by 50%.

The Taylor's series expansion for uxi-âˆ†x,yj about (xi,yj) is :

uxi+âˆ†x,yj=uxi,yj-âˆ†x1!uxxi,yj+âˆ†x22!uxxxi,yj-âˆ†x33!uxxxxi,yj+âˆ†x44!uxxxxxi,yj+â€¦ (2.4)

Â Suppose we decide to truncate the series on right hand side (RHS) of equation (2.4) beginning with the 3rd term. If âˆ†x is sufficiently small, the 4th and higher terms are much smaller than the 3rd term. So we write:

uxi-âˆ†x,yj=uxi,yj-âˆ†x1!uxxi,yj+0âˆ†x2 (2.5)

By dividing equation (2.5) by (âˆ†x) and rearrange the terms once we get:

uxxi,yj=uxi+âˆ†x,yj-uxi,yjâˆ†x+0âˆ†x (2.6)

uxi+âˆ†x,yj-uxi,yjâˆ†x Which is called a backward difference approximation to the partial derivative of âˆ‚uâˆ‚x at the points uxi,yj and it is known to be first order accurate or {0(âˆ†x)} accurate.

If we subtract equation (2.4) from equation (2.1) we get:

uxi+âˆ†x,yj- uxi-âˆ†x,yj=2âˆ†xuxxi+yj+0âˆ†x3 Â (2.7)

By dividing equation (2.7) by (2âˆ†x) and rearrange the terms once we get:

uxxi,yj=uxi+âˆ†x,yj-uxi-âˆ†x,yj2âˆ†x+0âˆ†x2 Â (2.8)

uxi+âˆ†x,yj-uxi-âˆ†x,yj2âˆ†x is called central difference approximation to the partial derivative âˆ‚uâˆ‚x at the points uxi,yj which is known as second order accurate or {0(âˆ†x2)} (Lecture note book MAT 518,2008).

Consider the figure (2.2) below, it is simple that the equation (2.8) approximates the slop of the tangent at point P by the slop of the chord AB, and is said to be central difference approximation** .** Moreover, we can also approximate the slope of tangent P by the slop of chord PB which gives a forward difference approximation formula

**Likewise, the slop of the chord AP gives a backward difference approximation formula (smith, 1985)**

*.*

*.*Adding equation (2.1) and equation (2.4) we get a similar approximation for âˆ‚2uâˆ‚x2Â at the points uxi,yj

uxi+âˆ†x,yj+ uxi-âˆ†x,yj=2uxi,yj+âˆ†x2uxxxi+yj+0(âˆ†x4) Â (2.9)

Dividing (2.9) by (âˆ†x2) Â and rearrange the terms to get :

uxxxi+yj=xi+âˆ†x,yj-2xi+yj+uxi-âˆ†x,yjâˆ†x2+0âˆ†x2 Â (2.10)

xi+âˆ†x,yj-2xi+yj+uxi-âˆ†x,yjâˆ†x2Â Which is called the central approximation to the second partial derivative âˆ‚2uâˆ‚x2Â and it is the second order {(âˆ†x2)} accurate.

Â Similarly as the previous steps, we can get the approximations to the partial derivates uy and uyy to obtain a *forward, backward *and* central difference* *approximation *âˆ†y

(Look at Figure 2.3) as follows:

uyxi,yj= uxi,yj+âˆ†y-uxi,yjâˆ†yÂ (Forward difference approximation) Â (2.11)

uyxi,yj= uxi,yj-uxi,yj-âˆ†yâˆ†y (Backward difference approximation) (2.12)

Â uyxi,yj= uxi,yj+âˆ†y-uxi,yj-âˆ†y2âˆ†y(Central difference approximation) (2.13)

uyyxi,yj= uxi,yj+âˆ†y-2xi,yj+uxi,yj-âˆ†y2âˆ†y2Â (Central difference approximation) Â (2.14)

Finite difference method is an effective method, but it also suffers from a drawback. A major drawback of the finite difference method appears in its incapability to deal efficiently with the solution of the problems over arbitrarily shaped complex geometries. This happens because of the interpolation difficulties between the interior pints and the boundaries conditions in order to construct the (FDM) expressions for nodes points next the boundaries conditions (Ozisik, 1994).

**CTCS Method**

We illustrate the CTCS method (Central Time Central Space) as an example of the finite difference method, this example is taken from (Burden; Faires, 2005).

Example 2.1:

Consider the IBVB below:

âˆ‚2uâˆ‚t2-âˆ‚2uâˆ‚x2=0, Â 0<x<1 , 0<tÂ ; Â (2.15)

With boundary conditions:

u0,t=u1,t=0Â ,Â 0<t;

And initial conditions:

ux,0=sinÏ€xÂ , 0â‰¤xâ‰¤1 ,

utx,0=0,Â 0â‰¤xâ‰¤1

Divide the domain into a set of uniform rectangles with grid lines parallel to the x axis and grid lines parallel to the t axis. Suppose we use âˆ†x=0.20 (i.e. we use 5 equal sub-intervals), consider the figure2.4 below:

The grid points are (xi,yi) where i = 1,2,3,4 and xi=iâˆ†x, and n=1,2,3,4,â€¦Â and tn=nâˆ†t. The initial conditions satisfies the values of u at points x0,x1,x2,x3,x4,x5 at t = 0. And The boundary conditions satisfies the values of u at points x=x0 and x=x5Â âˆ€t>0.

We use the central differences forÂ âˆ‚2uâˆ‚t2Â and âˆ‚2uâˆ‚x2Â once we get:

uin+1-2uin+uin-1âˆ†t2=ui+1n-2uin+ui-1nâˆ†x2 Â (2.16)

uin+1=21-r2uin+r2ui+1n+ui-1n-uin-1 (2.17)

Where r in our example is:

r= Â âˆ†tâˆ†x Â Â (2.18)

Equation (2.16) is called CTCS formula, which is an explicit formula.

If values of u at time level(n) and (n+1)

To implement the stability in CTCS method we make sure that the time step âˆ†t must be chosen such as râ‰¤1 which meansâˆ†tâ‰¤âˆ†x. We assume that âˆ†t=0.20 and âˆ†x=0.20 , so r=1.

Substituting the value of r in the formula (2.17) gives:

uin+1=ui+1n+ui-1n-uin-1 Â Â (2.19)

We now study the implementation of CTCS, so we start with:

i=1,n=0

u11=u20+u00-u1-1 (2.20)

Â All the values on right hand side are known exceptÂ u1-1, so we can compute the value of u1-1Â by using the central difference to the lift hand side of the initial condition utx,0=0, which gives u11=u1-1, substituting this value to equation (2.20)Â yields:

u11= u20+u002 Â (2.21)

Â Substituting the valuesÂ u00=0, u20=sin0.4Ï€=0.951 to the equation (2.21) we get, u11=0.476, which is the approximation to u(0.2,0.2).

i=2,n=0

u21=u30+u10-u2-1 Â Â (2.22)

Here also all values on right hand side are known exceptÂ u2-1, so we can compute the value of u2-1Â by using the central difference to the lift hand side of the initial condition utx,0=0, which gives similarlyÂ u21=u2-1, substituting this value to equation (2.22)Â gives:

u21=Â u30+u102 Â (2.23)

Using the values u10=sin0.2Ï€=0588, u30=sin0.6Ï€=0.951 in the equation (2.23) we get, u11=0.476, which is the approximation to u(0.4,0.2).

i=3,n=0

u31=u40+u20-u3-1 Â Â (2.24)

All the terms on right hand side are known exceptÂ u3-1, so we can compute the value of u2-1Â by using the central difference to the lift hand side of the initial condition utx,0=0, which gives similarlyÂ u31=u3-1, substituting this value to equation (2.24)Â gives:

u11=Â u40+u202 Â (2.25)

Â Using the values u20=sin0.4Ï€=0.951, u40=sin0.8Ï€=0.588 in the equation (2.25) we get, u31=0.769, which is the approximation to u(0.6,0.2).

i=4,n=0

u41=u50+u30-u4-1 Â (2.26)

All the terms on right hand side are known exceptÂ u4-1, so we can compute the value of u4-1Â by using the central difference to the lift hand side of the initial condition utx,0=0, which gives similarlyÂ u4-1=u41, substituting this value to equation (2.26)Â gives:

u41=u50+u302 Â (2.27)

Using the values u30=sin0.6Ï€=0.951, u50=0, in the equation (2.27)

we get, u41=0.476, which is the approximation to u(0.8,0.2).

By using the pervious values we can easily calculate the second row, so we start with:Â

i=1,n=1

u12=u21+u01-u10=0.769+0-0.88=0.181

i=2,n=1

u22=u31+u11-u20=0.769+0.476=0.951=0.294

i=3,n=1

u32=u41+u21-u30=0.476+0.764-0.951=0.289

i=4,n=1

u42=u51+u31-u40=0+0.769-0.588=0.181

So, the computational continues till final time level required.

**Crank â€“Nicolson method**

** Â **Crank â€“Nicolson method is an example of an implicit method, it can be use for solving diffusion equations (parabolic pde's), for general, consider the IBVP equation:

âˆ‚uâˆ‚t|ij+12=Âµ2âˆ‚u2âˆ‚x2|ij+12

We can approximate this equation by using the central difference approximation to âˆ‚u2âˆ‚x2 (space derivative) and central difference approximation to âˆ‚uâˆ‚tÂ (time derivative)

at pointsÂ [iâˆ†x, (j+12)âˆ†t], in other word, we can approximate the equation at non grid points (xi,tj+1/2) by:

Â uij+1-uijâˆ†t= Âµ2ui+1j+12-2uij+ui-1j+12âˆ†x2 (2.28)

RHS of equation (2.28) can be written as:

Â Âµ22ui+1j+12-2uij+ui-1j+12âˆ†x2+ui+1j-2uij+ui-1jâˆ†x2 Â (2.29)

Equation (2.28) can be written as:

-rui-1j+1+2-2uij+1-rui+1j+1=rui-1j+2-2uij-rui+1j Â (2.30)

Â This formula for Crank â€“Nicolson method and it is convergent and stable for all values of {r} (smith, 1985).

In Crank â€“Nicolson method, even we know the values of u at time level n from the IC and BC conditions, uij+1 cannot be calculated directly since ui-1j+1Â andÂ ui+1j+1 are known, that is why the Crank â€“Nicolson method is called *an implicit method*.

The following figure illustrates Crank â€“Nicolson method:

**The concept of stability**

Various of finite difference methods were applied for solving parabolic partial differential equations which forms the diffusion equations. All of these methods essentially reduce the solution of pde's to algebraic system of equations which can be done using the computer. It also known that the numerical computations carried out, even on a high-tech computer, it can be accurate just for a finite numbers of decimal places, which means, at any stage of calculations, someÂ round-offÂ error will be introduced. Let E be an error introduced into the calculations at all grid points{i,n,ui-n}. The resulting numerical solution of the differential equation will be:

âˆ‚uâˆ‚t =Î²2 âˆ‚2uâˆ‚x2 Â (2.31)

Let uin be the corresponding solution of the differential equation with finite difference method and suppose there is no round-off error in each calculation. Then we state:

Â Ein=uin-ui-n Â (2.32)

Stability is very important in getting the numerical solution of partial differential equations by using the finite difference methods. The solution of the finite difference equation is known to be stable, if the error does not increase exponentially as we proceed from one time step to another. In another word, Ein should be bounded as nâ†'âˆÅ¾ for a fixed âˆ†x and âˆ†t , or when both âˆ†x and âˆ†t tends to zero. There are two main methods to investigate whether a difference method is stable:Von-Neumann Fourier series expansion method and matrix method (Rao,2007).

**CHAPTER THREE**

**FTCS Method**

**Introduction **

Â FTCS method (Forward-Time Centered-Space) is one of some important and powerful finite difference methods used for solving parabolic partial differential equations that model many physics phenomenal, such that heat, mass transfer in fluids and unsteady behavior of fluid flow past bodies (Rao,2008).

Â For parabolic equations, FTCS approximates the partial derivative ut by the

(first-order Forward-Time) approximation and the partial derivative uxx by the

(second-order Centered-Space) approximation (Hoffman, 2001).

In this chapter, we will treat the FTCS method and other relevant subjects in more general details.

Â In our study and specifically in this chapter, FTCS method will be considered and used for solving one dimensional linear heat equations.

First of all and for more illustrations, let us consider a one-dimensional diffusion equation at the mesh grids points (i,n) as follow (Rao,2008):

âˆ‚uâˆ‚tÂ =Î±2Â âˆ‚u2âˆ‚x2 (3.1)

By introducing Forward difference approximation to the ut and Central difference approximation to the uxx , equation (3.1) becomes:

uxi,tn+âˆ†t-uxi,tnâˆ†t=Î±2uxi+âˆ†x,tn-2uxi,tn+uxi-âˆ†x,tnâˆ†x2 (3.2)

Notice that, uxi,tn+âˆ†t,uxi,tn, uxi+âˆ†x,tn and Â uxi-âˆ†x,tn in equation (3.2) are not actual values of u at Â uxi,tn+âˆ†t,uxi,tn,Â uxi+âˆ†x,tn

and uxi-âˆ†x,tn But approximations to u at these points.

So Eq(3.2) can be written as: Â

uin+1-uinâˆ†t=Î±2ui+1n-2uin+ui-1nâˆ†x2 (3.3)

So that;

uin+1=uin+Î±2âˆ†tâˆ†x2(ui+1n-2uin+ui-1n Â (3.4)

So we define

Â Î±2âˆ†tâˆ†x2=r Â Â Â (3.5)

Then, equation (3.5) becomes:

uin+1=uin+r(ui+1n-2uin+ui-1n) Â Â (3.6)

Eq(3.6) is the FTCS formula, since uin+1 can easily computed if all values of u at time level n are known (Rao, 2008).

**Stability analysis for FTCS: Fourier's method**

By formulating the Fourier forms component solutions in the form:

uin=Ï„ne-1Â iÎ¸ Â Â (3.7)

Into equation (3.6) and we assume that{Â -1 =I} , yields:

Ï„n+1eIÎ¸i=rÏ„neIÎ¸(i-1)+1-2rÏ„neIÎ¸i+rÏ„neIÎ¸(i+1) (3.8)

By dividing Eq (3.8) byÂ Ï„neIÎ¸i and rearrange the terms we get:

Â Ï„=1+reIÎ¸-2+e-IÎ¸ Â (3.9)

Ï„2=12-r+r(eIÎ¸-e-IÎ¸)2 Â (3.9.1)

Using the identity:

cosÎ¸=(eIÎ¸+e-IÎ¸)/2Â , Â 1-cosÎ¸=2sin2(Î¸/2) (3.10)

In Eq (3.9.1)Â and multiply Eq (3.9.1) by 2 gives:

Ï„=1-4rsin2(Î¸/2) Â Â (3.10)

For Fourier's method, a method is stable if and only if Ï„â‰¤1 for all values of Î¸.

FTCSÂ method to be stably,Â Fourier method requires:

1-4rsin2(Î¸/2)â‰¤1 âˆ€Î¸, so that -1â‰¤1-4rsin2(Î¸/2)â‰¤1. This inequality is satisfied if râ‰¤12 (Kincaid, Cheney, 1996; smith, 1985;Thomas,1999; LeVeqe, 2007).

**Stability analysis for FTCS: Matrix method**

**Concept of matrix**

In this field we consider five important concepts of matrix which accedes to prove the the stability requirement for the FTSC by using matrix method:

1st) the norm of the matrix A is a real positive number given a measure of the â€˜size'

Of the matrix and must satisfy the following requirements:

1) A>0 if Aâ‰ 0 and A=0 if A=0.

2) cA=cA for a real or complex scalar c.

3) ABâ‰¤AB.

4) A+Bâ‰¤A+B.

5) A1=maxiai,j.

6) AâˆÅ¾=maxjai,j.

7) A2=maxâ?¡{Î»,Î» an eigenvalue of ATA12}.

8) limxâ†'âˆÅ¾Ax=0 if Aâ‰¤1.

9) limxâ†'âˆÅ¾Ax=0 if and only if Î»i<1 for all eigenvalue.

2nd) let Î»i be an eigenvalue of the nÃ—n matrix A and Xi the corresponding eigenvector, hence

AXi=Î»iXi.

Therefore Î»iâ‰¤A for i=1â€¦,n hence Ï?(A)â‰¤A.

3rd) the eigenvalue of a common tridiagonal matrix:

The eigenvalue of nÃ—n matrix

a Â bÂ â‹¯â‹¯0Â c Â a Â b Â â‹® c â‹±â‹® â‹±â‹±â‹±â‹® Â â‹±â‹±b0 Â â‹¯ca

Are Î»i=a+2bc coskÏ€n+1for k=1,â€¦,n

Where a , b and c may be real or complex number.

4th) Gershgurian's first theorem:

In this theorem, Gershgurian proved that the largest of the moduli of the eigenvalue of a given square matrix A cannot exed the largest sum of the moduli of the elements along any row or any column, that is to say:

Ï?Aâ‰¤A1Â or Ï?(A)â‰¤AâˆÅ¾

5th) Gershgurian's second theorem:

Let pkbe the sum of the moduli of the elements along kth row excepting the diagonal elements ak,k , then each eigenvalue of matrix A lies inside or on the boundary of the as least one of the circles:

Î»-ak,k=pk (Smith, 1985).

So, we consider the following linear and non-homogeneous heat equation:

ut=uxx+f(x,t)

With the initial condition ux,0=g(x) and boundary conditions as follows:

Â u0,t=pt

ul,t=qt

The explicit FTCS method is used and obtains:

uin+1=sui-1n+1-2suin+sui+1n+sui+1n+kf(ih,nk)Â ,Â i=1,â€¦,M-1Â n=0,1,â€¦,N-1

Thus we have

u1n+1u2n+1â‹®â‹®â‹®uM-1n+1=1-2ss0â€¦â€¦0s1-2ss00â‹±â‹±â‹±â‹®â‹®â‹±â‹±â‹±â‹±â‹®â‹®â‹±â‹±s0â‹¯â‹¯â‹¯s1-2su1nu2nâ‹®â‹®â‹®uM-1n+su0n+kf(h,nk)kf(2h,nk)â‹®â‹®kf((M-1)h,nk)su0M+kf(Mh,nk)

Thus Un+1=AUn+B.

It is clear that the matrix A can be written as:

1-2ss0â€¦â€¦0s1-2ss00â‹±â‹±â‹±â‹®â‹®â‹±â‹±â‹±â‹±â‹®â‹®â‹±â‹±s0â‹¯â‹¯â‹¯s1-2s=

10â€¦â€¦00010â‹®â‹±â‹®â‹®â‹±â‹®01000â€¦â€¦01+s-21â€¦â€¦001-2â‹±0â‹®â‹±â‹±â‹±â‹®â‹®â‹±â‹±â‹±â‹®0â‹±-2100â€¦â€¦1-2

In other words,

A=IM-1+sTM-1

Where matrix IM-1is the identity matrix and matrix TM-1 is a tridiagonal matrix.

Thus, the eigenvalue of TM-1 can be obtained as the following:

Î»k=-4sin2kÏ€2M Â , Â k=0,1,â€¦,M-1.

So the equation will be stable where:

A2=max-4sin2kÏ€2Mâ‰¤1

Thus hâ†'0 , Mâ†'âˆÅ¾Â give usÂ 0<sâ‰¤12 (smith,1985;Rao, 2007).

**Consistency analysis for FTCS: Taylor's series**

To show that the FTCS method is consistent with the Eq (3.1) we substitute the exact solution of the equation ut=Î±2uxx into the formula (3.6), we abstain:

uxi,tn+1=uxi,tn+ruxi+1,tn-2uxi,tn+uxi-1,tn (3.10)

Expand uxi,tn+1 about xi,tn, we obtain:

uxi,tn+1=uxi,tn+âˆ†t=

uxi,tn+âˆ†tâˆ‚uâˆ‚tx=xi,t=tn+âˆ†t22!âˆ‚2uâˆ‚t2x=xi,t=tn+âˆ†t33!âˆ‚3uâˆ‚t3x=xi,t=tn+âˆ†t44!âˆ‚4uâˆ‚t4x=xi,t=tn+â€¦ Â (3.11)

Expand uxi+1,tn about xi,tn, we obtain:

uxi+1,tn=uxi+âˆ†x,tn=

uxi,tn+âˆ†xâˆ‚uâˆ‚xx=xi,t=tn+âˆ†x2!âˆ‚2uâˆ‚x2x=xi,t=tn+âˆ†x33!âˆ‚3uâˆ‚x3x=xi,t=tn+âˆ†x44!âˆ‚4uâˆ‚x4x=xi,t=tn+â€¦ Â (3.12)

Expand uxi-1,tn about xi,tn, we get:

uxi-1,tn=uxi-âˆ†x,tn=

uxi,tn-âˆ†xâˆ‚uâˆ‚xx=xi,t=tn+âˆ†x22!âˆ‚2uâˆ‚x2x=xi,t=tn-âˆ†x33!âˆ‚3uâˆ‚x3x=xi,t=tn+âˆ†x44!âˆ‚4uâˆ‚x4x=xi,t=tn-â€¦ Â (3.13)

If the expansions (3.11), (3.12) and (3.13) are substituted into Eq (3.10) and by considering that r=Î±2âˆ†tâˆ†x2, we obtain:

uxi,tn+âˆ†tâˆ‚uâˆ‚tx=xi,t=tn+âˆ†t22!âˆ‚2uâˆ‚t2x=xi,t=tn+âˆ†t33!âˆ‚3uâˆ‚t3x=xi,t=tn+âˆ†t44!âˆ‚4uâˆ‚t4x=xi,t=tn+â€¦=

uxi,tnÎ±2âˆ†tâˆ†x2uxi,tn+âˆ†xâˆ‚uâˆ‚xx=xi,t=tn+âˆ†x2!âˆ‚2uâˆ‚x2x=xi,t=tn+âˆ†x33!âˆ‚3uâˆ‚x3x=xi,t=tn+âˆ†x44!âˆ‚4uâˆ‚x4x=xi,t=tn+2uxi,tn+uxi,tn-âˆ†xâˆ‚uâˆ‚xx=xi,t=tn+âˆ†x22!âˆ‚2uâˆ‚x2x=xi,t=tn-âˆ†x33!âˆ‚3uâˆ‚x3x=xi,t=tn+âˆ†x44!âˆ‚4uâˆ‚x4x=xi,t=tn-â€¦

Then we get:

âˆ†tâˆ‚uâˆ‚tx=xi,t=tn+âˆ†t22!âˆ‚2uâˆ‚t2x=xi,t=tn+âˆ†t33!âˆ‚3uâˆ‚t3x=xi,t=tn+â€¦

=Î±2âˆ†tâˆ‚2uâˆ‚x2x=xi,t=tn+âˆ†x212âˆ‚4uâˆ‚x4x=xi,t=tn+â€¦ (3.14)

By dividing both sides by âˆ†t and rearrange the terms, then we get:

âˆ‚uâˆ‚tx=xi,t=tn+âˆ†t2!âˆ‚2uâˆ‚t2x=xi,t=tn+âˆ†t23!âˆ‚3uâˆ‚t3x=xi,t=tn+â€¦ =

Î±2âˆ‚2uâˆ‚x2x=xi,t=tn+âˆ†x212âˆ‚4uâˆ‚x4x=xi,t=tn+â€¦ Â (3.15)

As âˆ†tâ†'0 andÂ âˆ†xâ†'0 Eq (3.15)Â approaches ut=Î±2uxx , which is the diffusion equation, so the FTCS method for Eq(3.15) is consistent with ut=Î±2uxxÂ (Hoffman, 2001).

**Convergence**

The convergence of FTCS is related to the Consistency and stability, which is answered by study of both (Consistency & stability).

FTCS is said to be convergence if and only if FTCS is consistent and stable, otherwise it is not convergence (Hoffman, 2001;LeVeqe, 2007).

**Local truncation error**

Â By return to the equation (3.15), the FTCS formula for ut=Î±2uxx is:

Â Â uin+1-uinâˆ†t= Î±2ui+1n+2uin+ui-1nâˆ†x2 Â (3.16)

Let;

Fi,nu=uin+1-uinâˆ†t-Î±2ui+1n+2uin+ui-1nâˆ†x2=0 (3.17)

therefore;

Ti,n= Fi,nu=uin+1-uinâˆ†t-Î±2ui+1n+2uin+ui-1nâˆ†x2=0 Â (3.18)

By Taylor's expansions:

ui+1n=uxi+âˆ†x,tn= uxi,tn+âˆ†xâˆ‚uâˆ‚xx=xi,t=tn+âˆ†x2!âˆ‚2uâˆ‚x2x=xi,t=tn+âˆ†x33!âˆ‚3uâˆ‚x3x=xi,t=tn+âˆ†x44!âˆ‚4uâˆ‚x4x=xi,t=tn+â€¦ Â (3.19)

ui-1n=uxi-âˆ†x,tn=uxi,tn-âˆ†xâˆ‚uâˆ‚xx=xi,t=tn+âˆ†x22!âˆ‚2uâˆ‚x2x=xi,t=tn-âˆ†x33!âˆ‚3uâˆ‚x3x=xi,t=tn+âˆ†x44!âˆ‚4uâˆ‚x4x=xi,t=tn-â€¦ (3.20)

uin+1=uxi,tn+âˆ†t=uxi,tn+âˆ†tâˆ‚uâˆ‚tx=xi,t=tn+âˆ†t2!âˆ‚2uâˆ‚t2x=xi,t=tn+âˆ†t33!âˆ‚3uâˆ‚t3x=xi,t=tn+âˆ†t44!âˆ‚4uâˆ‚t4x=xi,t=tn+â€¦ (3.21)

Since ;

Â âˆ‚uâˆ‚t-Î±2âˆ‚u2âˆ‚x2= 0

Substituting Eq(s): (3.19),(3.20) and(3.21) in Eq(3.18), once we get:

Ti,n=âˆ‚uâˆ‚u-Î±2âˆ‚u2âˆ‚x2x=xi,t=tn+âˆ†t2âˆ‚2uâˆ‚t2x=xi,t=tn-Î±2âˆ†x212âˆ‚4uâˆ‚x4x=xi,t=tn+âˆ†t26âˆ‚3uâˆ‚t3x=xi,t=tn-Î±2âˆ†x4360âˆ‚6uâˆ‚x6x=xi,t=tn+â€¦ (3.22)

Therefore, the principal part of the local truncation error is:

âˆ†t2âˆ‚2uâˆ‚t2x=xi,t=tn-Î±2âˆ†x212âˆ‚4uâˆ‚x4x=xi,t=tn+â€¦ (3.23)

So, this is the local truncation error for the FTCS method for the diffusion equation at the points(iâˆ†x,nâˆ†t). And local truncation error for FTCS method for the diffusion equation is said to be 0âˆ†t+0(âˆ†x)2 (smith, 1985).

**Problem statement**

**First problem:**

Â Consider the second order one dimensional linear heat equation:

utx,t=Ïˆ(x,t)uxxx,t=x2tuxxx,t on 0,1Ã—0,T (3.24)

Subjected to the boundary conditions:

u0,t=0.

u1,t=expâ?¡(t2).

And the initial condition:

u(x,0)=fx=x2

The exact solution for this problem is given by (A.Soufyan; M.Boulmalf, 2005):

ux,t=x2expâ?¡(t2)

**F****TCS ****method for first problem:**

By applying the Forward difference approximation to ut(x,t) andÂ Central difference approximation toÂ uxx(x,t), so Eq(3.24) becomes:

uxi,tn+âˆ†t-u(xi,tn)âˆ†t=iâˆ†x2(nâˆ†t)uxi+âˆ†x,tn-2uxi,tn+uxi-âˆ†x,tnâˆ†x2

Multiplying both sides by âˆ†t(âˆ†x2) and rearrange the terms we obtain:

âˆ†x2uin+1-uin=iâˆ†x2nâˆ†t(âˆ†t)ui+1n-2uin+ui-1n

So:

âˆ†x2uin+1=âˆ†x2uin+iâˆ†x2nâˆ†t(âˆ†t)ui+1n-2uin+ui-1n

Dividing both sides by (âˆ†x2) yields:

uin+1=uin+âˆ†tâˆ†x2iâˆ†x2nâˆ†t(âˆ†t)ui+1n-2uin+ui-1n

Let r=âˆ†tâˆ†x2

So, the FTCS formula for equation 3.24 is:

uin+1=uin+riâˆ†x2nâˆ†tâˆ†tui+1n-2uin+ui-1n Â (3.25)

**Second ****problem:**

Â Consider the second order one dimensional linear heat equation:

utx,t=uxxx,t+Ï†x,t=uxxx,t+exp-x(cos t-sint) (3.26)

Subjected to the boundary conditions:

u0,t=sint

u1,t=sin(t/e)+1

And the initial condition:

ux,0=fx=x

The exact solution for this problem is given by (Bongsoo,2007):

ux,t=x+expâ?¡(-x)sint)

**F****TCS ****method for ****second**** problem:**

By applying the Forward difference approximation to ut(x,t) andÂ Central difference approximation toÂ uxx(x,t), ),, so Eq(3.26) becomes:

uxi,tn+âˆ†t-u(xi,tn)âˆ†t=uxi+âˆ†x,tn-2uxi,tn+uxi-âˆ†x,tnâˆ†x2+e-iâˆ†x(cosnâˆ†t)-sinnâˆ†t)

Multiplying both sides by âˆ†t(âˆ†x2) and rearrange the terms and by assuming that

Î¼=e-iâˆ†x(cosnâˆ†t)-sinnâˆ†t) so we get:

âˆ†x2uin+1-uin=âˆ†tui+1n-2uin+ui-1n+(âˆ†tâˆ†x2)Î¼

So it becomes:

âˆ†x2uin+1=âˆ†x2uin+âˆ†tui+1n-2uin+ui-1n+(âˆ†tâˆ†x2)Î¼

By solving the equivalent expression by uin+1 we get:

uin+1=uin+âˆ†tâˆ†x2ui+1n-2âˆ†tâˆ†x2uin+âˆ†tâˆ†x2ui-1n+(âˆ†t)Î¼

Let r=âˆ†tâˆ†x2Â so:

uin+1=uin+rui+1n-2ruin+rui-1n+(âˆ†t)Î¼

Substituting the value of Î¼ and rearrange the terms once we get:

uin+1=rui-1n+1-2ruin+rui+1n+âˆ†te-iâˆ†x(cosnâˆ†t-sinnâˆ†t) (3.27)

Which is the FTCS formula for equation 3.26.

**Third ****problem:**

Consider the second order one dimensional linear heat equation:

utx,t=uxxx,t+Ï•x,t=uxxx,t+exp-x(cosh t-sinht) Â (3.28)

Subjected to the boundary conditions:

u0,t=sinht

u1,t=e.sinht+t+1/6

And the initial condition:

ux,0=fx=x3/6

The exact solution for this problem is given by (Bongsoo,2007):

ux,t=expâ?¡(x)sinht+x3/6+xt

**F****TCS ****method for ****third**** problem:**

Â By applying the Forward difference approximation to ut(x,t) andÂ Central difference approximation toÂ uxx(x,t), so Eq(3.27) becomes:

uxi,tn+âˆ†t-u(xi,tn)âˆ†t=uxi+âˆ†x,tn-2uxi,tn+uxi-âˆ†x,tnâˆ†x2+e-iâˆ†x(coshnâˆ†t)-sinhnâˆ†t)

Let Î²=e-iâˆ†x(coshnâˆ†t)-sinhnâˆ†t)

Multiplying both sides by âˆ†t(âˆ†x2) and rearrange the terms and by assuming that

Î¼=e-iâˆ†x(cosnâˆ†t)-sinnâˆ†t) so we get:

âˆ†x2uin+1-uin=âˆ†tui+1n-2uin+ui-1n+(âˆ†tâˆ†x2)Î²

So it becomes:

âˆ†x2uin+1=âˆ†x2uin+âˆ†tui+1n-2uin+ui-1n+(âˆ†tâˆ†x2)Î²

By solving the equivalent expression by uin+1 we get:

uin+1=uin+âˆ†tâˆ†x2ui+1n-2âˆ†tâˆ†x2uin+âˆ†tâˆ†x2ui-1n+(âˆ†t)Î²

Let r=âˆ†tâˆ†x2Â so:

uin+1=uin+rui+1n-2ruin+rui-1n+(âˆ†t)Î²

Substituting the value of Î² and rearrange the terms ,the FTCS formula for Eq(3.28):

uin+1=rui-1n+1-2ruin+rui+1n+âˆ†te-iâˆ†x(coshnâˆ†t-sinhnâˆ†t) Â (3.29)

**CHAPTER FOUR**

**ADM Method**

**Introduction**

Â Over the last 10 years, many numerical methods are aimed to solve non-linear partial differential equations have appeared in many researcher's studies. However, most of these methods require a tedious analysis work which makes unnecessary linearization, perturbations, restrictive methods and assumptions which may change the problem itself being solved sometimes. Adomian decomposition method (ADM) overcame the major part of these difficulties were suffered and it was successfully applied in many applications in applied sciences (JAVIDI, M. G., A, 2007; LUO, X.-G, 2005).

In our study, Adomian decomposition method (ADM) will be considered and used for solving one dimensional linear heat equation.

**Literature survey**

Â In the early 1980s, George Adomian developed a fruitful powerful method for the computational of exact solution for linear and non-linear partial differential equations (pde's) (Adomian, G, 1984; Adomian, G, 1988; Adomian, G, 1990; Adomian, G, 1998).

A wide class of problems in mathematics, physics, engineering, biology, chemistry and in other sciences are presented and solved using Adomian decomposition method by many authors (Adomian, G, 1984; Adomian, G, 1998; Wazwaz, A.-M. E.-S., Salah, M, 2001; Wazwaz, A.-M, 2000; Chen, W. L., Zhengyi, 2004; Jang, B,2007; Amani, R., A.Sadeghi,J ,2007).

Furthermore, the method provides an effective tools for analytical solutions of a wide class of dynamical systems representing real physical and engineering problems, which generally converge very rapidly in order to reach an accurate solution, is needs less calculations to solve parabolic partial differential equations comparison with other methods (Adomian, G,1984; Adomian, G,1990; Soufyane, A. B., M,2005; Ray, S., S.Bera,K,R,2005; Bildik, N. B., Hatice,2005).

**Analysis of the method**

In this section and for convenience of the readers, we review the basic principles of the standard Adomian decomposition method by the following procedures:

Â We consider the following differential equation (G.Adomian, 1988):

Lu+Ru+Nu=g(x,t) Â (4.1)

Where L represents the highest order derivative which assumed to be invertible, R is the linear differential operator of order less than L, where g(x,t) is the source term, and finally, Nu represents the non-linear terms[ Â ].

Solving Lu from Eq(4.1), we obtain[ ]:

Lu=gx,t-Ru-Nu Â (4.2)

Since L is invertible, so by applying the inverse operator L-1 to both sides of equation (4,2) yields[ Â ]:

L-1Lu=L-1(gx,t)-L-1(Ru)-L-1(Nu) Â (4.3)

IfÂ L-1 is second order operator, so L-1 is a two-fold integration operator assumed to be stated from (to=0)Â toÂ t Â [ ].

Therefore, Lt , Lx are differential operators defined as[ ]:

Lt=âˆ‚2âˆ‚t2Â , Lx=âˆ‚2âˆ‚x2Â and L-1=0t0t.dtdt .

Solving Eq(4.3) and using the given conditions we obtain:

u=Ï‰+L-1(gx,t)-L-1(Ru)-L-1(Nu)Â Â (4.4)

Â Where the function (Ï‰=A+Bt) represents all the terms arising from the given conditions, all the conditions are assumed to be prescribed for a problem aimed to be solved[ Â ].

Â For non-linear equations, the non-linear term Nu=F(u) will be decomposed by an infinite series of the so-called Adomian polynomials[ ], so:

Fu=n=0âˆÅ¾An

Where An=Â 1n!dndÎ»nNn=0âˆÅ¾unÎ»nÎ»=0 Â [Â 5]

Where An are generated for all kinds of nonlinearity, so Ao depends on u0 only, A1 depends on u0 and u1 ,and so on until An depends on u0 , u1 , â€¦â€¦,un [ Â ].

For more details about formulating Adomian polynomials, specific algorithms were set in [ ].

The standard Adomian decomposition method decomposes the solution u(x,t) into an infinite series defined by[ ]:

ux,t=n=0âˆÅ¾un(x,t) [ ] (4.5)

Where u0 , u1 , â€¦â€¦,un are usually determined as[ Â ]:

u0=Ïˆ+L-1(gx,t)

u1=-L-1Ru0-L-1(A0)

:

:

un+1=-L-1Run-L-1(An) nâ‰¥0 (4.6)

Substituting the values of An in Eq(4.6) leads to the determination of the solution of u in series form defined by Eq(4.5) immediately[ Â ].

So, the approximate solution using n-terms approximation to u(x) is given by Ï†n such that:

ux,t=limnâ†'âˆÅ¾Ï†n

Where Ï†n is defined as:

Ï†n =k=0n-1uk(x,t)

Which can be written as:

Ï†1 =u0,

Ï†2 =u0+u1,

Ï†3 =u0+u1+u2,

Ï†4 =u0+u1+u2+u3,

:

:

Ï†n =u0+u1+u2+u3+â€¦+nn-1.

So , the absolute error is defined as:

absolute error =exact solution-approximate solution

[ 21Â 12 6 ]

Abdul-Majid Wazwaz have investigate the existence of the so-called noise terms 'selfâ€“cancelling' phenomenon, which if it exists, it gives a useful tool to fast convergence of the solution by the first two iterations only (u0 & u1). It is noted that the noise terms phenomenon may appear only in non-homogenous problems, whereas it is not exist in homogenous problems. The remaining â€œnon-cancelingâ€? terms may give

### Cite This Dissertation

To export a reference to this article please select a referencing stye below: