Skip to main content
Geosciences LibreTexts

2.2: Geological Applications of the Diffusion Equation

  • Page ID
    3527
  • General Diffusion Equation

    Diffusivity: Note that \(\kappa\) has units of \(\frac{m^2}{s}\), which does not depend on temperature or heat flow. Diffusivity, in general, describes the spreading out of some quantity (e.g. temperature) over space in some unit of time. Other diffusion processes will have other diffusivities, but they all will have the same unit. For each of the equations below, the solution is found by defining an initial condition (the initial temperature, topography, or concentration at every position) and two boundary conditions (the value of the temperature, topography, or concentration at two positions...they can be fixed or functions of time). We will go over these conditions in more detail in the next few sections.

    1-D vs 2-D/3-D: Heat spreads out in 2-D or 3-D because of the addition of the \(\frac{\partial^2 T}{\partial y^2}+\frac{\partial^2 T}{\partial z^2}\) term to the diffusion equation \(\frac{dT}{dt}=\kappa\frac{d^2T}{dz^2}\). The higher the degree of curvature on the graph indicates the direction of maximum cooling/heating.

    FIGURE cooling graph.

    Additionally, if \(\kappa\) is a constant, it will work the same in 1-D as it does in 2-D or 3-D, heat will just spread out in multiple directions.

    FIGURE T contours hot/cold.

    Cooling of the Oceanic Lithosphere

    \[\frac{dT}{dt}=\kappa\frac{d^2T}{dz^2}\]

    where \(\kappa\) is the thermal diffusivity. The solution for this equation is:

    \[T(x, t)=(T_s-T_m)erfc(\frac{y}{2\sqrt{\kappa t}})+T_m\]

    where \(T_s\) is the surface temperature, \(T_m\) is the mantle temperature, and y is the depth into the lithosphere. We will use this solution all the time, as it indicates the thermal structure of the oceanic lithosphere. It does not take into account reheating due to volcanism or mechanical feedback, but it is used to define the initial condition in models. We will use it to calculate the ridge push force. 

    As you can see above, in order to solve the partial differential diffusion equation, we took several large mathematical steps. Let's explore how we got there in more detail. Solving partial differential equations (PDE) might seem challenging at first, especially if you have never taken a class on PDE before, but by breaking it into a few simple steps, it becomes less daunting. Let's starting with the basic diffusion equation:

    \[\frac{dT}{dt}=\kappa\frac{d^2T}{dz^2}\]

    Now we will integrate to find \(T=f(x, t)\). \(\frac{dT}{dt}\) on takes the derivative in terms of t and holds x constant. For each derivative there will be a constant of integration. 

    \[T(x)=\int_{0}^{x_o}x dx=\frac{1}{2}x^2+C_o\]

    In order to find \(C_o\), we need a boundary condition. How many boundary conditions do we need? One for every derivative. If we have \(\frac{d}{dt}\), we need one initial condition. \(\int dt\rightarrow C_o \rightarrow time\)⇒T at all x locations at t=0. If we have \(\frac{d^2}{dx^2}\), we need two boundary conditions. \(\int(\int dx)dx \rightarrow C_1x+C_2\rightarrow\) space and distance ⇒T at all times at x=0.

    For almost every diffusion equation you come upon, the solution to the equation will include dependence on the error function.

    \[erf(\eta)=\frac{2}{\sqrt{\pi}} \int_{0}^{\eta} e^{-\eta^{\prime 2}} d\eta\prime\]

    The parameter \(\eta\) depends on the distance, time, and the diffusion coefficient, \(\eta=\frac{z}{2\sqrt{\kappa t}}\). Note that if you set your \(\eta=1\), you are back to the simple equation that we used to determine the length-scale and time-scale for diffusion. The value of the error function for any particular value of \(\eta\) can be calculated using a built in function or looked up in a table. There is also a related function called the complementary error function. 

    \[erfc(\eta)=1-erf(\eta)\]

    FIGURE error function graph

    Here is a plot of the error function and the complimentary error function (in the plot x is \(\eta\)). The final step to solve the PDE is to use the initial conditions and boundary conditions to determine the constants of integration in the general solution. 

    FIGURE Oceanic lithosphere and x,z figures. Now that we have a general idea of how to solve a PDE, let's return to the oceanic lithosphere example. Lithosphere forms by cooling of the mantle as the tectonic plate moves away from a spreading ridge. At the ridge, the mantle temperature of \(T_m\) almost reaches the surface (z=0). As the material formed at the ridge moves away from the spreading center at a rate equal to the speed of the tectonic plate (\(u\)), the mantle cools from above to become the lithosphere. The temperature at the surface of the lithosphere, \(T_s\), is essentially constant because the surface is at the bottom of the ocean (\(T=0^oC)\). The age of the plate at a distance \(x\) from the spreading center is determined by the speed of the plate: \(t=\frac{x}{u}\). A typical value for the thermal diffusion coefficient (thermal diffusivity) is \(\kappa_T=1.0x10^{-6}\frac{m^2}{s}\)

    at t=0  T(z=0)=\(T_m\)         and,       at z=0 T(t>0)=\(T_s\)

    t=\(0_+\) T(z=0)=\(T_s\)                            z=\(\infty\) T(t>0)=\(T_m\)

    \[\eta=\frac{z}{2\sqrt{\kappa t}}\]

    \[\frac{dT}{dt}=\kappa\frac{d^2T}{dz^2}\]

    For all times at z=0 \(T=T_s\)        and,     at t=0 T(z)=\(T_m\)

                          z=\(\infty\) \(T=T_m\)

    \[T(z, t)=(T_s-T_m) erfc(\frac{z}{2\sqrt{\kappa t}})+T_m\]

    at z=0, \(\eta\)=0, erfc=1 ⇒Ts-Tm+To⇒T=Ts

    at z=\(\infty\), \(\eta=\infty\), erfc=0⇒T=Tm

    at t=0, \(\eta=\infty\), erfc=0, ⇒T=Tm

    We can rewrite this using erf(\(eta\)).

    \[T=(T_s-T_m)erfc(\eta)+T_m\]

    Plug in \(erf(\eta)=(1-erf(\eta))\)

    \[=(T_s-T_m)(1-erf(\eta))+T_m)\]

    \[=(T_s-T_m)-(T_s-T_m)erf(\eta)+T_m\]

    \[=T_s-(T_s-T_m)erf(\eta)\]

    \[=T_s+(T_m-T_s)erf(\eta)\]

    \[T(z, t)=T_s+(T_m-T_s)erf(\frac{z}{2\sqrt{\kappa t}})\]

    FIGURE Plot T(z, t)

    Diffusion of Topography 

    For the case of erosion of a fault scarp with initial height, H, the diffusion equation can be written as:

    \[\frac{\partial H}{\partial t}=\kappa_H(\frac{\partial^2H}{\partial x^2})\]

    where \(\kappa_H\) if the diffusivity for fault scarp topography. The solution to this equation is:

    \[H(x, t)=erf(\frac{x}{2\sqrt{\kappa_H t}})\alpha+bx\]

    where \(\alpha\) is half the total fault offset and \(b\) is the background slope of topography.

    Diffusion of Chemical Species 

    For the case of diffusion of some chemical element denoted by C, the diffusion equation can be written as:

    \[\frac{\partial C}{\partial t}=\kappa_C(\frac{\partial^2C}{\partial x^2})\]

    where \(\kappa_C\) is the diffusivity of the chemical element and \(C\) is the concentration of the chemical element. The solution to this equation is:

    \[C(x, t)=(C_s-C_o)erf(\frac{x}{2\sqrt{\kappa_C t}})+C_o\]