1  Conservation laws and constitutive laws

In the first chapter of (Fowler 2011), there is a concise introduction to the different categories of techniques and approaches that you might use when doing mathematical modelling in the real world. Some of these ideas will be introduced to you in this course.

Here, we provide a brief intro to the highlights, involving the use of conservation laws (and PDEs) and also the concept of non-dimensionalisation (which you would have encountered previously), studied in Chapter 2.

Vectors and PDEs

Vectors and PDEs is not a prerequisite for this course, but naturally in studying anything related to the physical real world, we must discuss partial differential equations. The hope is that the necessary theory for PDEs will be presented to you as this course evolves, so that it can appreciated by both newcomers and experienced readers.

Conservation laws can be expressed as mathematical equations that represent the idea that some quantity is conserved. In processes governing the planet, these might correspond to conservation of heat, of water, of air, of momentum, etc.

In Chapter 3, we will develop the simplest possible model governing the temperature on the surface of the Earth. It is a conservation equation for energy and is zero-dimensional (does not involve time and does not involve spatial variation).

1.1 Derivation of the 1D heat equation

In order to demonstrate some of the basic principles of this course, let us demonstrate the derivation of the heat equation. We are interested in modelling the heat in a volume, \(V\), which, for the sake of concreteness is given by a long cylinder with its axis along \(x \in [0, L]\). We assume that the side walls of the cylinder are insulated and the temperature only varies along the \(x\) direction.

At any point along this rod, the internal heat is given by \(\rho c T(x, t)\), where \(\rho\) is the density of the material (kg/m3), \(c\) is the specific heat capacity (J/(kg K)), and \(T\) is the temperature (K). Therefore, the heat energy along any segment in the rod is calculated from

Internal heat energy

\[ \text{heat energy in $[a, b]$} = \int_a^b \rho c T \, \mathrm{d}x. \]

If the heat changes, then the rate of change of heat energy is given the time derivative of the above quantity. By conservation of energy, any change of the internal energy must be equal to the inflow or outflow of heat at the ends, \(x = a\) or \(x = b\). We therefore write \(q\) for the flux (or flow) of heat.

We need a constitutive law that dictates how energy is exchanged at the boundaries. Based on intuition, it is sensible to assume that the flow of heat proceeds from hot to cold. For example, hot air rises towards cool air; or heat from a hot mug of tea flows and diffuses outwards into a cold room. Therefore, we write this as

Fourier’s law

Fourier’s law in 1D specifies that the heat flux is given by \[ q(x, t) = -k \frac{{\partial}T}{{\partial}x}. \]

This is known as Fourier’s law. The quantity \(k\) is the thermal conductivity, and its units are W/(m K). Because a Watt is a Joule/s, you can also see that the units of \(k\) are J/(m K s). The quantity \(q\) is the flux, and you can verify that it is given in units of J/(m2 s).

Therefore by energy conservation, we have \[ \frac{\mathrm{d}}{\mathrm{d}t} \int_a^b \rho c T \, \mathrm{d}x = {\color{blue}q(x = a, t) - q(x = b, t)}, \] i.e. the change in internal heat is equal to the flow through the ends. Note that we have added-in the heat flux due to \(x = a\) assuming positive \(q\) refers to heat moving from left-to-right. Conversely, we subtract it away from \(x = b\). We can alternatively write this as \[ \frac{\mathrm{d}}{\mathrm{d}t} \int_a^b \rho c T \, \mathrm{d}x = {\color{blue} -\int_a^b \frac{\partial q}{\partial x} \, \mathrm{d}x}. \]

Substitution Fourier’s law, we can then write \[ \int_a^b \rho c \frac{{\partial}T}{{\partial}t} \, \mathrm{d}x = \int_a^b k \frac{{\partial}^2 T}{{\partial}x^2} \, \mathrm{d}x. \] Because the above integral identity needs to be true for all possible values of \(a\) and \(b\), then it must be true everywhere (this is sometimes referred to as the du Bois-Reymond lemma or the bump lemma. Therefore we are left with the classic heat equation.

Heat equation

\[ \rho c\frac{{\partial}T}{{\partial}t} = k \frac{{\partial}^2 T}{{\partial}x^2}. \tag{1.1}\]

In order to produce a sensible physical solution, partial differential equations are typically supplemented by initial conditions and boundary conditions. The initial condition prescribes the state of the function at some initial time, typically \(t = 0\). Boundary conditions prescribe how the function behaves on the boundary of its domain, which in this case is \(x = 0\) and \(x = L\). An example might be

Initial conditions (IC) and boundary conditions (BC)

\[ \begin{gathered} T(x, 0) = T_0 \\ T(0, 0) = T_a \\ T(L, 0) = T_b \end{gathered} \]

which expresses, respectively, that the temperature starts from a constant temperature, \(T_0\), and where the ends of the rod are kept at temperature \(T_a\) and \(T_b\).

1.1.1 Steady states and long-time behaviours

\(\nextSection\)
2024 note

This was a new addition in 2024.

When we refer to a steady-state solution, we are typically referring to a time-independent solution.

Definition 1.1 (Steady-state solutions) Given an evolving system described by a functio, say \(f(x, t)\), defined on some spatial domain and with \(t > 0\), the steady-state solution refers to time-independent solutions with \[ \frac{{\partial} f}{{\partial}t} = 0. \]

One can envisage that as the system evolves with \(t \to \infty\), it reaches a state that is independent of time. However, not all systems will approach a steady state. Moreover, not all steady states are stable or attractive (and might never be reacheable in a real-life experiment).

For the case of heat flow, such a steady-state solution would be \(T(x, t) = T(x)\). In this case, \[ \frac{{\partial}T}{{\partial}t} = 0 \Longrightarrow k \frac{{\partial}^2 T}{{\partial}x^2} = 0. \]

Therefore, for the heatflow along a segment of length \(L\) with left boundary held at \(T_a\) and right boundary held at \(T_b\) we have \[ T_\text{steady}(x) = \left(\frac{T_b - T_a}{L}\right) x + T_a. \]

1.2 Deriving the 1D transport (continuity) equation

\(\nextSection\)
2024 note

This was a new addition in 2024.

Consider the mass transport of some substance with density \(\rho(x, t)\), immersed in a fluid, along a one-dimensional line in \(x\). If desired, you may consider the substance as existing in three-dimensional space, and propagating along the \(x\)-direction, with its behaviour independent of \(y\) and \(z\). Initially, when considered at time \(t\), the mass of the substance between two points, \(a\), and \(b\), is:

\[ m_\text{blob}(t) = \int_{a(t)}^{b(t)} \rho (x, t) \, \mathrm{d}x. \]

As time increases the particles of the substance will move due to the fluid moving; at the same time, the fluid volume which is initially contained in \(x\in [a(0), b(0)]\) will also move. We want to find how the mass of the blob changes in time, and hence consider the quantity \[ \frac{\mathrm{d}m_\text{blob}}{\mathrm{d}t} = \int_a^b \frac{{\partial}\rho}{{\partial}t} \, \mathrm{d}x + \rho(b, t) \frac{\mathrm{d}b}{\mathrm{d}t} - \rho(a, t) \frac{\mathrm{d}a}{\mathrm{d}t}. \]

The above considers the intrinsic rate of change of the function within the integrand, but then adds the extra mass due to the right boundary shifting rightwards (\(b\)) and subtracts the mass due to the left boundary shifting rightwards (\(a\)). It is known as the Leibniz integral rule. We can thus write this within the integral as

\[ \frac{\mathrm{d}m_\text{blob}}{\mathrm{d}t} = \int_a^b \left[ \frac{{\partial}\rho}{{\partial}t} + \frac{{\partial}}{{\partial}x} \left(\rho(x, t) \frac{\mathrm{d}x}{\mathrm{d}t}\right) \right] \, \mathrm{d}x. \]

However, the quantity \[ \frac{\mathrm{d}x}{\mathrm{d}t} \equiv u(x, t) \]

represents the velocity of the fluid (which for the moment we assume to be a known and provided quantity). Therefore we can write the mass change as

\[ \frac{\mathrm{d}m_\text{blob}}{\mathrm{d}t} = \int_a^b \left[ \frac{{\partial}\rho}{{\partial}t} + \frac{{\partial}}{{\partial}x} \left(\rho u\right) \right] \, \mathrm{d}x. \]

This result, which explains how to pass a derivative through an integral express mass of a substance within a flow is known as the Reynolds Transport Theorem. We have just derived it in 1D.

Reynolds Transport Theorem

Let \(\rho = \rho(x, t)\) be some quantity (such as density) that is advected along a one-dimensional line in \(x\) due to a fluid with velocity \(u(x, t)\). Then \[ \frac{\mathrm{d}}{\mathrm{d}t} \int_{a(t)}^{b(t)} \rho \, \mathrm{d}x = \int_{a(t)}^{b(t)} \left[ \frac{{\partial}\rho}{{\partial}t} + \frac{{\partial}}{{\partial}x}(\rho u)\right] \mathrm{d}x. \]

We may now consider the substance being transported along the \(x\)-direction. If there is no interior creation or destruction of the source, then by conservation of mass, it must be the case that \[ \frac{\mathrm{d}m_\text{blob}}{\mathrm{d}t} = 0. \]

Thus, again since the above integral identity applies to all possible values of \(a\) and \(b\), it must be the case that the integrand is zero. Thus we conclude with the so-called tranport equation.

One-dimensional transport equation

The transport of a substance described by \(\rho(x, t)\) advected along a one-dimensional line in \(x\) due to a fluid with velocity \(u(x, t)\) is given by \[ \frac{{\partial}\rho}{{\partial}t} + \frac{{\partial}}{{\partial}x}(\rho u) = 0. \]

Again, we must consider the above problem in combination with potential initial conditions and boundary conditions. For instance, we might specify that the substance begins from some initial state, say \[ \rho(x, 0) = \rho_0(x). \]

The boundary conditions are more subtle. It is not always obvious what the boundary conditions should be on a problem.

Here is an example of a solution of such a problem.

Linear advection equation with constant speed

Solve the problem given by \[ \frac{{\partial}\rho}{{\partial}t} + c\frac{{\partial}\rho}{{\partial}x} = 0, \] where the velocity \(u = c\) is constant. You may assume the initial condition is given by \[ \rho(x, 0) = \rho_0(x) = e^{-x^2}. \]

Solution

The general solution is given by \[ \rho(x, t) = F(x - ct), \] where \(F\) is an arbitrary (differentiable) function. You can verify this via the chain rule, noting that \[ \frac{\partial F}{\partial t} = -c F'(x - ct) \quad \text{and} \quad \frac{\partial F}{\partial x} = F'(x - ct). \] Applying the initial condition at \(t = 0\) gives \(F(x) = \rho_0(x)\) and thus \[ \rho(x, t) = e^{-(x - ct)^2}. \]

You may plot this to discover that it is, as perhaps expected, a Gaussian profile moving to the right at constant speed (for \(c > 0\)).