2.3: Lithospheric Thickness Jupyter Notebook

Instructions

• To use this page, read the text between each of the python code windows, then press RUN to execute the code in the box.
• It is necessary to run each of the boxes in order.
• Don't click restart in the cells unless you come back to the top of the page and start over.
• If you modify the code, this modification will not remain after you leave or refresh this page.
• You must run simulations/cells that require user input to completion. Make sure you see the Message informing you a cell is complete before continuing to the next cell. Here's why you must do this:
  • These cells may not be re-run until the simulation is completed, otherwise no output will be generated.
  • All following cells will not run until the user input cell is completed.
• For Dropdown Menus: You do not need to rerun the cells with the menu to change your selection. However, you must rerun all of the following cells to implement the change.

Be patient, some times it can take 1-2 minutes for the juypter kernel to start.

An interactive example of how to calculate and plot the thickness of the lithosphere over time and distance from a spreading center.

Key Questions: Consider these as you work your way through this page.

• Describe how does the thickness of the lithosphere change as it moves away from the ridge in time and distance. Follow one isotherm and observe how it changes. Is the rate of thickening constant? Where/when is it faster or slower?
• Why does the thickness change with distance from the ridge?
• Why is the lithosphere thicker at a distance 1000 km for u = 2 cm/yr compared to u = 5 cm/yr?

First, import the necessary packages:

import matplotlib.pyplot as plt
import numpy as np
from scipy.special import erfcinv

print('The packages have been imported, move on the the next cell.')

Next, we define parameters and set up the necessary function. kappa is the thermal diffusion coefficient, age is the age of the lithosphere (it is created at the ridge), and Ts is the surface temperature. Tm is the temperature of the mantle, which is taken to be 1400C (ignoring the adiabatic depth dependence of temperature in the mantle). yr2my and sec2yr are conversion factors, and T is the range of isotherms that will be plotted. u_range is the range of spreading rates that we will consider. Then a diffusion profile function is created. This function does some unit conversions and solves for depth given the diffusion coefficient, age, spreading rate and the surface, mantle, and isotherm temperatures. The equation used in the function is the half-space cooling equation derived in 2.2: Geological Applications of the Diffusion Equation, but solved for depth.

kappa = 1e-6 # m^2 / s
age = np.arange(0.1,101,1) # my
Ts = 0 # C
Tm = 1400 # C
yr2my = 1e6   # yr/my
sec2yr = 3600*24*365   # s/yr 
T = np.arange(600,1300,200) # C
u_range = np.arange(1,6) #cm/yr

def diff_profile_z(Ts,Tm,kappa,age,T,u):
    u = u / 3.15e7 / 100 #Convert to m/s
    agesec = age*yr2my*sec2yr #Converts age to seconds
    zm = 2*np.sqrt(kappa*agesec)*erfcinv((T-Tm)/(Ts-Tm)) # Calculates Halfspace cooling equation, but solving for z in meters
    zkm = zm/1000 #Convert to km
    return zkm

print('The parameters and function have been defined. Move to the next cell.')

The following cell plots the thickness of the lithosphere vs time and distance away from the ridge. The for loop changes the spreading rate and updates the graph accordingly in each iteration. u_current is the current spreading rate in the for loop from the array of spreading rates defined above. The first set of code in the for loop calculates the distance traveled given u_current. The ax1.plot commands are dummy lines needed to provide a second x-axis for age to the plot. For each spreading rate in the for loop, the diffusion function is called to calculate the depth for each isotherm (this is done in the lines starting with line.set_ydata- line3.set_ydata). The horizontal distance traveled is also calculated here for each velocity, x_new.

handle = display(None, display_id=True)

y= np.zeros(len(age))

bbox = dict(boxstyle ="round", fc = '1') 
fig, ax = plt.subplots(figsize = (15,5))
ax1 = ax.twiny()
line, = ax.plot(age,y,color='purple')
line1, = ax.plot(age,y,color='blue')
line2, = ax.plot(age,y,color='green')
line3, = ax.plot(age,y,color='red')
line4, = ax1.plot(age,-np.ones(len(age)))
ax.set_ylim(130,0)
ax.set_xlim(0,5000)
ax.set_xlabel('Distance Traveled Away From Ridge (km)')
ax1.set_xlabel('Age (myr)')
ax.set_ylabel('Depth (km)')
ax.annotate('Ridge $ T_s = 0^\circ C $', (0,0),(500,10),arrowprops = {'arrowstyle':'->','color':'black'})
ax.set_title('Lithospheric Thickness versus Distance from the Ridge (Calculated to an age of 100 Myr)')
for u_current in u_range:
    print('Now Plotting:', u_current, 'cm/yr')
    agesec = age*yr2my*sec2yr #Convert age from myr to s
    x_new = agesec * (u_current/ 3.15e7 / 100) #m
    age_new = np.arange(0,5000)*1000 / (u_current/ 3.15e7 / 100) / yr2my / sec2yr #Myr
    x_new /= 1000 #Convert to km
    line4.set_ydata(-np.ones(len(age_new)))
    line4.set_xdata(age_new)
    ax1.set_xlim(0,np.max(age_new)+1)
    line.set_ydata(diff_profile_z(Ts,Tm,kappa,age,T[0],u_current))
    line1.set_ydata(diff_profile_z(Ts,Tm,kappa,age,T[1],u_current))
    line2.set_ydata(diff_profile_z(Ts,Tm,kappa,age,T[2],u_current))
    line3.set_ydata(diff_profile_z(Ts,Tm,kappa,age,T[3],u_current))
    line.set_xdata(x_new)
    line1.set_xdata(x_new)
    line2.set_xdata(x_new)
    line3.set_xdata(x_new)
    T1 = ax.annotate('T = '+str(T[0]) + '$ ^{\circ} C $',(np.max(x_new)-350,diff_profile_z(Ts,Tm,kappa,age,T[0],u_current)[-1]+5),color = 'purple')
    T2 = ax.annotate('T = '+str(T[1]) + '$ ^{\circ} C $',(np.max(x_new)-350,diff_profile_z(Ts,Tm,kappa,age,T[1],u_current)[-1]+5),color = 'blue')
    T3 = ax.annotate('T = '+str(T[2]) + '$ ^{\circ} C $',(np.max(x_new)-350,diff_profile_z(Ts,Tm,kappa,age,T[2],u_current)[-1]+5),color = 'green')
    T4 = ax.annotate('T = '+str(T[3]) + '$ ^{\circ} C $',(np.max(x_new)-350,diff_profile_z(Ts,Tm,kappa,age,T[3],u_current)[-1]+5),color = 'red')
    annotation = ax.annotate('u = ' + str(u_current) + ' cm/yr',(3500,20),bbox=bbox,fontsize=12)
    handle.update(fig)
    if u_current == u_range[0]:
     	print('Click enter or return to continue: ')
    wait = input()
    annotation.remove()
    T1.remove()
    T2.remove()
    T3.remove()
    T4.remove()
plt.close()

print('The simulation is complete. Consider the key questions.')

Key Questions: Consider these as you work your way through this page.

• Describe how does the thickness of the lithosphere change as it moves away from the ridge in time and distance. Follow one isotherm and observe how it changes. Is the rate of thickening constant? Where/when is it faster or slower?
• Why does the thickness change with distance from the ridge?
• Why is the lithosphere thicker at a distance 1000 km for u = 2 cm/yr compared to u = 5 cm/yr?