API Reference¶

Complete API documentation for gdrift. All classes and functions are organized by functionality.

Core Data Loading¶

load_dataset¶

load_dataset ¶

load_dataset(dataset_name: str, table_names=[], return_metadata=False)

Load a dataset from local cache, downloading if necessary.

Args: dataset_name (str): Dataset name (without .h5 extension) table_names (list, optional): Specific tables to load. Defaults to []. return_metadata (bool, optional): Whether to return file-level metadata.

Returns: dict: dictionary with all the datasets (and optionally metadata tuple)

DatasetRegistry¶

DatasetRegistry ¶

Registry for managing and querying datasets.

Central repository of all registered datasets in gdrift. Provides methods for filtering, searching, and validating dataset metadata. Constructed from the datasets.json manifest at module import time.

Parameters:

Name	Type	Description	Default
`datasets`	`list of Dataset`	List of Dataset objects to register. Names must be unique.	required

Attributes:

Name	Type	Description
`_datasets`	`dict`	Internal dictionary mapping dataset names to Dataset objects.

Examples:

>>> from gdrift.datasetnames import DATASET_REGISTRY, DatasetType
>>> # Get all dataset names
>>> print(DATASET_REGISTRY.get_dataset_names())
>>> # Filter by type
>>> solidus_profiles = DATASET_REGISTRY.filter_by_type(
...     DatasetType.SOLIDUS_PROFILE)
>>> # Search by pattern
>>> slb_datasets = DATASET_REGISTRY.search_by_name("slb")
>>> # Check existence
>>> if "1d_prem" in DATASET_REGISTRY:
...     prem = DATASET_REGISTRY.get_dataset("1d_prem")

Notes

The module-level DATASET_REGISTRY singleton is initialized automatically by reading gdrift/datasets.json. Users should not need to create DatasetRegistry instances manually.

get_dataset ¶

get_dataset(name: str) -> Optional[Dataset]

Get a dataset by name.

1D Radial Profiles¶

SplineProfile¶

SplineProfile ¶

Bases: AbstractProfile

A class to represent a spline profile.

Attributes: raw_depth (Number): Array of depths. raw_value (Number): Array of corresponding values. spline_type (str): Type of spline to use. Defaults to "linear". _is_spline_made (bool): Flag to indicate if the spline has been created.

init ¶

__init__(depth: Number, value: Number, name: Optional[str] = 'Profile', spline_type: str = 'linear', extrapolate: bool = False)

Initialize a radial profile with spline interpolation.

Creates a 1D profile by interpolating between (depth, value) pairs using scipy.interpolate.interp1d. The spline is constructed lazily on first query for efficiency.

Parameters:

Name	Type	Description	Default
`depth`	`array_like`	Array of depth values in meters from the surface. Must be monotonically increasing or decreasing.	required
`value`	`array_like`	Array of property values corresponding to each depth. Must have the same length as depth.	required
`name`	`str`	Name of the profile (e.g., "density", "vs", "temperature"). Default is "Profile".	`'Profile'`
`spline_type`	`str`	Interpolation method passed to scipy.interpolate.interp1d. Options: "linear", "cubic", "quadratic", etc. Default is "linear".	`'linear'`
`extrapolate`	`bool`	Whether to allow extrapolation outside the depth range. If False, queries outside the range raise ValueError. If True, uses linear extrapolation. Default is False.	`False`

Notes

The spline is not created until the first call to at_depth() to avoid unnecessary computation during initialization.

at_depth ¶

at_depth(depth: Number) -> Number

Query the profile value at a specified depth or depths.

Evaluates the spline at the requested depth(s). Constructs the spline on first call if not already created.

Parameters:

Name	Type	Description	Default
`depth`	`float or ndarray`	Depth(s) in meters from the surface at which to query the profile. Can be a scalar or array.	required

Returns:

Type	Description
`float or ndarray`	Profile value(s) at the specified depth(s). Shape matches the input depth array.

Raises:

Type	Description
`ValueError`	If extrapolate=False and the provided depth is outside the valid range [min_depth, max_depth].

Notes

The spline is created lazily on the first call to this method using scipy.interpolate.interp1d with the specified spline_type.

min_max_depth ¶

min_max_depth()

Calculate the minimum and maximum depth values of the profile to prevent extrapolation.

tuple: A tuple containing the minimum and maximum depth values (min, max).

RadialEarthModel¶

RadialEarthModel ¶

Class representing reference Earth Models such as PREM or AK135 Composite object containing multiple radial profiles representing different Earth properties, such as shear wave velocity (Vs), primary wave velocity (Vp), and density. T

Attributes: depth_profiles (dict): A dictionary of RadialProfile instances.

init ¶

__init__(profiles: Union[AbstractProfile, List[AbstractProfile]])

Initialise the RadialEarthModel with a dictionary of radial profiles instances.

Args: profiles (dict of RadialProfile): Profiles for different properties, keyed by property name.

get_profile ¶

get_profile(property_name: str) -> AbstractProfile

Retrieve a profile by its name.

Args: property_name (str): The name of the property profile to retrieve.

Returns: AbstractProfile: The profile corresponding to the specified property name.

Raises: ValueError: If the specified property name does not exist in the model.

RadialEarthModelFromFile¶

RadialEarthModelFromFile ¶

Bases: RadialEarthModel

A class for loading radial profiles from a dataset.

This class extends SplineProfile to specifically handle loading, and utilizing available profiles related to profiles in the mantle.

Attributes: model_name (str): The name of the model/dataset from which profiles are loaded. description (str, optional): A brief description of the profile's purpose or characteristics.

init ¶

__init__(model_name: str, description: str = None)

Initialize a radial Earth model by loading profiles from an HDF5 dataset.

Loads all property profiles from a registered dataset file. The HDF5 file must contain a "depth" array and one or more property arrays (e.g., "density", "vs", "vp"). Each property is wrapped in a SplineProfile for interpolation.

Parameters:

Name	Type	Description	Default
`model_name`	`str`	Name of the registered dataset (e.g., "1d_prem", "1d_solidus_Andrault_et_al_2011_EPSL"). Must exist in the dataset registry.	required
`description`	`str`	Human-readable description of the model. If None, no description is set. Default is None.	`None`

Raises:

Type	Description
`ValueError`	If model_name is not in the dataset registry.
`KeyError`	If the HDF5 file does not contain a "depth" array.

Examples:

>>> import gdrift
>>> andrault = gdrift.RadialEarthModelFromFile(
...     "1d_solidus_Andrault_et_al_2011_EPSL")
>>> print(andrault.get_profile_names())
['solidus temperature']
>>> solidus = andrault.get_profile("solidus temperature")
>>> T = solidus.at_depth(410e3)

PreliminaryRefEarthModel¶

PreliminaryRefEarthModel ¶

Bases: RadialEarthModelFromFile

Initialises the Preliminary Reference Earth Model (PREM). This model is based on the work by Dziewonski and Anderson (1981) and provides a reference Earth model that can be queried at specific depths for various profiles.

References: Dziewonski, Adam M., and Don L. Anderson. "Preliminary reference Earth model." Physics of the Earth and Planetary Interiors 25.4 (1981): 297-356.

The object is of type RadialEarthModel and is initialized by loading profiles from an existing dataset. Each profile is represented as a SplineProfile object.

Attributes: prem_profiles (list): A list of SplineProfile objects representing different profiles in the PREM dataset.

init ¶

__init__()

Initialize the Preliminary Reference Earth Model (PREM).

Loads the PREM dataset (Dziewonski & Anderson, 1981) containing reference profiles for density, seismic velocities, elastic moduli, pressure, and gravity as a function of radius/depth.

The model spans from Earth's center to the surface and includes discontinuities at major boundaries (e.g., core-mantle boundary, 410 km, 660 km discontinuities).

Notes

PREM is the standard 1D reference model for seismology and geodynamics. It represents a spherically symmetric, non-rotating, oceanless Earth.

Examples:

>>> import gdrift
>>> prem = gdrift.PreliminaryRefEarthModel()
>>> print(prem.get_profile_names())
['density', 'vs', 'vp', ...]
>>> rho_profile = prem.get_profile("density")
>>> rho_670km = rho_profile.at_depth(670e3)

References

Dziewonski, A. M., & Anderson, D. L. (1981). Preliminary reference Earth model. Physics of the Earth and Planetary Interiors, 25(4), 297-356. https://doi.org/10.1016/0031-9201(81)90046-7

HirschmannSolidus¶

HirschmannSolidus ¶

Bases: RadialEarthModel

A class representing the solidus model based on the work of Hirschmann (2000).

Attributes: solidus_profile (HirschmannSolidusProfile): The solidus profile based on Hirschmann (2000).

init ¶

__init__()

Initialize a Hirschmann solidus radial Earth model.

Creates a RadialEarthModel containing a single profile: the Hirschmann (2000) solidus temperature. This is a convenience wrapper around HirschmannSolidusProfile that conforms to the RadialEarthModel interface.

The solidus can be accessed via get_profile("solidus temperature") or directly through the solidus_profile attribute.

Examples:

>>> import gdrift
>>> solidus_model = gdrift.HirschmannSolidus()
>>> solidus_profile = solidus_model.get_profile("solidus temperature")
>>> T_410km = solidus_profile.at_depth(410e3)
>>> print(f"Solidus at 410 km: {T_410km:.1f} K")

See Also

HirschmannSolidusProfile : The underlying solidus implementation

3D Earth Models¶

EarthModel3D¶

EarthModel3D ¶

Bases: AbstractEarthModel

Three-dimensional Earth model with KD-tree spatial interpolation.

A container for 3D gridded geophysical data with efficient spatial queries using scipy's cKDTree. Supports multiple interpolation kernels (IDW, Gaussian, Wendland, linear, cubic, nearest neighbor) and handles both Cartesian and geographic coordinate systems.

Typical workflow: 1. Create model: model = EarthModel3D(nearest_neighbours=8) 2. Set coordinates: model.set_coordinates(x, y, z) 3. Add quantities: model.add_quantity("dvs", dvs_array, "label") 4. Query: value = model.at(lat=45, lon=0, depth=500e3, quantity="dvs")

Parameters:

Name	Type	Description	Default
`nearest_neighbours`	`int`	Number of nearest neighbors to use for interpolation. Higher values produce smoother interpolation but increase computational cost. Default is 8 (suitable for most 3D tomography models).	`8`
`default_max_distance`	`float`	Maximum search distance in meters for finding neighbors. Points farther than this from all data points return NaN. Default is 200e3 (200 km), suitable for typical mantle convection resolution.	`200000.0`

Attributes:

Name	Type	Description
`coordinates`	`ndarray or None`	Nx3 array of (x, y, z) coordinates in normalized units (scaled to R_earth). Set via `set_coordinates()`. None until coordinates are set.
`available_fields`	`dict`	Dictionary mapping quantity names to (data_array, label) tuples. Populated via `add_quantity()`.
`tree`	`cKDTree or None`	KD-tree for fast spatial queries. Built automatically when coordinates are set. None until `set_coordinates()` is called.
`nearest_neighbours`	`int`	Number of neighbors for interpolation (from initialization).
`default_max_distance`	`float`	Maximum search distance in meters (from initialization).

Examples:

>>> import gdrift
>>> import numpy as np
>>> # Create synthetic 3D model
>>> model = gdrift.EarthModel3D(nearest_neighbours=8)
>>> x = np.linspace(-0.5, 0.5, 10)  # normalized coords
>>> y = np.linspace(-0.5, 0.5, 10)
>>> z = np.linspace(-0.5, 0.5, 10)
>>> xx, yy, zz = np.meshgrid(x, y, z, indexing='ij')
>>> coords = np.column_stack([xx.ravel(), yy.ravel(), zz.ravel()])
>>> model.set_coordinates(coords[:, 0], coords[:, 1], coords[:, 2])
>>> # Add velocity perturbation
>>> dvs = np.random.randn(len(coords)) * 0.02
>>> model.add_quantity("dvs", dvs, "dln(Vs)")
>>> # Query using geographic coordinates
>>> value = model.at(lat=45.0, lon=10.0, depth=1000e3, quantity="dvs")

See Also

SeismicModel : Subclass for loading seismic tomography models interpolate_to_points : Underlying interpolation engine

init ¶

__init__(nearest_neighbours=8, default_max_distance=200000.0)

Initialize a 3D Earth model with interpolation parameters.

Parameters:

Name	Type	Description	Default
`nearest_neighbours`	`int`	Number of nearest neighbors for interpolation. Default is 8.	`8`
`default_max_distance`	`float`	Maximum search distance in meters. Default is 200e3 (200 km).	`200000.0`

SeismicModel¶

SeismicModel ¶

Bases: EarthModel3D

init ¶

__init__(model_name, nearest_neighbours: int = 8, default_max_distance: float = 200000.0, labels=[])

SeismicModel is a class for handling 3D seismic models.

This class inherits from EarthModel3D and is used to load and manage seismic models. It allows for the initialization of a seismic model with a specified name, number of nearest neighbours, and a default maximum distance. The model data is loaded from a dataset file, and quantities and coordinates are set accordingly.

Parameters:

model_name : str Name of the seismic model to load nearest_neighbours : int, optional Number of nearest neighbours to consider for interpolation. Defaults to 8. default_max_distance : float, optional The default maximum distance for the model in meters. Defaults to 200e3. labels : list, optional Specific labels to load from the dataset. If empty, loads all available fields.

Attributes: model_name (str): The name of the seismic model. nearest_neighbours (int): The number of nearest neighbours to consider. default_max_distance (float): The default maximum distance for the model in meters.

Available Seismic Models¶

AVAILABLE_SEISMIC_MODELS `module-attribute` ¶

AVAILABLE_SEISMIC_MODELS = [(replace('3d_seismic_', '')) for ds in (filter_by_type(TOMOGRAPHY_MODEL))]

Thermodynamics¶

ThermodynamicModel¶

ThermodynamicModel ¶

Bases: object

Thermodynamic lookup table for mantle mineral physics properties.

Provides 2D interpolation of pre-computed mineral physics properties (density, seismic velocities, elastic moduli) as a function of depth and temperature. Tables are based on self-consistent thermodynamic databases (Stixrude & Lithgow-Bertelloni 2011, 2021) and stored in HDF5 format with bivariate spline interpolation.

The model supports: - Forward queries: (depth, temperature) → property (e.g., Vs, Vp, rho) - Inverse queries: (depth, velocity) → temperature - Multiple compositions and database versions (see MODELS_AVAIL, COMPOSITIONS_AVAIL) - On-the-fly computation of Vs/Vp from elastic moduli

Parameters:

Name	Type	Description	Default
`model`	`str`	Thermodynamic database version (e.g. "SLB_21"). See MODELS_AVAIL for all available versions, derived from the datasets.json manifest.	required
`composition`	`str`	Mantle composition. Available options depend on the model. See COMPOSITIONS_AVAIL for all available compositions (e.g. "pyroliteCFMAS", "pyroliteNCMAS", "bulk-oceanic-crustCFMAS").	required
`temps`	`array_like`	Temperature grid for subsampling (Kelvin). If None, uses full temperature range from dataset. Default is None.	`None`
`depths`	`array_like`	Depth grid for subsampling (meters). If None, uses full depth range from dataset. Default is None.	`None`
`extrapolate`	`bool`	Whether to allow extrapolation outside the table bounds. If False, queries outside the range raise ValueError. Default is False.	`False`

Attributes:

Name	Type	Description
`model`	`str`	Database version (e.g., "SLB_21").
`composition`	`str`	Mantle composition (e.g., "pyroliteCFMAS").
`extrapolate`	`bool`	Extrapolation flag.
`_tables`	`dict`	Dictionary mapping property names to Table objects. Keys are the property names from the HDF5 file (e.g., "rho", "bulk_mod", "shear_mod", "v_s", "v_p").

Examples:

>>> import gdrift
>>> # Load SLB 2021 pyrolite model
>>> tm = gdrift.ThermodynamicModel("SLB_21", "pyroliteCFMAS")
>>> # List available properties
>>> print(tm.available_tables())
['rho', 'bulk_mod', 'shear_mod', 'v_s', 'v_p', ...]
>>> # Forward query: temperature to Vs
>>> vs = tm.temperature_to_vs(temperature=1600, depth=500e3)
>>> print(f"Vs = {vs:.1f} m/s")
>>> # Inverse query: Vs to temperature
>>> T = tm.vs_to_temperature(vs=4500, depth=500e3)
>>> print(f"Temperature = {T:.1f} K")
>>> # Generic property lookup
>>> rho = tm.temperature_to_property("rho", temperature=1600, depth=500e3)

Notes

All depths are in meters from the surface
All temperatures are in Kelvin
Velocities are in m/s, densities in kg/m³, moduli in Pa
Not all model/composition combinations are available
Use available_tables() to see which properties are loaded
Vs and Vp are computed from moduli if not present in HDF5 file

See Also

regularise_thermodynamic_table : Smooth phase transitions apply_anelastic_correction : Convert to seismic-frequency velocities compute_swave_speed : Calculate Vs from shear modulus and density compute_pwave_speed : Calculate Vp from bulk/shear moduli and density

References

Stixrude, L., & Lithgow-Bertelloni, C. (2011). Thermodynamics of mantle minerals—II. Phase equilibria. Geophysical Journal International, 184(3), 1180-1213. https://doi.org/10.1111/j.1365-246X.2010.04890.x

Stixrude, L., & Lithgow-Bertelloni, C. (2021). Thermal expansivity, heat capacity and bulk modulus of the mantle. Geophysical Journal International. (In preparation for SLB_21 parameters)

init ¶

__init__(model: str, composition: str, temps=None, depths=None, extrapolate=False)

available_tables ¶

available_tables()

Return a list of available table names in this model.

temperature_to_property ¶

temperature_to_property(property_name, temperature, depth)

Convert temperature and depth to a material property value.

Args: property_name (str): Name of the property (e.g. 'rho', 'alpha', 'Cp', 'vs', 'vp'). temperature: Temperature value(s). depth: Depth value(s).

Returns: Interpolated property value(s) at the given temperature and depth. Returns NaN if extrapolate=False and inputs are out of bounds.

temperature_to_vs ¶

temperature_to_vs(temperature, depth)

Convert temperature and depth to shear wave velocity (Vs).

Convenience wrapper for temperature_to_property("vs", ...). Computes Vs from shear modulus and density using the relation: Vs = sqrt(shear_modulus / density).

Parameters:

Name	Type	Description	Default
`temperature`	`float or array_like`	Temperature in Kelvin.	required
`depth`	`float or array_like`	Depth in meters from the surface.	required

Returns:

Type	Description
`float or ndarray`	Shear wave velocity in m/s. Returns NaN for out-of-bounds queries if extrapolate=False.

Examples:

>>> import gdrift
>>> tm = gdrift.ThermodynamicModel("SLB_21", "pyroliteCFMAS")
>>> vs = tm.temperature_to_vs(temperature=1600, depth=500e3)
>>> print(f"Vs = {vs:.1f} m/s")

temperature_to_vp ¶

temperature_to_vp(temperature, depth)

Convert temperature and depth to compressional wave velocity (Vp).

Convenience wrapper for temperature_to_property("vp", ...). Computes Vp from bulk modulus, shear modulus, and density using: Vp = sqrt((bulk_modulus + 4/3 * shear_modulus) / density).

Parameters:

Name	Type	Description	Default
`temperature`	`float or array_like`	Temperature in Kelvin.	required
`depth`	`float or array_like`	Depth in meters from the surface.	required

Returns:

Type	Description
`float or ndarray`	Compressional wave velocity in m/s. Returns NaN for out-of-bounds queries if extrapolate=False.

Examples:

>>> import gdrift
>>> tm = gdrift.ThermodynamicModel("SLB_21", "pyroliteCFMAS")
>>> vp = tm.temperature_to_vp(temperature=1600, depth=500e3)
>>> print(f"Vp = {vp:.1f} m/s")

temperature_to_rho ¶

temperature_to_rho(temperature, depth)

Convert temperature and depth to density.

Convenience wrapper for temperature_to_property("rho", ...). Queries the pre-computed density lookup table.

Parameters:

Name	Type	Description	Default
`temperature`	`float or array_like`	Temperature in Kelvin.	required
`depth`	`float or array_like`	Depth in meters from the surface.	required

Returns:

Type	Description
`float or ndarray`	Density in kg/m³. Returns NaN for out-of-bounds queries if extrapolate=False.

Examples:

>>> import gdrift
>>> tm = gdrift.ThermodynamicModel("SLB_21", "pyroliteCFMAS")
>>> rho = tm.temperature_to_rho(temperature=1600, depth=500e3)
>>> print(f"Density = {rho:.1f} kg/m³")

vs_to_temperature ¶

vs_to_temperature(vs: Number, depth: Number, bounds: Optional[Union[Tuple[float, float], Tuple[ndarray, ndarray]]] = (300, 7000)) -> Number

Convert S-wave velocity (vs) to temperature at a given depth.

Parameters:

vs : Number The S-wave velocity. depth : Number The depth at which the temperature is to be calculated. bounds : Optional[Union[Tuple[float, float], Tuple[numpy.ndarray, numpy.ndarray]]], default=(300, 7000) The bounds for the temperature calculation. It can be a tuple of floats or numpy arrays.

Returns:

Number The temperature corresponding to the given S-wave velocity and depth.

vp_to_temperature ¶

vp_to_temperature(vp: Number, depth: Number, bounds: Optional[Union[Tuple[float, float], Tuple[ndarray, ndarray]]] = (300, 7000)) -> Number

Convert P-wave velocity (vp) to temperature at a given depth.

Parameters:

vp : Number The P-wave velocity. depth : Number The depth at which the temperature is to be calculated. bounds : Optional[Union[Tuple[float, float], Tuple[numpy.ndarray, numpy.ndarray]]], default=(300, 7000) The bounds for the temperature calculation. It can be a tuple of floats or numpy arrays.

Returns:

Number The temperature corresponding to the given P-wave velocity and depth.

Phase Transition Regularization¶

regularise_thermodynamic_table ¶

regularise_thermodynamic_table(slb_pyrolite: ThermodynamicModel, temperature_profile: AbstractProfile, regular_range: Dict[str, Tuple] = default_regular_range)

Regularises the thermodynamic table by creating a regularised thermodynamic model that uses precomputed regular tables for S-wave and P-wave speeds.

Args: slb_pyrolite (ThermodynamicModel): The original thermodynamic model. temperature_profile (AbstractProfile): The temperature profile to be used for regularisation. This is supposed to be a 1D profile of average temperature profiles. regular_range (Dict[str, Tuple], optional): Dictionary specifying the regularisation range for each parameter. Defaults to gdrift.mineralogy.default_regular_range.

Returns: RegularisedThermodynamicModel: A regularised thermodynamic model with precomputed tables for S-wave and P-wave speeds.

Anelasticity¶

CammaranoAnelasticityModel¶

CammaranoAnelasticityModel ¶

Bases: BaseAnelasticityModel

A specific implementation of an anelasticity model following the approach by Cammarano et al.

init ¶

__init__(B: Callable, g: Callable, a: Callable, solidus: SplineProfile, Q_bulk: Callable = lambda x: 10000, omega: Callable = lambda x: 1.0)

Initialize the model with the given parameters.

Args: B (Callable): Scaling factor for the Q model. g (Callable): Activation energy parameter. a (Callable): Frequency dependency parameter. solidus (SplineProfile): Solidus temperature profile for mantle. Q_bulk (Callable): Bulk quality factor (default is 10000). omega (Callable): Seismic frequency (default is 1).

from_q_profile `classmethod` ¶

from_q_profile(q_profile: str) -> CammaranoAnelasticityModel

Create a CammaranoAnelasticityModel from a predefined Q-profile.

Uses the six Q-profiles (Q1-Q6) from Cammarano et al. (2003), with upper/lower mantle parameter switching at 660 km depth. The Q factor is computed as:

Q = B * omega^a * exp(a * g * T_solidus / T)

where B controls the overall attenuation amplitude (related to grain size), g is a dimensionless activation energy parameter, a is the frequency exponent, and T_solidus is the solidus temperature at the given depth.

Q-profiles represent different assumptions about grain size, water content, and attenuation mechanisms in the mantle:

======= ============= ============= ====================================== Profile B (UM / LM) g (UM / LM) Physical Interpretation ======= ============= ============= ====================================== Q1 0.5 / 10 20 / 10 Fine grain size, low water content Q2 0.8 / 15 20 / 10 Medium grain size Q3 1.1 / 20 20 / 10 Coarse grain, higher water content Q4 0.035 / 2.25 30 / 15 SLB activation energy, fine grain Q5 0.056 / 3.6 30 / 15 SLB activation energy, medium grain Q6 0.077 / 4.95 30 / 15 SLB activation energy, coarse grain ======= ============= ============= ======================================

UM = upper mantle (< 660 km), LM = lower mantle (>= 660 km). Q1-Q3 use Cammarano et al. (2003) activation energy parameterization. Q4-Q6 use Stixrude & Lithgow-Bertelloni activation energy values. All profiles use a = 0.2 (frequency exponent) and omega = 1.0 Hz.

Args: q_profile (str): One of "Q1" through "Q6".

Returns: CammaranoAnelasticityModel: The configured model.

Raises: ValueError: If q_profile is not one of Q1-Q6.

GoesAnelasticityModel¶

GoesAnelasticityModel ¶

Bases: BaseAnelasticityModel

Anelasticity model following the approach by Goes et al. (2004). Uses parameters Q0 and xi instead of B and g.

init ¶

__init__(Q0: Callable, xi: Callable, a: Callable, omega: Callable, solidus: SplineProfile, Q_bulk: Callable = lambda x: 10000)

Initialize the Goes anelasticity model.

Args: Q0 (Callable): Reference quality factor as a function of depth. xi (Callable): Activation energy parameter as a function of depth. a (Callable): Frequency dependency parameter as a function of depth. omega (Callable): Seismic frequency as a function of depth. solidus (SplineProfile): Solidus temperature profile for mantle. Q_bulk (Callable): Bulk quality factor (default is 10000).

from_q_profile `classmethod` ¶

from_q_profile(q_profile: str) -> GoesAnelasticityModel

Create a GoesAnelasticityModel from a predefined Q-profile.

Uses the Goes et al. (2000) parameterization where the Q factor is computed as:

Q = Q0 * omega^a * exp(a * xi * T_solidus / T)

where Q0 is a reference quality factor, xi is a dimensionless activation energy parameter, a is the frequency exponent, and T_solidus is the solidus temperature at the given depth.

======= ============== ============= ====================================== Profile Q0 (UM / LM) xi (UM / LM) Physical Interpretation ======= ============== ============= ====================================== Q4 4.9 / 22.2 26 / 14 Standard Goes parameterization Q6 1.91 / 48.84 26 / 14 Alternative parameterization ======= ============== ============= ======================================

UM = upper mantle (< 660 km), LM = lower mantle (>= 660 km). All profiles use a = 0.15 (frequency exponent) and omega = 1.0 Hz.

Args: q_profile (str): One of "Q4" or "Q6".

Returns: GoesAnelasticityModel: The configured model.

Raises: ValueError: If q_profile is not Q4 or Q6.

apply_anelastic_correction¶

apply_anelastic_correction ¶

apply_anelastic_correction(thermo_model, anelastic_model: BaseAnelasticityModel)

Apply anelastic corrections to seismic velocity data using the provided "anelastic_model" within the low attenuation limit. The corrections are based on the equation: \(1 - \frac{V(anelastic)}{V(elastic)} = \frac{1}{2} \cot(\frac{\alpha \pi}{2}) Q^{-1}\) as described by Stixrude & Lithgow-Bertelloni (doi:10.1029/2004JB002965, Eq-10).

thermo_model (ThermodynamicModel): The thermodynamic model containing temperature and depth data.

ThermodynamicModel: A new thermodynamic model with anelastically corrected seismic velocities.

The returned model includes the following methods with anelastic corrections:

compute_swave_speed: Calculates the anelastic effect on shear wave speed.
compute_pwave_speed: Calculates the anelastic effect on compressional wave speed.

The compute_swave_speed method applies the anelastic correction to the shear wave speed using the provided anelastic model. The compute_pwave_speed method applies the anelastic correction to the compressional wave speed using the quality factor derived from the equations provided by Don L. Anderson & R. S. Hart (1978, PEPI eq 1-3).

The corrections are applied by meshing the depths and temperatures from the thermodynamic model and computing the quality factor matrices for shear and bulk moduli. The corrected seismic velocities are then calculated and returned as new tables with the corrected values.

Utility Functions¶

Coordinate Transformations¶

geodetic_to_cartesian ¶

geodetic_to_cartesian(lat, lon, depth, earth_radius=R_earth)

Convert geographic coordinates to Cartesian coordinates.

Parameters: lat (float or numpy.ndarray): Latitude in degrees. lon (float or numpy.ndarray): Longitude in degrees. depth (float or numpy.ndarray): Depth below Earth's surface in km. earth_radius (float): Radius of the Earth in km. Default is 6371 km.

Returns: tuple: Cartesian coordinates (x, y, z).

cartesian_to_geodetic ¶

cartesian_to_geodetic(x, y, z, earth_radius=R_earth)

Convert Cartesian coordinates to geographic coordinates.

Parameters: x (float or numpy.ndarray): x coordinate in km. y (float or numpy.ndarray): y coordinate in km. z (float or numpy.ndarray): z coordinate in km. earth_radius (float): Radius of the Earth in km. Default is 6371e3 m.

Returns: tuple: Geographic coordinates (lat, lon, depth).

nondimensionalise_coords ¶

nondimensionalise_coords(x, y, z, R_nd_earth=2.22, R_nd_cmb=1.22)

Convert dimensional Cartesian coordinates to nondimensional form.

Applies linear radial scaling to map physical coordinates (meters) to a nondimensional reference frame. Commonly used in geodynamic simulations to improve numerical conditioning.

Parameters:

Name	Type	Description	Default
`x`	`float or array_like`	Cartesian coordinates in meters (origin at Earth's center).	required
`y`	`float or array_like`	Cartesian coordinates in meters (origin at Earth's center).	required
`z`	`float or array_like`	Cartesian coordinates in meters (origin at Earth's center).	required
`R_nd_earth`	`float`	Nondimensional radius for Earth's surface. Default is 2.22.	`2.22`
`R_nd_cmb`	`float`	Nondimensional radius for core-mantle boundary. Default is 1.22.	`1.22`

Returns:

Type	Description
`tuple of (float or ndarray)`	Nondimensionalized (x', y', z') coordinates.

Notes

Uses linear scaling: r' = a*r + b, where a and b are determined by mapping R_earth → R_nd_earth and R_cmb → R_nd_cmb. Angular coordinates (theta, phi) are preserved.

See Also

dimensionalise_coords : Inverse transformation

dimensionalise_coords ¶

dimensionalise_coords(x, y, z, R_nd_cmb=1.22, R_nd_earth=2.22)

Convert nondimensional Cartesian coordinates back to meters.

Inverse of nondimensionalise_coords. Maps nondimensional coordinates back to physical units (meters) using linear radial scaling.

Parameters:

Name	Type	Description	Default
`x`	`float or array_like`	Nondimensional Cartesian coordinates.	required
`y`	`float or array_like`	Nondimensional Cartesian coordinates.	required
`z`	`float or array_like`	Nondimensional Cartesian coordinates.	required
`R_nd_cmb`	`float`	Nondimensional radius for core-mantle boundary. Must match the value used in nondimensionalisation. Default is 1.22.	`1.22`
`R_nd_earth`	`float`	Nondimensional radius for Earth's surface. Must match the value used in nondimensionalisation. Default is 2.22.	`2.22`

Returns:

Type	Description
`tuple of (float or ndarray)`	Dimensional (x, y, z) coordinates in meters.

Notes

Uses inverse linear scaling: r = a*r' + b, where a and b are determined by mapping R_nd_earth → R_earth and R_nd_cmb → R_cmb.

See Also

nondimensionalise_coords : Forward transformation

Gravity and Pressure¶

compute_gravity ¶

compute_gravity(radius, mass_enclosed)

Compute gravitational acceleration at each radius based on the enclosed mass.

Args: radius (numpy.ndarray): Array of radii from the center. mass_enclosed (numpy.ndarray): Array of cumulative mass enclosed up to each radius.

Returns: numpy.ndarray: Array of gravitational acceleration at each radius.

compute_pressure ¶

compute_pressure(radius, density, gravity)

Calculate the hydrostatic pressure at each radius based on the density and gravitational acceleration.

Args: radius (numpy.ndarray): Array of radii from the center to the surface. density (numpy.ndarray): Array of densities at each radius. gravity (numpy.ndarray): Array of gravitational accelerations at each radius.

Returns: numpy.ndarray: Array of pressures calculated from the surface inward to each radius.

compute_mass ¶

compute_mass(radius, density)

Compute the mass enclosed within each radius using the cumulative trapezoidal rule.

Args: radius (numpy.ndarray): Array of radii from the center of the Earth or other celestial body. density (numpy.ndarray): Array of densities corresponding to each radius.

Returns: numpy.ndarray: Array of cumulative mass enclosed up to each radius.

Mesh Generation¶

fibonacci_sphere ¶

fibonacci_sphere(n)

Generates points on a sphere using the Fibonacci sphere algorithm, which distributes points approximately evenly over the surface of a sphere.

This method calculates coordinates for each point using the golden angle, ensuring that each point is equidistant from its neighbors. The algorithm is particularly useful for creating evenly spaced points on a sphere's surface without clustering at the poles, a common issue in other spherical point distribution methods.

Args: n (int): The number of points to generate on the sphere's surface.

Returns: numpy.ndarray: A 2D array of shape (n, 3), where each row contains the [x, y, z] coordinates of a point on the sphere.

Example: >>> sphere = _fibonacci_sphere(100) >>> print(sphere.shape) (100, 3)

Constants¶

R_earth `module-attribute` ¶

R_earth = 6371000.0

R_cmb `module-attribute` ¶

R_cmb = R_earth - 2890000.0

Dataset Information¶

Dataset Types¶

All datasets in gdrift are categorized by type:

Reference Earth Model: 1D radial profiles (e.g., PREM)
Solidus Profile: Mantle solidus temperature vs. depth
Seismic Tomography Model: 3D seismic velocity anomalies
Thermodynamic Model: Temperature-to-property lookup tables
Geodynamic Profile: Pre-computed adiabatic profiles

Registered Datasets¶

Reference Earth Models¶

Name	Utility Class	Citation
`1d_prem`	`PreliminaryRefEarthModel`	Dziewonski & Anderson (1981)

Solidus Profiles¶

Name	Utility Class	Citation
`1d_solidus_Andrault_et_al_2011_EPSL`	`RadialEarthModelFromFile`	Andrault et al. (2011) EPSL
`1d_solidus_Fiquet_et_al_2010_SCIENCE`	`RadialEarthModelFromFile`	Fiquet et al. (2010) Science
`1d_solidus_Ghelichkhan_et_al_2021_GJI`	`RadialEarthModelFromFile`	Ghelichkhan et al. (2021) GJI
`1d_solidus_Nomura_et_al_2014_SCIENCE`	`RadialEarthModelFromFile`	Nomura et al. (2014) Science
`1d_solidus_Zerr_et_al_1988_SCIENCE`	`RadialEarthModelFromFile`	Zerr & Boehler (1994) Nature

Seismic Tomography Models¶

See AVAILABLE_SEISMIC_MODELS for the complete list of 25 tomography models including:

S40RTS (Ritsema et al. 2011)
GLAD-M25 (Lei et al. 2020)
SEMUCB-WM1 (French & Romanowicz 2014)
REVEAL (Thrastarson et al. 2024)
And 21 more...

Load with:

model = gdrift.SeismicModel("3d_seismic_<model_name>")

Thermodynamic Models¶

Name	Composition	Citation
`SLB_21_pyroliteCFMAS`	Pyrolite CFMAS	Stixrude & Lithgow-Bertelloni (2021)
`SLB_21_pyroliteNCMAS`	Pyrolite NCMAS	Stixrude & Lithgow-Bertelloni (2021)

Additional thermodynamic models (SLB_08, SLB_11, SLB_24) with various compositions are available. Use gdrift.mineralogy.MODELS_AVAIL and gdrift.mineralogy.COMPOSITIONS_AVAIL to see the full list.

Geodynamic Profiles¶

Name	Description	Citation
`1d_geodynamic_SLB21_pyroliteCFMAS`	Adiabatic profile for pyrolite	Computed from SLB_21

Complete Citations¶

Andrault et al. (2011)¶

Andrault, Denis, et al. "Solidus and liquidus profiles of chondritic mantle: Implication for melting of the Earth across its history." Earth and Planetary Science Letters 304.1-2 (2011): 251-259.

Dziewonski & Anderson (1981)¶

Dziewonski, Adam M., and Don L. Anderson. "Preliminary reference Earth model." Physics of the Earth and Planetary Interiors 25.4 (1981): 297-356.

Fiquet et al. (2010)¶

Fiquet, G., et al. "Melting of peridotite to 140 gigapascals." Science 329.5998 (2010): 1516-1518.

French & Romanowicz (2014)¶

French, Scott W., and Barbara Romanowicz. "Whole-mantle radially anisotropic shear velocity structure from spectral-element waveform tomography." Geophysical Journal International 199.3 (2014): 1303-1327.

Ghelichkhan et al. (2021)¶

Ghelichkhan, S., et al. "The adjoint method applied to time-dependent mantle convection." Geophysical Journal International (2021).

Lei et al. (2020)¶

Lei, W., et al. "Global adjoint tomography—model GLAD-M25." Geophysical Journal International 223.1 (2020): 1-21.

Nomura et al. (2014)¶

Nomura, Ryuichi, et al. "Low core-mantle boundary temperature inferred from the solidus of pyrolite." Science 343.6170 (2014): 522-525.

Ritsema et al. (2011)¶

Ritsema, Jeroen, et al. "S40RTS: a degree-40 shear-velocity model for the mantle from new Rayleigh wave dispersion, teleseismic traveltime and normal-mode splitting function measurements." Geophysical Journal International 184.3 (2011): 1223-1236.

Stixrude & Lithgow-Bertelloni (2011, 2016, 2021)¶

Stixrude, L., and C. Lithgow-Bertelloni. "Thermodynamics of mantle minerals - I. Physical properties." Geophysical Journal International 162.2 (2005): 610-632.
Stixrude, L., and C. Lithgow-Bertelloni. "Thermodynamics of mantle minerals - II. Phase equilibria." Geophysical Journal International 184.3 (2011): 1180-1213.
Updates in 2016 and 2021 (see SLB website for details).

Thrastarson et al. (2024)¶

Thrastarson, Solvi, et al. "REVEAL: A global full‐waveform inversion model." Bulletin of the Seismological Society of America 114.3 (2024): 1392-1406.

Zerr & Boehler (1994)¶

Zerr, A., and R. Boehler. "Constraints on the melting temperature of the lower mantle from high-pressure experiments on MgO and magnesioüstite." Nature 371.6497 (1994): 506-508.

Usage Examples¶

See the Examples section for Jupyter notebooks demonstrating these APIs in action.

API Reference¶

Core Data Loading¶

load_dataset¶

load_dataset ¶

DatasetRegistry¶

DatasetRegistry ¶

get_dataset ¶

1D Radial Profiles¶

SplineProfile¶

SplineProfile ¶

__init__ ¶

at_depth ¶

min_max_depth ¶

RadialEarthModel¶

RadialEarthModel ¶

__init__ ¶

get_profile ¶

RadialEarthModelFromFile¶

RadialEarthModelFromFile ¶

__init__ ¶

PreliminaryRefEarthModel¶

PreliminaryRefEarthModel ¶

__init__ ¶

HirschmannSolidus¶

HirschmannSolidus ¶

__init__ ¶

3D Earth Models¶

EarthModel3D¶

EarthModel3D ¶

__init__ ¶

SeismicModel¶

SeismicModel ¶

__init__ ¶

Available Seismic Models¶

AVAILABLE_SEISMIC_MODELS module-attribute ¶

Thermodynamics¶

ThermodynamicModel¶

ThermodynamicModel ¶

__init__ ¶

available_tables ¶

temperature_to_property ¶

temperature_to_vs ¶

temperature_to_vp ¶

temperature_to_rho ¶

vs_to_temperature ¶

vp_to_temperature ¶

Phase Transition Regularization¶

regularise_thermodynamic_table ¶

Anelasticity¶

CammaranoAnelasticityModel¶

CammaranoAnelasticityModel ¶

__init__ ¶

from_q_profile classmethod ¶

GoesAnelasticityModel¶

GoesAnelasticityModel ¶

__init__ ¶

from_q_profile classmethod ¶

apply_anelastic_correction¶

apply_anelastic_correction ¶

Utility Functions¶

Coordinate Transformations¶

geodetic_to_cartesian ¶

cartesian_to_geodetic ¶

nondimensionalise_coords ¶

dimensionalise_coords ¶

Gravity and Pressure¶

compute_gravity ¶

compute_pressure ¶

compute_mass ¶

Mesh Generation¶

fibonacci_sphere ¶

Constants¶

R_earth module-attribute ¶

R_cmb module-attribute ¶

Dataset Information¶

Dataset Types¶

Registered Datasets¶

Reference Earth Models¶

Solidus Profiles¶

Seismic Tomography Models¶

init ¶

init ¶

init ¶

init ¶

init ¶

init ¶

init ¶

AVAILABLE_SEISMIC_MODELS `module-attribute` ¶

init ¶

init ¶

from_q_profile `classmethod` ¶

init ¶

from_q_profile `classmethod` ¶

R_earth `module-attribute` ¶

R_cmb `module-attribute` ¶