MATH3203 Lecture 3 PDEs in Biology and Classification of PDEs

Document Sample
MATH3203 Lecture 3 PDEs in Biology and Classification of PDEs Powered By Docstoc
					                  MATH3203 Lecture 3
        PDEs in Biology and Classification of PDEs
                                  Dion Weatherley
                      Earth Systems Science Computational Centre,
                                University of Queensland

                                    March 6, 2006


1 Recap of previous lecture                                                                                                2
  1.1 PDE form of the fundamental conservation law . . . . . . . . . . . . . . .                                           2
  1.2 Advection and diffusion equations . . . . . . . . . . . . . . . . . . . . . . .                                       2

2 PDE models in Biology                                                                                                    3
  2.1 Refinement of the population model . . . . . . . . . . . . . . . . . . . . . .                                        3
  2.2 Other PDE models in the life sciences . . . . . . . . . . . . . . . . . . . . .                                      4

3 Types of BCs and PDE classification                                                                                       4
  3.1 Types of PDE Boundary Conditions . . . . .           .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   4
      3.1.1 Dirichlet Boundary Conditions . . . . .        .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   4
      3.1.2 Neumann Boundary Conditions . . . .            .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   5
      3.1.3 Radiation Boundary Conditions . . . .          .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   5
  3.2 Classification of PDEs . . . . . . . . . . . . .      .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   5
  3.3 Characteristic coordinates and canonical forms       .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   6
      3.3.1 Hyperbolic case, D > 0 . . . . . . . . .       .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   7
      3.3.2 Parabolic case, D = 0 . . . . . . . . . .      .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   7
      3.3.3 Elliptic case, D < 0 . . . . . . . . . . .     .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   .   7

1     Recap of previous lecture
1.1     PDE form of the fundamental conservation law
Commencing with the principle of conservation of a quantity within a given domain, we
derived the PDE formulation of the fundamental conservation law:

                                         εt = −         .φ + σ.                         (1)
where ε(x, t) represents the density of the conserved quantity for each point in the domain
(x ∈ Ω) at any time, t. The flux of the conserved quantity through any given point is given
by the flux vector, φ(x, t) and the density of sources or sinks of the conserved quantity is
given by σ(x, t).
    Although we derived the fundamental conservation law assuming the conserved quan-
tity was energy, the same form applies for any conserved quantity. It is useful when
analysing a particular PDE, to consider it as a conservation law and identify the terms
playing the roles of density, flux and sources.

1.2     Advection and diffusion equations
The fundamental conservation law is a single equation in three unknowns, a so-called
underconstrained problem. In order to specify a properly constrained PDE model, we
must introduce constitutive relations relating the flux or sources to the density. Constitu-
tive relations arise from assumptions about the physical properties of the medium within
the domain. We examined two constitutive laws for the flux: advection and diffusion
constitutive laws.
    The advection constitutive law assumes that the conserved quantity is carried by an
underlying velocity field, v(x) that is independant of the density of the conserved quantity.
The flux term for advection is given by:

                                        φ(x, t) = ε(x, t)v(x)                           (2)
It is common to further assume that the flow is incompressible which is expressed as
  .v = 0, leading to the incompressible flow advection equation:

                                          εt + v.           ε=0                         (3)
    By solving the 1D advection equation, we demonstrated that advection models unde-
formed transport of the conserved quantity through the domain.
    The second constitutive law we examined was the diffusion equation, modelling the
flux of the conserved quantity from regions of high density to regions of lower density. In
this case the flux is proportional to the spatial gradient of the density:

                                   φ(x, t) = −D              ε(x, t)                    (4)
   When the constitutive law is substituted into the fundamental conservation law, we
obtain the diffusion equation:

                                          εt − D            ε=0                         (5)
where        is the Laplace operator:

                                    2        δ2  δ2  δ2
                                        =       + 2+ 2                                 (6)
                                            δx2 δy  δz
    Constitutive laws may be combined to model more complex situations. To model
both advection and diffusion, we can sum the appropriate constitutive laws for the flux,
to obtain the advection-diffusion equation:

                           εt − D           ε+ε       .v + v.   ε=0                    (7)

2     PDE models in Biology
In the physical sciences and engineering, mathematical models are usually expressed in
terms of nearly exactly measurable quantities and the goal is to obtain precise quanti-
tative results to be compared with experimental observations. In the biological sciences
however, it is rare that the phenomenon to be studied can be so precisely quantified,
due to the complexity of the systems to be considered. It is difficult to account for all
the intracacies of predator-prey interactions, tumor growth or the spread of infectious
diseases for example. Typically biological models are more qualitative, aiming to predict
only quanitative features rather than yielding detailed quantitative results. This section
aims to introduce some standard PDE models in the biological sciences, as a complement
to our derivation of advection and diffusion equations from physical principles.

2.1    Refinement of the population model
In the first lecture, we studied the Malthus and Logistic models for population growth.
Both such models were ODE models, involving only one independant variable, time.
We introduced the physical variable U (t) to represent the size of the population. It is
implicitly understood however, that this variable represents the population density (or
population per unit volume).
    While both these models capture some aspects of population growth, both ignore the
spatial distribution of an organism within its ecosystem. Even for an ecosystem with
only a single species, population pressure will cause the organisms to migrate apart, to
minimise the impact of competition for resources. To attempt to model this, we must
extend our definition of the population density to U (x, t), representing the population
density within a small volume surrounding x at time t. Implicit in the model is the
assumption that there are sufficient organisms that defining a local population density is
    To model population changes in both space and time, we must introduce a PDE
population growth model that includes terms governing the motion of organisms within
the ecosystem. One of the fundamental PDE models in biology is Fisher’s equation.
This equation assumes that the population will spread out spatially within a 1D domain.
Fisher’s equation is:
                       Ut = DUxx + rU 1 −        , U = U (x, t)                     (8)
By comparison with the fundamental conservation law, we can identify this as a diffusion
equation, with an additional source term modelling population growth due to reproduc-

2.2     Other PDE models in the life sciences
In the biological sciences, it is not always necessary to model the spatial structure however
one might be interested in modelling situations with other types of structure. For example,
we might be interested in modelling the ages of the organisms. Then the population
density will depend on both the age (a) and the time (t) i.e. U = U (a, t). In life history
theory, models are developed that predict the age structure of a population at any time
t given the initial age structure U (a, 0) and some conditions on birth and death rates,
related to age. For example, the birth rate might be higher for times when the average
age is lower and conversely the death rate might be higher for older populations.
    Other structured models may include independant variables such as size, weight etc.;
such models are known as physiologically structured models. Other common PDE models
in the life sciences include models for the biochemical reactions and motion of chemical
species in the blood, cells or other media. Advection is most commonly used to model
transport of particles, chemicals or organisms, via bulk motion of the supporting medium
e.g. wind, water, blood etc. Diffusion models the random motions of particles etc. that
cause these entities to disperse from high concentrations to lower concentrations.
    Sources and sinks depend upon the particular model. For population models the
source term represents birth or growth rates, and a sink represents the death rate either
by natural causes or by predation. For chemical problems, the source/sink term is called
a reaction term and represents the rate at which a chemical species is created and/or
consumed in chemical reactions, or lost by absorption (e.g. across a cell boundary).

3       Types of BCs and PDE classification
3.1     Types of PDE Boundary Conditions
As previously mentioned, to completely define a PDE model, the PDE equation itself is
insufficient. We also need to specify initial conditions that give the density at t = 0 for
all points x in the domain Ω i.e. ε(x, 0) and boundary conditions specifying conditions
on the density or flux at the boundaries of the domain (δΩ) for all time i.e. ε(xb , t) or
φ(xb , t) for xb ∈ δΩ.
    Typically the initial condition is simply some given function of space:

                                  ε(x, 0) = f (x),    x∈Ω                                (9)
   There are three types of boundary conditions that often occur in physical problems:
Dirichlet, Neumann and radiation boundary conditions. It is important to note that
the boundary conditions for a given problem must be consistent with the specified PDE
equation within the domain. For example, if we have a diffusion equation within the
domain, it is invalid to specify an advection-like boundary condition and vice versa.

3.1.1    Dirichlet Boundary Conditions
The first type of boundary conditions, Dirichlet BCs, provide a function of time (g(t)) as
the constraint on the density at a specified set of boundary points, xb ∈ (D ⊂ δΩ):

                                  ε(xb , t) = g(t),   t>0                               (10)

Suppose we are modelling heat flow within the domain. The Dirichlet BC asserts that
the temperature is constant along the specified set of boundary points xb .

3.1.2   Neumann Boundary Conditions
The second type of BC is Neumann BCs, in which the flux vector is specified by a given
vector function (h(t)) along a subset of the domain boundary (N ):

                         φ(xb , t) = h(t),   t > 0, xb ∈ (N ⊂ δΩ)                     (11)
The Neumann BC specifies the influx/outflux of the conserved quantity into/out of the
domain from the surrounding environment. If h(t) = 0 then we say that the boundary
is insulated ; there is no flow across the boundary and the domain is isolated from the
surrounding environment.

3.1.3   Radiation Boundary Conditions
A third type of boundary condition of use for diffusive PDEs is the radiation or Robin

                          φ(xb , t) = β (ε(xb , t) − ψ(t)) ,   t>0                    (12)
In heat flow, for example, this law expresses Newton’s law of cooling, which states that the
heat flux is proportional to the temperature difference between the domain boundary and
the temperature of the surrounding environment, ψ(t). β is a heat-loss factor defining the
thermal coupling of the domain with the surrounding environment. In this case, one may
wish to specify temperature profile for the environment which evolves in time or which is
a constant (ψ(t) =const.); in such situations we refer to the environment as a heat bath.

3.2     Classification of PDEs
We classify ODEs in terms of their order and whether they are linear or nonlinear. PDE
models are more difficult to classify, because of the greater variety of basic forms (or
structures) the PDE may take. Not only are the order and linearity important, but also
the PDE structure. It is the PDE structure that dictates what types of boundary and
initial conditions can be imposed and ultimately what types of physical processes the
PDE models.
    Let us consider a general second-order PDE in two independant variables. We can
write such an equation in the form:

                       Auxx + Buxt + Cutt + F (x, t, u, ux , ut ) = 0                 (13)
where A, B and C are constants. The classification of this PDE is based upon the principal
part of the equation, Lu ≡ Auxx + Buxt + Cutt . Here we have used x and t but this
discussion applies equally well for two spatial variables also e.g. x and y.
    The classification is based upon the sign of the quantity, D ≡ B 2 − 4AC which is
called the discriminant. If D > 0, the equation is classified as hyperbolic; if D < 0 the
equation is elliptic; and if D = 0 the equation is parabolic. This terminology comes from
the classification of plane curves; for example Ax2 + Ct2 = 1, where A, C > 0 and B = 0
(i.e. D < 0) graphs an ellipse in the xt-plane. Similarly, Ax2 − Ct2 = 1 graphs as a

    Under this classification, the diffusion equation is parabolic, Laplace’s equation is
elliptic and the 1D wave equation (uxx = c2 utt ) is hyperbolic. Note that the advection
equation does not have a classification under this scheme; advection is a first-order PDE.
As it turns out, all parabolic equations are diffusion-like, all hyperbolic equations are
wave-like, and all elliptic equations are static.

3.3    Characteristic coordinates and canonical forms
Now we show that the principal part Lu can be simplified for all three types of PDEs
by introducing a new set of independant variables. We seek a linear transformation of
coordinates that simplifies Lu:

                                 ξ = ax + bt,       τ = cx + dt                        (14)
where ξ and τ are new independant variables, and a, b, c, d are to be determined. So that
we can transform back to the original coordinates, we require the transformation to be
invertible i.e. ad − bc = 0.
    The dependent function u in the new variables will be denoted by U = U (ξ, τ ); that
is, u(x, t) = U (ax + bt, cx + dt). Then using the chain rule for derivatives we can write
the first-order partial derivatives of u as:

                              ux = Uξ ξx + Uτ τx = aUξ + cUτ
                              ut = Uξ ξt + Uτ τt = bUξ + dUτ

Applying the chain rule again, we obtain relations for the second-order partials:

                          uxx = a2 Uξξ + 2acUξτ + c2 Uτ τ
                           utt = b2 Uξξ + 2bdUξτ + d2 Uτ τ
                          uxt = abUξξ + (ad + cb)Uξτ + cdUτ τ

Substituting these quantities into the principal part and collecting like terms, we get

               Auxx + Buxt + Cutt = (Aa2 + Bab + Cb2 )Uξξ
                                    +(2acA + B(ad + bc) + 2Cbd)Uξτ
                                    +(Ac2 + Bcd + Cd2 )Uτ τ

Now we can select a, b, c, d so that some of the partials in the new variables disappear. We
have some flexibility so we choose a = c = 1. Then the coefficients of Uξξ and Uτ τ become
quadratic expressions in b and d respectively and we can set these to zero to obtain values
for b and d that eliminate Uξξ and Uτ τ from the equation e.g.

                                     A + Bb + Cb2 = 0                                  (15)
has two possible solutions:
                                         −B ± D
                                      b=                                               (16)

To ensure the transformation is invertible, we need (d − b) = 0 so we can select the
following values:
                                  √                      √
                           −B + D                  −B − D
                       b=                and d =                                (17)
                               2C                     2C

3.3.1   Hyperbolic case, D > 0
For the hyperbolic equation where D > 0, we can use the following transformation to
characteristic coordinates:
                                        √                  √
                                  −B + D              −B − D
                         ξ =x+               tτ = x +        t                 (18)
                                     2C                 2C
to simplify the equation to its canonical form:

                                Uξτ + G(ξ, τ, U, Uξ , Uτ ) = 0                        (19)
where only the mixed second-order partial derivative appears. This is a significant simpli-
fication over the original where all second-order derivatives occurred. Note that if C = 0
then these coordinates are not valid. We can however find a different set of characteristic
coordinates by asserting that b = d = 1 instead.

3.3.2   Parabolic case, D = 0
For the parabolic case, Equations 17 give b = d which results in a non-invertible transfor-
mation. After some trial and error, we can find a suitable transformation that eliminates
Uτ τ and Uξτ , namely:
                                   ξ = x,   τ =x−   t                       (20)
i.e. a = c = 1 and b = 0, leaving d = B/2C. Therefore, the canonical form in the
parabolic case becomes a diffusion-like equation:

                                Uξξ + H(ξ, τ, U, Uξ , Uτ ) = 0                        (21)

3.3.3   Elliptic case, D < 0
Equations 17 now yield values for b and d that are complex conjugates i.e. d = b. Selecting
a = c = 1 again, we get a complex transformation:

                                 ξ = x + bt,       τ = x + bt.                        (22)
This is a little nasty but it can be avoided by making a second transformation to real
variables α and β:
                                   1                 1
                              α = (ξ + τ ), β = (ξ − τ ).                             (23)
                                   2                 2i
After a bit of routine algebra, the elliptic equation reduces to:

                           Uαα + Uββ + K(α, β, U, Uα , Uβ ) = 0                       (24)

where only the mixed partial is missing. Recognise that the combination of second-order
partials is the Laplace Operator, so we have a nonhomogeneous Laplace equation for
K = 0 or simply the homogeneous Laplace equation for K = 0. Generally transformation
to characteristic coordinates is only of practical use for hyperbolic equations.



Shared By:
Description: MATH3203 Lecture 3 PDEs in Biology and Classification of PDEs