Lithospheric strength profile

Lithospheric strength profile#

For rifts to form, the extensional stresses must exceed the lithospheric strength profile, expressed by the brittle and ductile deformation types.

Following the Mohr-Coulomb failure criterion, the shear stress \(\tau\) required to cause failure is given by:

(1)#\[\begin{equation} \tau = P + C \tan(\theta) \end{equation}\]

where \(C\) is the cohesion, \(P\) is the pressure (i.e., \(\rho g h\)) , and \(\theta\) is the internal angle of friction.

At deeper depths, the above profile does not hold true as the ductile deformation becomes dominant. The ductile deformation is governed by the power law creep, which can be expressed as:

(2)#\[\begin{equation} \sigma = A^{-1/n}\dot{\varepsilon}^{1/n}\exp(\frac{E}{nRT}) \end{equation}\]

where \(\sigma\) is the stress, \(A\) is a prefactor, \(\dot{\varepsilon}\) is the strain rate, \(n\) is the stress exponent, \(E\) is the activation energy, \(R\) is the gas constant, and \(T\) is the temperature.

../../../_images/strength_envelope.png

Conceptual strength envelope assuming a simplified temperature and pressure distribution. Image: modified from Brune et al. (2023)

Plotting the strength profile#

We can plot the above set of equations to visualize the strength profile of the lithosphere that needs to be overcome for rifting to occur.

import matplotlib.pyplot as plt
from scipy import special
import numpy as np
import matplotlib.cm as cm
import matplotlib.colors as colors

# define the constants appropriate for the problem
g = 9.8              # Gravitational acceleration
R = 8.31451          # Gas constant
strainrate = 1e-15   # s^-1

# ------------------- Temperatures------------------------
# we define continental geotherm and oceanic geotherm as functions of depth
def continental_geotherm(depth):
    'Return temperatures based on a constant linear geothermal gradient'
    T_surface  = 273.15; # K
    T_gradient = 0.02; # K/m
    return T_surface + T_gradient * depth

def oceanic_geotherm(depth):
    'Return temperature based on oceanic cooling model'
    T_surface  = 273.15; # K
    T_mantle   = 1673.15; # K
    diffusivity= 1e-6; # m^2/s
    plate_age  = 80e6 * 365.25 * 24 * 3600; # plate of 80 Ma

    return T_surface + (T_mantle - T_surface) * special.erf(depth / (np.sqrt(diffusivity * plate_age)))

# ---------------- Rheology------------------------
# Upper crust rheology assuming Wet Quartz (Rutter & Brodie, 2004)
A_wq = 8.57e-28 # prefactor Pa^-n s^-1
n_wq = 4.       # stress exponent
E_wq = 223.e3   # activation energy J/mol

# Lower crust rheology assuming Wet Anorthite (Rybacki et al., 2006)
A_wa = 7.13e-18 # prefactor Pa^-n s^-1
n_wa = 3        # stress exponent
E_wa = 345.e3   # activation energy J/mol

# Upper mantle rheology assuming dry Olivine (Hirth & Kohlstedt, 2003)
A_do = 6.52e-16 # prefactor MPa^-n s^-1
n_do = 3.5      # stress exponent
E_do = 530.e3   # activation energy J/mol

# Upper mantle rheology assuming wet Olivine (Hirth & Kohlstedt, 2003)
A_wo = 5.33e-19
n_wo = 3.5
E_wo = 480.e3

# ---------------Brittle strength envelope---------------------
density  = 3250 # kg/m^3
depths   = np.linspace(0, 100e3, 100) # depth in m
pressure = g*density*depths

phi=20./180.*3.1416
cohesion = 2e6
yield_stress=pressure*0.6+cohesion # theta = 30 degrees

def ductile_stress(A, n, E, T):
    'Calculate the ductile stress for given strain rate and rheological parameters'
    return (strainrate / A * np.exp(E / (R * T))) ** (1/n)

# plot lithosphere strength profile for continental and oceanic lithosphere
strength_continental = np.zeros_like(depths)
strength_oceanic     = np.zeros_like(depths)

# Let's assume the following structure for continental (left) and oceanic
# lithosphere (right)

# ------------------------------------|---------------------------------
# upper crust (20km) wet quartz       |  upper crust (10km) wet quartz
# ------------------------------------| --------------------------------
# lower crust (20-40km) wet anorthite |  mantle (40-150km) wet olivine
# ------------------------------------|---------------------------------
# mantle (40-150km)  dry olivine      |
# ------------------------------------|

# Calculate stress profile for continental lithosphere

for d in range(len(depths)):
    if depths[d] <= 20e3:
        strength_continental[d] = ductile_stress(A_wq, n_wq, E_wq, continental_geotherm(depths[d]))
    elif 20e3 < depths[d] <= 40e3:
        strength_continental[d] = ductile_stress(A_wa, n_wa, E_wa, continental_geotherm(depths[d]))
    else:
        strength_continental[d] = ductile_stress(A_do, n_do, E_do, continental_geotherm(depths[d]))


for d in range(len(depths)):
    if depths[d] <= 10e3:
        strength_oceanic[d] = ductile_stress(A_wq, n_wq, E_wq, oceanic_geotherm(depths[d]))
    else:
        strength_oceanic[d] = ductile_stress(A_wo, n_wo, E_wo, oceanic_geotherm(depths[d]))

# compute the effective strength as the minimum of brittle and ductile strength curves
effective_strength_continental = np.minimum(strength_continental, yield_stress)
effective_strength_oceanic     = np.minimum(strength_oceanic, yield_stress)

# plot the strength profiles
fig, ax = plt.subplots(1, 2, figsize=(8, 5), sharey=True)
ax[0].plot(effective_strength_continental / 1e6, depths / 1e3, label='Continental lithosphere', color='green')
ax[1].plot(effective_strength_oceanic / 1e6, depths / 1e3, label='Oceanic lithosphere', color='orange')
ax[0].invert_yaxis()
ax[0].set_xlabel('Strength (MPa)')
ax[1].set_xlabel('Strength (MPa)')
ax[0].set_ylabel('Depth (km)')
ax[0].set_title('Continental Lithosphere')
ax[1].set_title('Oceanic Lithosphere')
../../../_images/strength_profiles.png

In the above plot, any applied force exceeding the strength of the lithosphere would create cracks (brittle deformation) or shear zones (ductile deformation).

Things to try#

  • How does the strength profile change with increasing strain rate?

  • Investigate the variations in the strength profile with different temperature profiles. Would it require more or less total stress to initiate rifting in tectonic regions with colder geotherms.

References#

  • Hirth, G., & Kohlstedt, D. (2003). Rheology of the upper mantle and the mantle wedge: A view from the experimentalists. Geophysical monograph series, 138, 83-105.

  • Rutter, E. H., & Brodie, K. H. (2004). Experimental intracrystalline plastic flow in hot-pressed synthetic quartzite prepared from Brazilian quartz crystals. Journal of Structural Geology, 26(2), 259-270.

  • Rybacki, E., Gottschalk, M., Wirth, R., & Dresen, G. (2006). Influence of water fugacity and activation volume on the flow properties of fine‐grained anorthite aggregates. Journal of Geophysical Research: Solid Earth, 111(B3).

  • Brune, S., Kolawole, F., Olive, J. A., Stamps, D. S., Buck, W. R., Buiter, S. J., … & Shillington, D. J. (2023). Geodynamics of continental rift initiation and evolution. Nature Reviews Earth & Environment, 4(4), 235-253.

 


→ NEXT: ASPECT.prm file description    icon