Layers and Planets#

Layer#

class burnman.Layer(name=None, radii=None, verbose=False)[source]#

Bases: object

The base class for a planetary layer. The user needs to set the following before properties can be computed:

  • set_material(), which sets the material of the layer, e.g. a mineral, solid_solution, or composite

  • set_temperature_mode(), either predefine, or set to an adiabatic profile

  • set_pressure_mode(), to set the self-consistent pressure (with user-defined option the pressures can be overwritten). To set the self-consistent pressure the pressure at the top and the gravity at the bottom of the layer need to be set.

  • make(), computes the self-consistent part of the layer and starts the settings to compute properties within the layer

Note that the entire planet this layer sits in is not necessarily self-consistent, as the pressure at the top of the layer is a function of the density within the layer (through the gravity). Entire planets can be computed self-consistently with the planet class. Properties will be returned at the pre-defined radius array, although the evaluate() function can take a newly defined depthlist and values are interpolated between these (sufficient sampling of the layer is needed for this to be accurate).

reset()[source]#

Resets all cached material properties. It is typically not required for the user to call this function.

set_material(material)[source]#

Set the material of a Layer with a Material

set_temperature_mode(temperature_mode='adiabatic', temperatures=None, temperature_top=None)[source]#

Sets temperatures within the layer as user-defined values or as a (potentially perturbed) adiabat.

Parameters:
  • temperature_mode (string) – This can be set to ‘user-defined’, ‘adiabatic’, or ‘perturbed-adiabatic’. ‘user-defined’ fixes the temperature with the profile input by the user. ‘adiabatic’ self-consistently computes the adiabat when setting the state of the layer. ‘perturbed-adiabatic’ adds the user input array to the adiabat. This allows the user to apply boundary layers (for example).

  • temperatures (array of float) – The desired fixed temperatures in [K]. Should have same length as defined radii in layer.

  • temperature_top (float) – Temperature at the top of the layer. Used if the temperature mode is chosen to be ‘adiabatic’ or ‘perturbed-adiabatic’. If ‘perturbed-adiabatic’ is chosen as the temperature mode, temperature_top corresponds to the true temperature at the top of the layer, and the reference isentrope at this radius is defined to lie at a temperature of temperature_top - temperatures[-1].

set_pressure_mode(pressure_mode='self-consistent', pressures=None, gravity_bottom=None, pressure_top=None, n_max_iterations=50, max_delta=1e-05)[source]#

Sets the pressure mode of the layer, which can either be ‘user-defined’, or ‘self-consistent’.

Parameters:
  • pressure_mode (string) – This can be set to ‘user-defined’ or ‘self-consistent’. ‘user-defined’ fixes the pressures with the profile input by the user in the ‘pressures’ argument. ‘self-consistent’ forces Layer to calculate pressures self-consistently. If this is selected, the user will need to supply values for the gravity_bottom [m/s^2] and pressure_top [Pa] arguments.

  • pressures (array of floats) – Pressures [Pa] to set layer to (if the ‘user-defined’ pressure_mode has been selected). The array should be the same length as the layers user-defined radii array.

  • pressure_top (float) – Pressure [Pa] at the top of the layer.

  • gravity_bottom (float) – Gravity [m/s^2] at the bottom of the layer.

  • n_max_iterations (integer) – Maximum number of iterations to reach self-consistent pressures.

  • max_delta (float) – Relative update to the highest pressure in the layer between iterations to stop iterations.

make()[source]#

This routine needs to be called before evaluating any properties. If pressures and temperatures are not user-defined, they are computed here. This method also initializes an array of copied materials from which properties can be computed.

evaluate(properties, radlist=None, radius_planet=None)[source]#

Function that is used to evaluate properties across the layer. If radlist is not defined, values are returned at the internal radlist. If asking for different radii than the internal radlist, pressure and temperature values are interpolated and the layer material evaluated at those pressures and temperatures.

Parameters:
  • properties (list of strings) – List of properties to evaluate.

  • radlist (array of floats) – Radii to evaluate properties at. If left empty, internal radii list is used.

  • planet_radius (float) – Planet outer radius. Used only to calculate depth.

Returns:

1D or 2D array of requested properties (1D if only one property was requested)

Return type:

numpy.array

property mass#

Calculates the mass of the layer [kg]

property moment_of_inertia#

Returns the moment of inertia of the layer [kg m^2]

property gravity#

Returns gravity profile of the layer [m s^(-2)]

property bullen#

Returns the Bullen parameter across the layer. The Bullen parameter assess if compression as a function of pressure is like homogeneous, adiabatic compression. Bullen parameter =1 , homogeneous, adiabatic compression Bullen parameter > 1 , more compressed with pressure, e.g. across phase transitions Bullen parameter < 1, less compressed with pressure, e.g. across a boundary layer.

property brunt_vasala#

Returns the brunt-vasala (or buoyancy) frequency, N, across the layer. This frequency assess the stabilty of the layer: N < 0, fluid will convect N= 0, fluid is neutral N > 0, fluid is stabily stratified.

property pressure#

Returns current pressures across the layer that was set with set_state().

Aliased with P().

Returns:

Pressures in [Pa] at the predefined radii.

Return type:

numpy.array

property temperature#

Returns current temperature across the layer that was set with set_state().

  • Aliased with T().

Returns:

Temperatures in [K] at the predefined radii.

Return type:

numpy.array

property molar_internal_energy#

Returns the molar internal energies across the layer.

Notes

  • Needs to be implemented in derived classes.

  • Aliased with energy().

Returns:

The internal energies in [J/mol] at the predefined radii.

Return type:

numpy.array

property molar_gibbs#

Returns the molar Gibbs free energies across the layer.

Needs to be implemented in derived classes. Aliased with gibbs().

Returns:

Gibbs energies in [J/mol] at the predefined radii.

Return type:

numpy.array

property molar_helmholtz#

Returns the molar Helmholtz free energies across the layer.

Needs to be implemented in derived classes. Aliased with helmholtz().

Returns:

Helmholtz energies in [J/mol] at the predefined radii.

Return type:

numpy.array

property molar_mass#

Returns molar mass of the layer.

Needs to be implemented in derived classes.

Returns:

Molar mass in [kg/mol].

Return type:

numpy.array

property molar_volume#

Returns molar volumes across the layer.

Needs to be implemented in derived classes. Aliased with V().

Returns:

Molar volumes in [m^3/mol] at the predefined radii.

Return type:

numpy.array

property density#

Returns the densities across this layer.

Needs to be implemented in derived classes. Aliased with rho().

Returns:

The densities of this material in [kg/m^3] at the predefined radii.

Return type:

numpy.array

property molar_entropy#

Returns molar entropies acroos the layer.

Needs to be implemented in derived classes. Aliased with S().

Returns:

Entropies in [J/K/mol] at the predefined radii.

Return type:

numpy.array

property molar_enthalpy#

Returns molar enthalpies across the layer.

Needs to be implemented in derived classes. Aliased with H().

Returns:

Enthalpies in [J/mol] at the predefined radii.

Return type:

numpy.array

property isothermal_bulk_modulus#

Returns isothermal bulk moduli across the layer.

Notes

  • Needs to be implemented in derived classes.

  • Aliased with K_T().

Returns:

Bulk moduli in [Pa] at the predefined radii.

Return type:

numpy.array

property adiabatic_bulk_modulus#

Returns the adiabatic bulk moduli across the layer.

Needs to be implemented in derived classes. Aliased with K_S().

Returns:

Adiabatic bulk modulus in [Pa] at the predefined radii.

Return type:

numpy.array

property isothermal_compressibility#

Returns isothermal compressibilities across the layer (or inverse isothermal bulk moduli).

Needs to be implemented in derived classes. Aliased with beta_T().

Returns:

Isothermal compressibilities in [1/Pa] at the predefined radii.

Return type:

numpy.array

property adiabatic_compressibility#

Returns adiabatic compressibilities across the layer (or inverse adiabatic bulk moduli).

Needs to be implemented in derived classes. Aliased with beta_S().

Returns:

Adiabatic compressibilities in [1/Pa] at the predefined radii.

Return type:

numpy.array

property shear_modulus#

Returns shear moduli across the layer.

Needs to be implemented in derived classes. Aliased with beta_G().

Returns:

Shear moduli in [Pa] at the predefined radii.

Return type:

numpy.array

property p_wave_velocity#

Returns P wave speeds across the layer.

Needs to be implemented in derived classes. Aliased with v_p().

Returns:

P wave speeds in [m/s] at the predefined radii.

Return type:

numpy.array

property bulk_sound_velocity#

Returns bulk sound speeds across the layer.

Needs to be implemented in derived classes. Aliased with v_phi().

Returns:

Bulk sound velocities in [m/s] at the predefined radii.

Return type:

numpy.array

property shear_wave_velocity#

Returns shear wave speeds across the layer.

Needs to be implemented in derived classes. Aliased with v_s().

Returns:

Shear wave speeds in [m/s] at the predefined radii.

Return type:

numpy.array

property grueneisen_parameter#

Returns the grueneisen parameters across the layer.

Needs to be implemented in derived classes. Aliased with gr().

Returns:

Grueneisen parameters [unitless] at the predefined radii.

Return type:

numpy.array

property thermal_expansivity#

Returns thermal expansion coefficients across the layer.

Needs to be implemented in derived classes. Aliased with alpha().

Returns:

Thermal expansivities in [1/K] at the predefined radii.

Return type:

numpy.array

property molar_heat_capacity_v#

Returns molar heat capacity at constant volumes across the layer.

Needs to be implemented in derived classes. Aliased with C_v().

Returns:

Heat capacities in [J/K/mol] at the predefined radii.

Return type:

numpy.array

property molar_heat_capacity_p#

Returns molar_heat capacity at constant pressures across the layer.

Needs to be implemented in derived classes. Aliased with C_p().

Returns:

Heat capacities in [J/K/mol] at the predefined radii.

Return type:

numpy.array

property P#

Alias for pressure()

property T#

Alias for temperature()

property energy#

Alias for molar_internal_energy()

property helmholtz#

Alias for molar_helmholtz()

property gibbs#

Alias for molar_gibbs()

property V#

Alias for molar_volume()

property rho#

Alias for density()

property S#

Alias for molar_entropy()

property H#

Alias for molar_enthalpy()

property K_T#

Alias for isothermal_bulk_modulus()

property K_S#

Alias for adiabatic_bulk_modulus()

property beta_T#

Alias for isothermal_compressibility()

property beta_S#

Alias for adiabatic_compressibility()

property G#

Alias for shear_modulus()

property v_p#

Alias for p_wave_velocity()

property v_phi#

Alias for bulk_sound_velocity()

property v_s#

Alias for shear_wave_velocity()

property gr#

Alias for grueneisen_parameter()

property alpha#

Alias for thermal_expansivity()

property C_v#

Alias for molar_heat_capacity_v()

property C_p#

Alias for molar_heat_capacity_p()

class burnman.BoundaryLayerPerturbation(radius_bottom, radius_top, rayleigh_number, temperature_change, boundary_layer_ratio)[source]#

Bases: object

A class that implements a temperature perturbation model corresponding to a simple thermal boundary layer. The model takes the following form: T = a*exp((r - r1)/(r0 - r1)*c) + b*exp((r - r0)/(r1 - r0)*c) The relationships between the input parameters and a, b and c are given below.

This model is a simpler version of that proposed by [RM81].

Parameters:
  • radius_bottom (float) – The radius at the bottom of the layer (r0) [m].

  • radius_top (float) – The radius at the top of the layer (r1) [m].

  • rayleigh_number (float) – The Rayleigh number of convection within the layer. The exponential scale factor is this number to the power of 1/4 (Ra = c^4).

  • temperature_change (float) – The total difference in potential temperature across the layer [K]. temperature_change = (a + b)*exp(c).

  • boundary_layer_ratio (float) – The ratio of the linear scale factors (a/b) corresponding to the thermal boundary layers at the top and bottom of the layer. A number greater than 1 implies a larger change in temperature across the top boundary than the bottom boundary.

temperature(radii)[source]#

Returns the temperature at one or more radii [K].

Parameters:

radii (float or numpy.array) – The radii at which to evaluate the temperature.

Returns:

The temperatures at the requested radii.

Return type:

float or numpy.array

dTdr(radii)[source]#

Returns the thermal gradient at one or more radii [K/m].

Parameters:

radii (float or numpy.array) – The radii at which to evaluate the thermal gradients.

Returns:

The thermal gradient at the requested radii.

Return type:

float or numpy.array

set_model_thermal_gradients(dTdr_bottom, dTdr_top)[source]#

Reparameterizes the model based on the thermal gradients at the bottom and top of the model.

Parameters:
  • dTdr_bottom (float) – The thermal gradient at the bottom of the model [K/m]. Typically negative for a cooling planet.

  • dTdr_top (float) – The thermal gradient at the top of the model [K/m]. Typically negative for a cooling planet.

Planet#

class burnman.Planet(name, layers, n_max_iterations=50, max_delta=1e-05, verbose=False)[source]#

Bases: object

A class to build (self-consistent) Planets made out of Layers (burnman.Layer). By default the planet is set to be self-consistent (with zero pressure at the surface and zero gravity at the center), but this can be overwritte using the set_pressure_mode(). Pressure_modes defined in the individual layers will be ignored. If temperature modes are already set for each of the layers, when the planet is initialized, the planet will be built immediately.

reset()[source]#

Resets all cached material properties. It is typically not required for the user to call this function.

get_layer(name)[source]#

Returns a layer with a given name

Parameters:

name (str) – Given name of a layer

Returns:

Layer with the given name.

Return type:

burnman.Layer

get_layer_by_radius(radius)[source]#

Returns a layer in which this radius lies

Parameters:

radius (float) – Radius at which to evaluate the layer.

Returns:

Layer in which the radius lies.

Return type:

burnman.Layer

evaluate(properties, radlist=None)[source]#

Function that is generally used to evaluate properties of the different layers and stitch them together. If asking for different radii than the internal radlist, pressure and temperature values are interpolated and the layer material evaluated at those pressures and temperatures.

Parameters:
  • properties (list of strings) – List of properties to evaluate

  • radlist (array of floats) – Radii to evaluate properties at. If left empty, internal radius lists are used.

Returns:

1D or 2D array of requested properties (1D if only one property was requested)

Return type:

numpy.array

set_pressure_mode(pressure_mode='self-consistent', pressures=None, pressure_top=0.0, gravity_bottom=0.0, n_max_iterations=50, max_delta=1e-05)[source]#

Sets the pressure mode of the planet by user-defined values are in a self-consistent fashion. pressure_mode is ‘user-defined’ or ‘self-consistent’. The default for the planet is self-consistent, with zero pressure at the surface and zero pressure at the center.

Parameters:
  • pressure_mode (str) – This can be set to ‘user-defined’ or ‘self-consistent’.

  • pressures (array of floats) – Pressures (Pa) to set layer to (‘user-defined’). This should be the same length as defined radius array for the layer.

  • pressure_top (float) – Pressure (Pa) at the top of the layer.

  • gravity_bottom (float) – Gravity (m/s^2) at the bottom the layer

  • n_max_iterations (int) – Maximum number of iterations to reach self-consistent pressures.

make()[source]#

This routine needs to be called before evaluating any properties. If pressures and temperatures are self-consistent, they are computed across the planet here. Also initializes an array of materials in each Layer to compute properties from.

property mass#

calculates the mass of the entire planet [kg]

property average_density#

calculates the average density of the entire planet [kg/m^3]

property moment_of_inertia#

#Returns the moment of inertia of the planet [kg m^2]

property moment_of_inertia_factor#

#Returns the moment of inertia of the planet [kg m^2]

property depth#

Returns depth of the layer [m]

property gravity#

Returns gravity of the layer [m s^(-2)]

property bullen#

Returns the Bullen parameter

property brunt_vasala#
property pressure#

Returns current pressure that was set with set_state().

Aliased with P().

Returns:

Pressure in [Pa].

Return type:

array of floats

property temperature#

Returns current temperature that was set with set_state().

Aliased with T().

Returns:

Temperature in [K].

Return type:

array of floats

property molar_internal_energy#

Returns the molar internal energy of the planet.

Needs to be implemented in derived classes. Aliased with energy().

Returns:

The internal energy in [J/mol].

Return type:

array of floats

property molar_gibbs#

Returns the molar Gibbs free energy of the planet.

Needs to be implemented in derived classes. Aliased with gibbs().

Returns:

Gibbs energy in [J/mol].

Return type:

array of floats

property molar_helmholtz#

Returns the molar Helmholtz free energy of the planet.

Needs to be implemented in derived classes. Aliased with helmholtz().

Returns:

Helmholtz energy in [J/mol].

Return type:

array of floats

property molar_mass#

Returns molar mass of the planet.

Needs to be implemented in derived classes.

Returns:

Molar mass in [kg/mol].

Return type:

array of floats

property molar_volume#

Returns molar volume of the planet.

Needs to be implemented in derived classes. Aliased with V().

Returns:

Molar volume in [m^3/mol].

Return type:

array of floats

property density#

Returns the density of this planet.

Needs to be implemented in derived classes. Aliased with rho().

Returns:

The density of this material in [kg/m^3].

Return type:

array of floats

property molar_entropy#

Returns molar entropy of the planet.

Needs to be implemented in derived classes. Aliased with S().

Returns:

Entropy in [J/K/mol].

Return type:

array of floats

property molar_enthalpy#

Returns molar enthalpy of the planet.

Needs to be implemented in derived classes. Aliased with H().

Returns:

Enthalpy in [J/mol].

Return type:

array of floats

property isothermal_bulk_modulus#

Returns isothermal bulk modulus of the planet.

Needs to be implemented in derived classes. Aliased with K_T().

Returns:

Isothermal bulk modulus in [Pa].

Return type:

array of floats

property adiabatic_bulk_modulus#

Returns the adiabatic bulk modulus of the planet.

Needs to be implemented in derived classes. Aliased with K_S().

Returns:

Adiabatic bulk modulus in [Pa].

Return type:

array of floats

property isothermal_compressibility#

Returns isothermal compressibility of the planet (or inverse isothermal bulk modulus).

Needs to be implemented in derived classes. Aliased with beta_T().

Returns:

Isothermal compressibility in [1/Pa].

Return type:

array of floats

property adiabatic_compressibility#

Returns adiabatic compressibility of the planet (or inverse adiabatic bulk modulus).

Needs to be implemented in derived classes. Aliased with beta_S().

Returns:

Adiabatic compressibility in [1/Pa].

Return type:

array of floats

property shear_modulus#

Returns shear modulus of the planet.

Needs to be implemented in derived classes. Aliased with beta_G().

Returns:

Shear modulus in [Pa].

Return type:

array of floats

property p_wave_velocity#

Returns P wave speed of the planet.

Needs to be implemented in derived classes. Aliased with v_p().

Returns:

P wave speed in [m/s].

Return type:

array of floats

property bulk_sound_velocity#

Returns bulk sound speed of the planet.

Needs to be implemented in derived classes. Aliased with v_phi().

Returns:

Bulk sound velocity in [m/s].

Return type:

array of floats

property shear_wave_velocity#

Returns shear wave speed of the planet.

Needs to be implemented in derived classes. Aliased with v_s().

Returns:

Shear wave speed in [m/s].

Return type:

array of floats

property grueneisen_parameter#

Returns the grueneisen parameter of the planet.

Needs to be implemented in derived classes. Aliased with gr().

Returns:

Grueneisen parameters [unitless].

Return type:

array of floats

property thermal_expansivity#

Returns thermal expansion coefficient of the planet.

Needs to be implemented in derived classes. Aliased with alpha().

Returns:

Thermal expansivity in [1/K].

Return type:

array of floats

property molar_heat_capacity_v#

Returns molar heat capacity at constant volume of the planet.

Needs to be implemented in derived classes. Aliased with C_v().

Returns:

Isochoric heat capacity in [J/K/mol].

Return type:

array of floats

property molar_heat_capacity_p#

Returns molar heat capacity at constant pressure of the planet.

Needs to be implemented in derived classes. Aliased with C_p().

Returns:

Isobaric heat capacity in [J/K/mol].

Return type:

array of floats

property P#

Alias for pressure()

property T#

Alias for temperature()

property energy#

Alias for molar_internal_energy()

property helmholtz#

Alias for molar_helmholtz()

property gibbs#

Alias for molar_gibbs()

property V#

Alias for molar_volume()

property rho#

Alias for density()

property S#

Alias for molar_entropy()

property H#

Alias for molar_enthalpy()

property K_T#

Alias for isothermal_bulk_modulus()

property K_S#

Alias for adiabatic_bulk_modulus()

property beta_T#

Alias for isothermal_compressibility()

property beta_S#

Alias for adiabatic_compressibility()

property G#

Alias for shear_modulus()

property v_p#

Alias for p_wave_velocity()

property v_phi#

Alias for bulk_sound_velocity()

property v_s#

Alias for shear_wave_velocity()

property gr#

Alias for grueneisen_parameter()

property alpha#

Alias for thermal_expansivity()

property C_v#

Alias for molar_heat_capacity_v()

property C_p#

Alias for molar_heat_capacity_p()