example_lookup_tables

Here we show how to read in and use lookup tables using terratools.

Import required modules.

In [1]:

import numpy as np
import matplotlib.pyplot as plt
from terratools import lookup_tables
from terratools import example_data as examples

import numpy as np import matplotlib.pyplot as plt from terratools import lookup_tables from terratools import example_data as examples

Define a function to plot the multitables

In [2]:

def plot_multi(ax, xdats, ydat, xlabs, ylab):
    for xdat, xlab in zip(xdats, xlabs):
        ax.plot(xdat, ydat, label=xlab)
    ax.set_xlabel("V$_p$ (km/s)")
    ax.set_ylabel("Pressure (Pa)")
    ax.legend()
    plt.gca().invert_yaxis()

    return ax

def plot_multi(ax, xdats, ydat, xlabs, ylab): for xdat, xlab in zip(xdats, xlabs): ax.plot(xdat, ydat, label=xlab) ax.set_xlabel("V$_p$ (km/s)") ax.set_ylabel("Pressure (Pa)") ax.legend() plt.gca().invert_yaxis() return ax

Load in three example sesimic lokoup tables This creates interpolation objects that can be queried for the sesimic property at a temperature and pressure at each point.

In [3]:

hzb_table = lookup_tables.SeismicLookupTable(examples.example_hzb_table())
lhz_table = lookup_tables.SeismicLookupTable(examples.example_lhz_table())
bas_table = lookup_tables.SeismicLookupTable(examples.example_bas_table())

hzb_table = lookup_tables.SeismicLookupTable(examples.example_hzb_table()) lhz_table = lookup_tables.SeismicLookupTable(examples.example_lhz_table()) bas_table = lookup_tables.SeismicLookupTable(examples.example_bas_table())

Downloading file 'example_hzb_table.dat' from 'doi:10.6084/m9.figshare.24100362.v3/example_hzb_table.dat' to '/home/docs/.cache/terratools'.

Downloading file 'example_lhz_table.dat' from 'doi:10.6084/m9.figshare.24100362.v3/example_lhz_table.dat' to '/home/docs/.cache/terratools'.

Downloading file 'example_bas_table.dat' from 'doi:10.6084/m9.figshare.24100362.v3/example_bas_table.dat' to '/home/docs/.cache/terratools'.

In this example we will find Vp at a simple range of temperatures and pressures for two compostions (hzb (C=0) and lhz(C=0.2)). We then find the seismic properties for a compostion of C=0.15.

In [4]:

# A range of temperatures and pressures
temperatures = np.linspace(500, 3000, 32)
pressures = np.linspace(1e9, 9e10, 32)

# Vp at temperature/pressures
vp_hzb = hzb_table.interp_points(pressures, temperatures, "vp")
vp_lhz = lhz_table.interp_points(pressures, temperatures, "vp")

# Define new and known bulk compositions
Cnew = 0.15
Chzb = 0.0
Clhz = 0.2
# vp_Cnew is vp at the defined temperatures and pressures for C=0.15 mterial.
vp_Cnew = lookup_tables.linear_interp_1d(vp_hzb, vp_lhz, 0.0, 0.2, Cnew)

print(f"vp_Cnew = {vp_Cnew}")

# Plot
fig, ax = plt.subplots(figsize=(4, 7))
plot_multi(
    ax,
    [vp_hzb, vp_lhz, vp_Cnew],
    pressures,
    ["C=0", "C=0.2", "C=0.15"],
    "Pressure (Pa)",
)
plt.show()
plt.close()

# A range of temperatures and pressures temperatures = np.linspace(500, 3000, 32) pressures = np.linspace(1e9, 9e10, 32) # Vp at temperature/pressures vp_hzb = hzb_table.interp_points(pressures, temperatures, "vp") vp_lhz = lhz_table.interp_points(pressures, temperatures, "vp") # Define new and known bulk compositions Cnew = 0.15 Chzb = 0.0 Clhz = 0.2 # vp_Cnew is vp at the defined temperatures and pressures for C=0.15 mterial. vp_Cnew = lookup_tables.linear_interp_1d(vp_hzb, vp_lhz, 0.0, 0.2, Cnew) print(f"vp_Cnew = {vp_Cnew}") # Plot fig, ax = plt.subplots(figsize=(4, 7)) plot_multi( ax, [vp_hzb, vp_lhz, vp_Cnew], pressures, ["C=0", "C=0.2", "C=0.15"], "Pressure (Pa)", ) plt.show() plt.close()

vp_Cnew = [ 8.121325    8.51601475  8.73708843  8.94551319  9.75129172 10.22358645
 10.37688767 10.48186604 10.77025111 11.19774592 11.31918111 11.41340528
 11.50548894 11.59486526 11.68134804 11.76573668 11.84709592 11.92708421
 12.00501479 12.08146158 12.15644967 12.22973859 12.30186487 12.37235592
 12.44183714 12.50991387 12.57681491 12.64257477 12.70735105 12.77104765
 12.8337061  12.89545   ]

In this second example, we will add a third composition (bas (C=1.0)) and find the Seismic properties for a mechanical mixure of known proportions of hzb, lhz and bas.

In [5]:

# Fraction of each compostion that makes up the mechanical mixture
hzbfrac = 0.3
lhzfrac = 0.5
basfrac = 0.2

# Vp at same temps/pressures for basalt
vp_bas = bas_table.interp_points(pressures, temperatures, "vp")

# Calculate Vp of the mechanical mixture using harmonic mean
vp_mm = lookup_tables._harmonic_mean(
    [vp_bas, vp_hzb, vp_lhz], [basfrac, lhzfrac, hzbfrac]
)

fig, ax2 = plt.subplots(figsize=(4, 7))
plot_multi(
    ax2,
    [vp_hzb, vp_lhz, vp_bas, vp_mm],
    pressures,
    [
        f"hzb (frac={hzbfrac})",
        f"lhz (frac={lhzfrac})",
        f"bas (frac={basfrac})",
        "Mechanical Mixture",
    ],
    "Pressure (Pa)",
)
plt.show()
plt.close()

# Fraction of each compostion that makes up the mechanical mixture hzbfrac = 0.3 lhzfrac = 0.5 basfrac = 0.2 # Vp at same temps/pressures for basalt vp_bas = bas_table.interp_points(pressures, temperatures, "vp") # Calculate Vp of the mechanical mixture using harmonic mean vp_mm = lookup_tables._harmonic_mean( [vp_bas, vp_hzb, vp_lhz], [basfrac, lhzfrac, hzbfrac] ) fig, ax2 = plt.subplots(figsize=(4, 7)) plot_multi( ax2, [vp_hzb, vp_lhz, vp_bas, vp_mm], pressures, [ f"hzb (frac={hzbfrac})", f"lhz (frac={lhzfrac})", f"bas (frac={basfrac})", "Mechanical Mixture", ], "Pressure (Pa)", ) plt.show() plt.close()