Solving equations

Goals of this lesson

• Provide a quick overview of what we need to know in order to solve equations

• Establish the necessary background for why we need ot use numerical solutions to equations

Finding a solution

We have now seen the general form of several different equations we will be using to model geodynamic processes, such as the heat transfer equation below

\[ \rho c_{\mathrm{p}} \left( \frac{\partial T}{\partial t} + V \cdot \nabla T \right) = \nabla \cdot k \nabla T + A \]

but how do we go from this form of the equation to a solution in the form

\[ T(x, y, z, t) = \mathrm{...?} \]

A simpler example

Let’s start simpler.

Steady-state heat conduction in one dimension (1D)

\[ \rho c_{\mathrm{p}} \left( \frac{\partial T}{\partial t} + v_{z} \frac{\partial T}{\partial z} \right) = \frac{\partial}{\partial z} \left(k \frac{\partial T}{\partial z} \right) + A \]

where

• Time dependence is represented by \(\frac{\partial T}{\partial t}\)

• Advection is represented by \(v_{z} \frac{\partial T}{\partial z}\)

• Conduction is represented by \(\frac{\partial}{\partial z} \left(k \frac{\partial T}{\partial z} \right)\)

• Heat production is represented by \(A\)

In the case of steady-state heat conduction without advection or heat production, our equation simplifies considerably:

\[ \rho c_{\mathrm{p}} \left( \cancel{\frac{\partial T}{\partial t}} + \cancel{v_{z} \frac{\partial T}{\partial z}} \right) = \frac{\partial}{\partial z} \left(k \frac{\partial T}{\partial z} \right) + \cancel{A} \]

resulting in

\[ \frac{\partial}{\partial z} \left(k \frac{\partial T}{\partial z} \right) = 0 \]

In fact, we can even factor out the thermal conductivity \(k\) to yield an even simpler relationship:

\[\begin{split} \begin{aligned} \frac{\partial}{\partial z} \left(k \frac{\partial T}{\partial z} \right) &= 0 \\ \frac{\partial}{\partial z} \left(\frac{\partial T}{\partial z} \right) &= 0 \\ \frac{\partial^{2} T}{\partial z^{2}} &= 0 \end{aligned} \end{split}\]

Here we are thus left with simply the second derivative of temperature \(T\) with respect to depth \(z\). We can solve this! We just need to integrate twice with respect to \(z\).

A simpler example?

We can integrate with respect to \(z\) once to find

\[\begin{split} \begin{aligned} \int \frac{d^{2} T}{d z^{2}} dz &= \int 0 dz \\ \frac{d T}{d z} &= c_{1} \\ \frac{d T}{d z} + c_{1} &= 0 \end{aligned} \end{split}\]

noting both integration constants would be combined as \(c_{1}\).

A second integration yields

\[\begin{split} \begin{aligned} \int \left(\frac{d T}{d z} + c_{1} \right) dz &= \int 0 dz \\ \int \frac{d T}{d z} dz + \int c_{1} dz &= \int 0 dz \\ T(z) + c_{1}z &= c_{2} \end{aligned} \end{split}\]

which can be rearranged to

\[ T(z) + c_{1}z + c_{2} = 0 \]

\[ T(z) = -c_{1}z - c_{2} \]

OK, now what?

\[ T(z) = -c_{1}z - c_{2} \]

OK, so we now have a solution, but we cannot calculate temperatures yet because we don’t know the values of \(c_{1}\) and \(c_{2}\), the constants of integration.

• What do we need to be able to find \(c_{1}\) and \(c_{2}\)?

Boundary conditions

By integrating twice we end up with two unknown constants.

In order to find the values of the constants, we need to use boundary conditions, known (or assumed) values for certain variables in our equation.

There are two natural choices for our heat transfer problem:

• Known temperatures at certain depths

• Known temperature gradients at certain depths

\[ T(z) = -c_{1}z - c_{2} \]

Let’s make some assumptions:

• We know the surface temperature at \(z = 0\) is \(T = 0\)

• We know the temperature at some depth \(z = L\) is \(T = 1000 ^{\circ}\mathrm{C}\)

If we substitute in assumption 1 we find \(c_{2} = 0\). If we substitute in assumption 2 we find \(c_{1} = -1000 / L\).

Now we can finally see our exciting result with \(c_{1} = -1000 / L\) and \(c_{2} = 0\)…

\[\begin{split} \begin{aligned} T(z) &= -c_{1}z - c_{2} \\ T(z) &= \frac{1000}{L} z \end{aligned} \end{split}\]

…a straight line.

Other considerations

So that was an almost embarrassingly simple example

• Heat conducted between two known temperatures at steady state will simply produce a linear temperature increase across the thickness of the material

What happens, though, if we consider a heat transfer problem that is not at a steady state?

• What other information might we need to know for that case?

Initial conditions

In order to consider the time evolution of the temperatures in the Earth, we would also need to know the starting distribution of temperature

• This is known as an initial condition

Clearly, the temperatures we expect after some time can strongly depend on the initial conditions we assume

Why don’t we need an initial condition for the steady-state case?

What about other scenarios?

What if we want to solve slightly more complicated problems?

• 2D heat transfer with surface topography

• Spatially variable material properties

• Temperature-dependent material properties

In most of these cases we cannot simply directly integrate our heat transfer equation because we cannot find the constants of integration.

• Consider the example above with topography

• We likely do not have a function that can give us the required known values of temperature and elevation at the surface

The need for numerical integration

In cases where the geometry of the problem or material property distributions/behaviors are more complex, we need to use numerical methods to integrate our equations of interest

Much of days 2 and 3 in this course will focus on the finite difference method for solving equations, how it can be used, and its limitations.

The finite difference method is one of the most popular approaches used for solving equations in geodynamic models

• We will start learning the details of the finite difference method tomorrow (Tuesday) morning

Summary

To solve equations used in geodynamic modelling (or anything really) we need to know

• Boundary conditions

• Initial conditions (in some cases)

We can solve these equations directly in some cases, but most of the time we will need to use numerical methods such as the finite difference method to find solutions

References

• Stüwe, K. (2007). Geodynamics of the lithosphere: an introduction. Springer Science & Business Media.