MATH661

Lecture 27: Differential Conservation Laws

1.The relevance of physics for scientific computation

Efficient algorithms often arise from the specifities of an underlying application domain, perhaps none more so than those inspired from physics. Classical physics can be derived from a remarkably small set of experimentally verified postulates.

The least action principle asserts that a physical system can be described by a function $L (t, q, \dot{q})$ of the system generalized coordinates $q (t)$ and velocities $\dot{q} (t) = d q / d t$ , known as the Lagrangian, itself the difference of the system's kinetic and potential energy $L = K - U$ . The time evolution of the system is known as the system's trajectory $(q (t), \dot{q} (t))$ , and of all possible trajectories consistent with system constraints the trajectory actually followed by the system from initial time $t_{0}$ to final time $t_{1}$ minimizes a functional known as the action $S$
$S (q, \dot{q}) = \int_{t_{0}}^{t_{1}} L (t, q (t), \dot{q} (t)) d t .$
Example. However complex a physical system might be, application of the least action principle follows the procedure exemplified here for a simple mass-spring system. A point mass $m$ attached to a spring of stiffness $k$ is at distance $q (t)$ away from the equilibrium position $q = 0$ . For constant $m$ , this harmonic oscillator motion is described by the differential system
$\frac{d}{d t} (m v) = - k q \Rightarrow \frac{d v}{d t} = - \frac{k}{m} q, \frac{d q}{d t} = v .$
A state of this system is given by the values for position and velocity $(q, v)$ , and the above equations specify the time evolution of the system. Denoting the velocity as $v = \dot{q}$ , $d v / d t = \ddot{q}$ , and eliminating $v$ gives the familiar

$m \ddot{q} + k q = 0 .$ (1)

The same equation also results from the minimization of the action $S (q, \dot{q})$ of the Lagrangian

$L (q, \dot{q}) = \frac{1}{2} m {\dot{q}}^{2} - \frac{1}{2} k q^{2},$ (2)

with $K = m {\dot{q}}^{2} / 2$ , $U = k q^{2} / 2$ . The minimization is performed over all trajectories $(q (t), \dot{q} (t))$ with the same end-point values at $t_{0}, t_{1}$ . Let the $δ$ operator denote a small change in a trajectory. Since all trajectories have the same endpoints $δ q (t_{0}) = δ q (t_{1}) = 0$ . The change in the action is
$δ S = \int_{t_{0}}^{t_{1}} δ L (t, q (t), \dot{q} (t)) d t = \int_{t_{0}}^{t_{1}} (\frac{\partial L}{\partial q} δ q + \frac{\partial L}{\partial \dot{q}} δ \dot{q}) d t .$
Consider changes in overall trajectory to be independent of time such that the $δ$ and $d / d t$ operators commute, and apply integration by parts
$\int_{t_{0}}^{t_{1}} [\frac{\partial L}{\partial \dot{q}} δ (\frac{d q}{d t})] d t = \int_{t_{0}}^{t_{1}} [\frac{\partial L}{\partial \dot{q}} \frac{d}{d t} (δ q)] d t = {[\frac{\partial L}{\partial \dot{q}} δ q]}_{t_{0}}^{t_{1}} - \int_{t_{0}}^{t_{1}} [δ q \frac{d}{d t} (\frac{\partial L}{\partial \dot{q}})] d t = - \int_{t_{0}}^{t_{1}} [δ q \frac{d}{d t} (\frac{\partial L}{\partial \dot{q}})] d t .$
For $S$ to be at a minimum the change in the action must be stationary $δ S = 0$ ,
$\int_{t_{0}}^{t_{1}} [\frac{\partial L}{\partial q} - \frac{d}{d t} (\frac{\partial L}{\partial \dot{q}})] δ q d t = 0 .$
For the above to be valid for all $δ q$ the equation
$\frac{\partial L}{\partial q} - \frac{d}{d t} (\frac{\partial L}{\partial \dot{q}}) = 0$
must hold, and replacing (2) recovers (1).
One class of constraints are the conservation laws, the experimental observation that certain quantities remain constant during the system's evolution. Classical mechanics identifies three conserved quantities: mass, momentum, and energy. It is a matter of personal preference to also consider conservation of angular momentum as fundamental or as a consequence of conservation of linear momentum. Classical electrodynamics adds conservation of electric charge, while quantum mechanics also defines conservation of certain microscopic quantities known as quantum numbers such as the baryon or lepton numbers.

Other constraints often refer to allowed spatial positions and are known generically as geometric constraints. Note that this is an idealization: in reality some other physical system $M$ is interacting with the one being considered $P$ and it is assumed that the system $M$ is so much larger that its position does not change. Such idealization or modeling assumptions are often encountered. As another important example, the system $P$ may exhibit energy dissipation, such as the decrease of an object's momentum due to friction. Energy is indeed lost from system $P$ to the surrounding medium $M$ , but the overall energy of the combined system $M + P$ is conserved.

Scientific computation uses many concepts and terms from physics, such as the characterization of a numerical scheme for differential equations as being “conservative”, in the sense of maintaining the conserved physical quantities. Also, a remarkably large number of efficient algorithms arise from the desire to mimic physical properties. Conservation laws can be stated for a small enough spatial domain that it can be considered to be infinitesimal in the sense of calculus. In this case differential conservation laws are obtained. Alternatively, consideration of a finite-sized spatial domain leads to integral formulation of the conservation laws. In a large class of physical systems of current research interest models are constructed in which the evolution of the system $P$ depends on the history of interactions with the surrounding medium $M$ . Such systems are described by integro-differential laws, elegantly expressed through fractional derivatives.

2.Conservation laws

Banking example.

Conservation of some physical quantity is stated for a hypothesized isolated system. In reality no system is truly isolated and the most interesting applications come about from the study of interaction between two or more systems. This leads to the question of how one can follow the changes in physical quantities of the separate systems. An extremely useful procedure is to set up an accounting procedure. A mundane but illuminating example is quantity of Euros

E

in a building

B

. If the building is a commonplace one, it is to be expected that when completely isolated, the amount of currency in the building is fixed

E = E_{0} .

(3)

$E_{0}$ is some constant. Equation (3) is self-evident but not particularly illuminating – of course the amount of money is constant if nothing goes in or out! Similar physics statements such as “the total mass-energy of the universe is constant” are again not terribly useful, though one should note this particular statement is not obviously true. Things get more interesting when we consider a more realistic scenario in which the system is not isolated. People might be coming and going from building $B$ and some might actually have money in their pockets. In more leisurely economic times, one might be interested just in the amount of money in the building at the end of the day. Just a bit of thought leads to

E_{n} = E_{n - 1} + Δ E_{n - 1, n}

where $E_{n}$ is the amount of money at the end of day $n$ , $E_{n - 1}$ that from the previous day and $Δ E_{n - 1, n}$ the difference between money received and that payed in the building during day $n$

Δ E_{n - 1, n} = R_{n - 1, n} - P_{n - 1, n} .

As economic activity picks up and we take building $B$ to mean “bank” it becomes important to keep track of the money at all times, not just at the end of the day. It then makes sense to think of the rate at which money is moving in or out of the building so we can not only track the amount of currency at any given time, but also be able to make future predictions. Instead of separate receipts $R$ and payments $P$ , use a single quantity $F$ to denote the amount of money leaving or entering building $B$ during time interval $Δ t$ with the understanding that positive values of $F$ represent incomes and negative ones expenditures. Such understandings go by the name of sign conventions. They're not especially meaningful but it aids communication if a single convention is adopted. The amount of euros in the building then changes in accordance to

E (t + Δ t) = E (t) + F Δ t,

(4)

and $F$ is known as a flux, the Latin term for flow.

While (4) is a good approximation for small intervals, errors arise when $Δ t$ is large since economic acitvity might change from hour to hour. Better accounting is obtained by considering $F$ as defined at any given time $t$ , such that $F (t)$ is the instantaneous flux of euros at time $t$ . The fundamental theorem of calculus then states

E (t + Δ t) = E (t) + \int_{t}^{t + Δ t} F (τ) d τ,

(5)

with the same significance as (4).

In a large bank one keeps track of the amount of money in individual rooms and the inflows and outflows through individual doors. A room or door can be identified by its spatial position $𝒙 = (x_{1}, x_{2}, x_{3})$ , but $𝒙$ refers to a single point and physical currency occupies some space. The conceptual difficulty is overcome by introducing a fictitious density of currency at time $t$ denoted by $e (x, t)$ . The only real meaning associated with this density is that the sum of all values of $e (x, t)$ in some volume $ω$ is the amount of currency in that volume

E (ω, t) = \int_{ω} e (𝒙, t) d 𝒙 .

(6)

On afterthought, the same sort of question should have arisen when $E (t)$ was defined at one instant in time. Ingrained psychological perspectives make $E (t)$ more plausible, but were we to live our lives such that quantum fluctuations are observable, $E (t)$ would be much more questionable.

By an analogous procedure, define $𝒇 (𝒙, τ)$ as the instantaneous flux density of euros in a small region around $(𝒙, τ)$ . This flux is a vector quantity to distinguish fluxes along different spatial directions. The flux density along direction $𝒏 (𝒙)$ is given by $𝒇 (𝒙, τ) \cdot 𝒏 (x)$ . Consider $𝒏 (𝒙)$ as the inward pointing unit vector normal to the surface $\partial B$ that bounds the bank. The total flux is again obtained by integrating flux densities

F (τ) = \int_{\partial B} 𝒇 (𝒙, τ) \cdot 𝒏 (𝒙) d 𝒙 .

(7)

Gathering the above leads to re-expressing (4) or (5) gives

E (B, t + Δ t) = E (B, t) + \int_{t}^{t + Δ t} \int_{\partial B} 𝒇 (𝒙, τ) \cdot 𝒏 (𝒙) d 𝒙 d τ

(8)

Using (6) leads to the statement,

\int_{B} e (𝒙, t + Δ t) d 𝒙 = \int_{B} e (𝒙, t) d 𝒙 + \int_{t}^{t + Δ t} \int_{\partial B} 𝒇 (𝒙, τ) \cdot 𝒏 (𝒙) d 𝒙 d τ .

There are special cases in which additional events affecting the balance of $E$ can occur. When $B$ is a reserve bank money might be (legally) printed and destroyed in the building. Again by analogy with fluid dynamics, such events are said to be sources of $E$ within $B$ , much like a underground spring is a source of surface water. Let $Σ (t)$ be the total sources at time $t$ . As before, $Σ (t)$ might actually be obtained by summing over several sources placed in a number of positions, for instance the separate printing presses and furnaces that exist in $B$ . It is useful to introduce a spatial density of sources $σ (x, t)$ . The conservation statement now becomes

\int_{B} e (𝒙, t + Δ t) d 𝒙 - \int_{B} e (𝒙, t) d 𝒙 = \int_{t}^{t + Δ t} \int_{\partial B} 𝒇 (𝒙, τ) \cdot 𝒏 (𝒙) d 𝒙 d τ + \int_{t}^{t + Δ t} \int_{\partial B} σ (𝒙, τ) d 𝒙 d τ .

(9)

The above encompasses all physical conservation laws, and is quite straightforward in interpretation:

change in Euros in $B$ = net Euros coming in or going out of $B$ + net Euros produced or destroyed in $B$ .

It should be emphasized that the above statement has true physical meaning and is referred to as an integral formulation of a conservation law. The key term is “integral” and refers to the integration over some spatial domain.

Local formulations.

Equation (9) is useful and often applied directly in the analysis of physical systems. From an operational point of view it does have some inconveniences though. These have mainly to do with the integration domains

B

, sometimes difficult to describe and to perform integrations over. Avoid this by considering

𝒇 (𝒙, t)

defined everywhere, not only on

\partial B

(the doors and windows of

B

). These internal fluxes can be shown to have a proper physical interpretation. Assuming that

𝒇 (𝒙, t)

is smooth allows use of the Gauss theorem to transform the surface integral over

\partial B

into a volume integral over

B

\int_{\partial B} 𝒇 (𝒙, τ) \cdot 𝒏 (𝒙) d 𝒙 = - \int_{B} \nabla \cdot 𝒇 (𝒙, τ) d 𝒙

(10)

The minus assign arises from the convention of an inward pointing normal. Applying (10) to (9) leads to

\int_{B} [e (x, t + Δ t) - e (x, t) + \int_{t}^{t + Δ t} \nabla \cdot 𝒇 (𝒙, τ) d τ] d x = \int_{t}^{t + Δ t} \int_{\partial B} σ (𝒙, τ) d 𝒙 d τ .

(11)

There was nothing special about the shape of the building $B$ or the length of the time interval $Δ t$ we used in deriving (11), hence the equality should hold for infinitesimal domains

\frac{\partial e}{\partial t} + \nabla \cdot 𝒇 = σ,

(12)

where, as is customary, the dependence of $e, 𝒇, σ$ on space and time is understood but not written out explicitly. Equation (12) is known as the local or differential form of the conservation law for $E$ .

3.Special forms of conservation laws

Second law of dynamics.

The full general form (12) often arises in applications, but simplifications can arise from specific system properties. As a simple example, the dynamics of a point mass

m

which has no internal structure is described by the conservation of momentum statement

\frac{d}{d t} (m 𝒗) = \sum 𝑭 .

(13)

The correspondence with (12) is given by $e ⟷ (m v)$ , $σ ⟷ \sum 𝑭$ , hence the statement: “external forces are sources of momentum”. Instead of a PDE, the lack of internal structure has led to an ODE.

Advection equation.

Other special forms of (12) are not quite so trivial. Often

𝒇, σ

depend on

e

, with the specific form of this dependence is given by physical analysis. Accounting for all physical effects is so difficult that simple approximations are often used. For instance if

𝒇 (e)

is sufficiently smooth Taylor series expansion gives

𝒇 (e) = 𝒇_{0} + 𝒇^{'} (e_{0}) (e - e_{0}) + \dots .

(14)

Choosing the origin such that $𝒇_{0} = 0$ and $e_{0} = 0$ , the simplest truncation is

𝒇 (e) ≅ 𝒇^{'} (0) e = 𝒖 e,

(15)

and the $σ = 0$ form of (12) is

\frac{\partial e}{\partial t} + \nabla \cdot (𝒖 e) = 0 .

(16)

In this approximation $𝒖$ is a constant giving

\frac{\partial e}{\partial t} + 𝒖 \cdot \nabla e = 0

(17)

known as the constant velocity advection equation. Its one-dimensional form is the basis of much development in numerical methods for PDE's

\frac{\partial e}{\partial t} + u \frac{\partial e}{\partial x} = 0

(18)

Diffusion equation.

Another widely encountered dependence of

𝒇

e

is of the form

𝒇 (e) = - α \nabla e

(19)

and this leads to

\frac{\partial e}{\partial t} - \nabla \cdot (α \nabla e) = σ (e) .

(20)

If there are no sources and $α$ is a constant we have

\frac{\partial e}{\partial t} = α \nabla^{2} e

(21)

the heat or diffusion equation.

Combined effects.

Both above flux types can appear in which case the associated conservation law is

\frac{\partial e}{\partial t} + \nabla \cdot (𝒖 e) = α \nabla^{2} e,

(22)

known as the advection-diffusion equation, and is a linear PDE. If sources $σ$ exist the above becomes

\frac{\partial e}{\partial t} + \nabla \cdot (𝒖 e) = α \nabla^{2} e + σ,

(23)

or for constant advection velocity $𝒖$

\frac{\partial e}{\partial t} + 𝒖 \cdot \nabla e = α \nabla^{2} e + σ .

(24)

It is often the case that the flux depends on the conserved quantity itself, $𝒇 (e) = 𝒖 (e) e$ , in which case (22) becomes a non-linear PDE.

Steady-state transport.

Various effects can balance leading to no observable time dependence,

\partial e / \partial t = 0

. If there is no overall diffusive flux

𝒇 (e) = - α \nabla e

within an infinitesimal volume, then

\nabla \cdot 𝒇 = - α \nabla^{2} e = 0

leads to the Laplace equation

\nabla^{2} e = 0 .

If the infinitesimal volume contains sources the Poisson equation

\nabla^{2} e = σ,

is obtained.

Separation of variables.

Often, the time dependence can be isolated from the spatial dependence,

e (𝒙, t) = X (𝒙) T (t)

, in which case the diffusion equation for constant

α

leads to

\frac{\dot{T}}{T} = α \frac{\nabla^{2} X}{X} = - λ,

with $λ$ a positive constant to avoid unphysical exponential growth. The spatial part of the solution statisfies the Helmholtz equation

\nabla^{2} X = - κ^{2} X,

with $κ^{2} = λ / α$ . The above is interpreted as an eigenproblem for the Laplacian operator $Δ = \nabla^{2}$ .

The above special forms of differential conservation laws play an important role in scientific computation. Numerical techniques have been developed to capture the underlying physical behavior expressed in say the diffusion equation or the Helmholtz equation. These equations were first studied within physics, but they reflect universal behavior. Consider the Black-Scholes financial model for the price of an option $V (S, t)$ on an asset $S (t)$

\frac{\partial V}{\partial t} + r S \frac{\partial V}{\partial S} = r V - \frac{1}{2} σ^{2} S^{2} \frac{\partial^{2} V}{\partial S^{2}},

with $σ$ the standard deviation of stock market returns and $r$ the annualized risk-free interest rate. The terminology might be totally different, but the same patterns emerge and the Black Scholes model can be interpreted as an advection diffusion equation with non-constant advection velocity $r S$ , negative diffusion coefficient $- σ^{2} S^{2} / 2$ and source term $r V$ . The similarity to the physics advection and diffusion equations arises from the same type of modeling assumptions relating fluxes to state variables.