## Explicit Methods for Reaction AdvectionDiffusion Equations

The equation whose solution we will focus on is u — D(u, x,y)q2xu + g(u, x, t)dxu + f (u, x, t), u(x, to) — «o(x) [1]

where t and x are the time and space variables and u = 8t«. (A generalization of [1] to two spatial dimensions will be considered later.) The three terms on the RHS of [1] are associated with diffusion, advection, and reaction, respectively. To illustrate the main issues ofnumerical solution of [1], we first consider its simpler version - the heat equation - where we retain only the diffusion term and take the diffusion coefficient D to be a constant:

The last two equations are called boundary conditions; such (or similar) conditions must always be supplied for a partial differential equation (PDE). Boundary conditions that fix the value of the unknown, as in [3], are called Dirichlet boundary conditions. Boundary condition of other types

 ( n +1 r A r J n m ) V -1 m J V m ) + 1

Figure 1 The stencil of method [5].

will be considered later. Note that setting the boundary conditions [3] to zero does not restrict the generality of our treatment. Indeed, for the more general conditions, «(0) = a(t), «(X) = b(t), one defines a new variable w = « — [(X — x)a(t) + xb(t)]/X, which satisfies [2] and [3].

We discretize tn = t0 + nr and xm = 0 + mh, where r and h are the time step and the mesh size, and n = 0, 1, 2,... and m = 0, 1, ..., X/h=M. Denote «nmx «(xm, tn) and Um to be, respectively, the exact and the numerical solutions at node (xm, tn). Replacing

one obtains the finite-difference approximation to [2]:

Un+1 = rUl+i + (1 - 2r )un + rUl -1, un = 0, UM = 0 [5]

where

The meaning of the notation O(rp) is explained in Numerical Methods for Local Models. Equation [5] is the counterpart of the explicit Euler method for ordinary differential equations (ODEs). As follows from [4], its accuracy is 0(t) + O(h2). (However, for an initial condition uo(x) that is either discontinuous or has a discontinuous slope, the spatial accuracy is reduced: the smaller the r in [6], the closer the accuracy to the maximum possible value O(h ).) The numerical solution at node (xm, tn + j) can thus be found from [5] and the initial condition lfm — u0(xm) if one knows the solution at nodes

±!, tn). These four nodes form a so-called stencil for scheme [5], as shown schematically in Figure 1.