Utilities¶
BurnMan has a number of low-level utilities to help achieve common goals. Several of these have already been described in previous sections.
Unit cell¶
- burnman.utils.unitcell.molar_volume_from_unit_cell_volume(unit_cell_v, z)[source]¶
Converts a unit cell volume from Angstroms^3 per unitcell, to m^3/mol.
- Parameters
- unit_cell_vfloat
Unit cell volumes [A^3/unit cell]
- zfloat
Number of formula units per unit cell
- Returns
- Vfloat
Volume [m^3/mol]
- burnman.utils.unitcell.cell_parameters_to_vectors(cell_parameters)[source]¶
Converts cell parameters to unit cell vectors.
- Parameters
- cell_parameters1D numpy array
An array containing the three lengths of the unit cell vectors [m], and the three angles [degrees]. The first angle (\(\alpha\)) corresponds to the angle between the second and the third cell vectors, the second (\(\beta\)) to the angle between the first and third cell vectors, and the third (\(\gamma\)) to the angle between the first and second vectors.
- Returns
- M2D numpy array
The three vectors defining the parallelopiped cell [m]. This function assumes that the first cell vector is colinear with the x-axis, and the second is perpendicular to the z-axis, and the third is defined in a right-handed sense.
- burnman.utils.unitcell.cell_vectors_to_parameters(M)[source]¶
Converts unit cell vectors to cell parameters.
- Parameters
- M2D numpy array
The three vectors defining the parallelopiped cell [m]. This function assumes that the first cell vector is colinear with the x-axis, the second is perpendicular to the z-axis, and the third is defined in a right-handed sense.
- Returns
- cell_parameters1D numpy array
An array containing the three lengths of the unit cell vectors [m], and the three angles [degrees]. The first angle (\(\alpha\)) corresponds to the angle between the second and the third cell vectors, the second (\(\beta\)) to the angle between the first and third cell vectors, and the third (\(\gamma\)) to the angle between the first and second vectors.
Mathematical¶
- burnman.utils.math.unit_normalize(a, order=2, axis=-1)[source]¶
Calculates the L2 normalized array of numpy array a of a given order and along a given axis.
- burnman.utils.math.linear_interpol(x, x1, x2, y1, y2)[source]¶
Linearly interpolate to point x, between the points (x1,y1), (x2,y2)
- burnman.utils.math.bracket(fn, x0, dx, args=(), ratio=1.618, maxiter=100)[source]¶
Given a function and a starting guess, find two inputs for the function that bracket a root.
- Parameters
- fnfunction
The function to bracket
- x0float
The starting guess
- dxfloat
Small step for starting the search
- argsparameter list
Additional arguments to give to fn
- ratio
The step size increases by this ratio every step in the search. Defaults to the golden ratio.
- maxiterint
The maximum number of steps before giving up.
- Returns
- xa, xb, fa, fb: floats
xa and xb are the inputs which bracket a root of fn. fa and fb are the values of the function at those points. If the bracket function takes more than maxiter steps, it raises a ValueError.
- burnman.utils.math.smooth_array(array, grid_spacing, gaussian_rms_widths, truncate=4.0, mode='inverse_mirror')[source]¶
Creates a smoothed array by convolving it with a gaussian filter. Grid resolutions and gaussian RMS widths are required for each of the axes of the numpy array. The smoothing is truncated at a user-defined number of standard deviations. The edges of the array can be padded in a number of different ways given by the ‘mode’ parameter.
- Parameters
- arraynumpy ndarray
The array to smooth
- grid_spacingnumpy array of floats
The spacing of points along each axis
- gaussian_rms_widthsnumpy array of floats
The Gaussian RMS widths/standard deviations for the Gaussian convolution.
- truncatefloat (default=4.)
The number of standard deviations at which to truncate the smoothing.
- mode{‘reflect’, ‘constant’, ‘nearest’, ‘mirror’, ‘wrap’, ‘inverse_mirror’}
The mode parameter determines how the array borders are handled either by scipy.ndimage.filters.gaussian_filter. Default is ‘inverse_mirror’, which uses
burnman.tools.math._pad_ndarray_inverse_mirror()
.
- Returns
- smoothed_array: numpy ndarray
The smoothed array
- burnman.utils.math.interp_smoothed_array_and_derivatives(array, x_values, y_values, x_stdev=0.0, y_stdev=0.0, truncate=4.0, mode='inverse_mirror', indexing='xy')[source]¶
Creates a smoothed array on a regular 2D grid. Smoothing is achieved using
burnman.tools.math.smooth_array()
. Outputs scipy.interpolate.interp2d() interpolators which can be used to query the array, or its derivatives in the x- and y- directions.- Parameters
- array2D numpy array
The array to smooth. Each element array[i][j] corresponds to the position x_values[i], y_values[j]
- x_values1D numpy array
The gridded x values over which to create the smoothed grid
- y_values1D numpy array
The gridded y_values over which to create the smoothed grid
- x_stdevfloat
The standard deviation for the Gaussian filter along the x axis
- y_stdevfloat
The standard deviation for the Gaussian filter along the x axis
- truncatefloat (optional)
The number of standard deviations at which to truncate the smoothing (default = 4.).
- mode{‘reflect’, ‘constant’, ‘nearest’, ‘mirror’, ‘wrap’, ‘inverse_mirror’}
The mode parameter determines how the array borders are handled either by scipy.ndimage.filters.gaussian_filter. Default is ‘inverse_mirror’, which uses
burnman.tools.math._pad_ndarray_inverse_mirror()
.- indexing{‘xy’, ‘ij’}, optional
Cartesian (‘xy’, default) or matrix (‘ij’) indexing of output. See numpy.meshgrid for more details.
- Returns
- interps: tuple of three interp2d functors
interpolation functions for the smoothed property and the first derivatives with respect to x and y.
- burnman.utils.math.compare_l2(depth, calc, obs)[source]¶
Computes the L2 norm for N profiles at a time (assumed to be linear between points).
- Parameters
depths (array of float) – depths. \([m]\)
calc (list of arrays of float) – N arrays calculated values, e.g. [mat_vs,mat_vphi]
obs (list of arrays of float) – N arrays of values (observed or calculated) to compare to , e.g. [seis_vs, seis_vphi]
- Returns
array of L2 norms of length N
- Return type
array of floats
- burnman.utils.math.compare_chifactor(calc, obs)[source]¶
Computes the chi factor for N profiles at a time. Assumes a 1% a priori uncertainty on the seismic model.
- Parameters
calc (list of arrays of float) – N arrays calculated values, e.g. [mat_vs,mat_vphi]
obs (list of arrays of float) – N arrays of values (observed or calculated) to compare to , e.g. [seis_vs, seis_vphi]
- Returns
error array of length N
- Return type
array of floats
- burnman.utils.math.l2(x, funca, funcb)[source]¶
Computes the L2 norm for one profile(assumed to be linear between points).
- Parameters
x (array of float) – depths \([m]\).
funca (list of arrays of float) – array calculated values
funcb (list of arrays of float) – array of values (observed or calculated) to compare to
- Returns
L2 norm
- Return type
array of floats
- burnman.utils.math.nrmse(x, funca, funcb)[source]¶
Normalized root mean square error for one profile :type x: array of float :param x: depths in m. :type funca: list of arrays of float :param funca: array calculated values :type funcb: list of arrays of float :param funcb: array of values (observed or calculated) to compare to
- Returns
RMS error
- Return type
array of floats
- burnman.utils.math.chi_factor(calc, obs)[source]¶
\(\chi\) factor for one profile assuming 1% uncertainty on the reference model (obs) :type calc: list of arrays of float :param calc: array calculated values :type obs: list of arrays of float :param obs: array of reference values to compare to
- Returns
\(\chi\) factor
- Return type
array of floats
- burnman.utils.math.independent_row_indices(array)[source]¶
Returns the indices corresponding to an independent set of rows for a given array. The independent rows are determined from the pivots used during row reduction/Gaussian elimination.
- Parameters
- array2D numpy array of floats
The input array
- Returns
- indices1D numpy array of integers
The indices corresponding to a set of independent rows of the input array.
- burnman.utils.math.array_to_rational_matrix(array)[source]¶
Converts a numpy array into a sympy matrix filled with rationals
- burnman.utils.math.generate_complete_basis(incomplete_basis, array)[source]¶
Given a 2D array with independent rows and a second 2D array that spans a larger space, creates a complete basis for the combined array using all the rows of the first array, followed by any required rows of the second array. So, for example, if the first array is: [[1, 0, 0], [1, 1, 0]] and the second array is: [[1, 0, 0], [0, 1, 0], [0, 0, 1]], the complete basis will be: [[1, 0, 0], [1, 1, 0], [0, 0, 1]].
- Parameters
- incomplete_basis2D numpy array
An array containing the basis to be completed.
- array2D numpy array
An array spanning the full space for which a basis is required.
- Returns
- complete_basis2D numpy array
An array containing the basis vectors spanning both of the input arrays.
Miscellaneous¶
- class burnman.utils.misc.OrderedCounter(**kwds)[source]¶
Bases:
Counter
,OrderedDict
Counter that remembers the order elements are first encountered
- clear() None. Remove all items from od. ¶
- copy()¶
Return a shallow copy.
- elements()¶
Iterator over elements repeating each as many times as its count.
>>> c = Counter('ABCABC') >>> sorted(c.elements()) ['A', 'A', 'B', 'B', 'C', 'C']
# Knuth’s example for prime factors of 1836: 2**2 * 3**3 * 17**1 >>> prime_factors = Counter({2: 2, 3: 3, 17: 1}) >>> product = 1 >>> for factor in prime_factors.elements(): # loop over factors … product *= factor # and multiply them >>> product 1836
Note, if an element’s count has been set to zero or is a negative number, elements() will ignore it.
- classmethod fromkeys(iterable, v=None)¶
Create a new ordered dictionary with keys from iterable and values set to value.
- get(key, default=None, /)¶
Return the value for key if key is in the dictionary, else default.
- items() a set-like object providing a view on D's items ¶
- keys() a set-like object providing a view on D's keys ¶
- most_common(n=None)¶
List the n most common elements and their counts from the most common to the least. If n is None, then list all element counts.
>>> Counter('abcdeabcdabcaba').most_common(3) [('a', 5), ('b', 4), ('c', 3)]
- move_to_end(/, key, last=True)¶
Move an existing element to the end (or beginning if last is false).
Raise KeyError if the element does not exist.
- pop(k[, d]) v, remove specified key and return the corresponding ¶
value. If key is not found, d is returned if given, otherwise KeyError is raised.
- popitem(/, last=True)¶
Remove and return a (key, value) pair from the dictionary.
Pairs are returned in LIFO order if last is true or FIFO order if false.
- setdefault(/, key, default=None)¶
Insert key with a value of default if key is not in the dictionary.
Return the value for key if key is in the dictionary, else default.
- subtract(**kwds)¶
Like dict.update() but subtracts counts instead of replacing them. Counts can be reduced below zero. Both the inputs and outputs are allowed to contain zero and negative counts.
Source can be an iterable, a dictionary, or another Counter instance.
>>> c = Counter('which') >>> c.subtract('witch') # subtract elements from another iterable >>> c.subtract(Counter('watch')) # subtract elements from another counter >>> c['h'] # 2 in which, minus 1 in witch, minus 1 in watch 0 >>> c['w'] # 1 in which, minus 1 in witch, minus 1 in watch -1
- update(**kwds)¶
Like dict.update() but add counts instead of replacing them.
Source can be an iterable, a dictionary, or another Counter instance.
>>> c = Counter('which') >>> c.update('witch') # add elements from another iterable >>> d = Counter('watch') >>> c.update(d) # add elements from another counter >>> c['h'] # four 'h' in which, witch, and watch 4
- values() an object providing a view on D's values ¶
- burnman.utils.misc.copy_documentation(copy_from)[source]¶
Decorator @copy_documentation(another_function) will copy the documentation found in a different function (for example from a base class). The docstring applied to some function a() will be
(copied from BaseClass.some_function): <documentation from BaseClass.some_function> <optionally the documentation found in a()>
- burnman.utils.misc.merge_two_dicts(x, y)[source]¶
Given two dicts, merge them into a new dict as a shallow copy.
- burnman.utils.misc.pretty_print_values(popt, pcov, params)[source]¶
Takes a numpy array of parameters, the corresponding covariance matrix and a set of parameter names and prints the parameters and principal 1-s.d.uncertainties (np.sqrt(pcov[i][i])) in a nice text based format.
- burnman.utils.misc.pretty_print_table(table, use_tabs=False)[source]¶
Takes a 2d table and prints it in a nice text based format. If use_tabs=True then only is used as a separator. This is useful for importing the data into other apps (Excel, …). The default is to pad the columns with spaces to make them look neat. The first column is left aligned, while the remainder is right aligned.
- burnman.utils.misc.array_from_file(filename)[source]¶
Generic function to read a file containing floats and commented lines into a 2D numpy array.
Commented lines are prefixed by the characters # or %.
- burnman.utils.misc.attribute_function(m, attributes, powers=[])[source]¶
Function which returns a function which can be used to evaluate material properties at a point. This function allows the user to define the property returned as a string. The function can itself be passed to another function (such as nonlinear_fitting.confidence_prediction_bands()).
Properties can either be simple attributes (e.g. K_T) or a product of attributes, each raised to some power.
- Parameters
- mMaterial
The material instance evaluated by the output function.
- attributeslist of strings
The list of material attributes / properties to be evaluated in the product
- powerslist of floats
The powers to which each attribute should be raised during evaluation
- Returns
- ——-
- ffunction(x)
Function which returns the value of product(a_i**p_i) as a function of condition (x = [P, T, V])
Tools¶
Burnman has a number of high-level tools to help achieve common goals. Several of these have already been described in previous sections.
Plotting¶
- burnman.tools.plot.pretty_plot()[source]¶
Makes pretty plots. Overwrites the matplotlib default settings to allow for better fonts. Slows down plotting
- burnman.tools.plot.plot_projected_elastic_properties(mineral, plot_types, axes, n_zenith=31, n_azimuth=91, n_divs=100)[source]¶
Plot types must be one of: ‘vp’: V$_{P}$ (km/s) ‘vs1’: ‘V$_{S1}$ (km/s) ‘vs2’: V$_{S2}$ (km/s) ‘vp/vs1’: V$_{P}$/V$_{S1}$ ‘vp/vs2’: V$_{P}$/V$_{S2}$ ‘s anisotropy’: S-wave anisotropy (%), 200(vs1s - vs2s)/(vs1s + vs2s)) ‘linear compressibility’ Linear compressibility (GPa$^{-1}$) ‘youngs modulus’: Youngs Modulus (GPa)
axes objects must be initialised with projection=’polar’
Output for seismology¶
- burnman.tools.output_seismo.write_tvel_file(planet_or_layer, modelname='burnmanmodel', background_model=None)[source]¶
Function to write input file for obspy travel time calculations. Note: Because density isn’t defined for most 1D seismic models, densities are output as zeroes. The tvel format has a column for density, but this column is not used by obspy for travel time calculations.
- Parameters
- planet_or_layerburnman.Planet() or burnman.Layer()
Planet or layer to write out to tvel file
- filenamestring
Filename to read to
- background_modelburnman.seismic.Seismic1DModel()
1D seismic model to fill in parts of planet (likely to be an earth model) that aren’t defined by layer (only need when using Layer())
- burnman.tools.output_seismo.write_axisem_input(layers, modelname='burnmanmodel_foraxisem', axisem_ref='axisem_prem_ani_noocean.txt', plotting=False)[source]¶
Writing velocities and densities to AXISEM (www.axisem.info) input file. The input can be a single layer, or a list of layers taken from a planet (planet.layers). Currently this function will implement explicit discontinuities between layers in the seismic model. Currently this function is only set for Earth.
- Parameters
- layerslist of one or more burnman.Layer()
List of layers to put in axisem file
- modelnamestring
Name of model, appears in name of output file
- axisem_refstring
Reference file, used to copy the header and for the rest of the planet, in the case of a Layer(), default = ‘axisem_prem_ani_noocean.txt’
- plotting: Boolean
True means plot of the old model and replaced model will be shown (default = False)
- burnman.tools.output_seismo.write_mineos_input(layers, modelname='burnmanmodel_formineos', mineos_ref='mineos_prem_noocean.txt', plotting=False)[source]¶
Writing velocities and densities to Mineos (https://geodynamics.org/cig/software/mineos/) input file Note: currently, this function only honors the discontinuities already in the synthetic input file, so it is best to only replace certain layers with burnman values
- Parameters
- layerslist of one or more burnman.Layer()
List of layers to put in axisem file
- modelnamestring
Name of model, appears in name of output file
- mineos_refstring
Reference file, used to copy the header and for the rest of the planet, in the case of a Layer(), default = ‘mineos_prem_noocean.txt’
- plotting: Boolean
True means plot of the old model and replaced model will be shown (default = False)
Equations of state¶
- burnman.tools.eos.check_eos_consistency(m, P=1000000000.0, T=300.0, tol=0.0001, verbose=False, including_shear_properties=True)[source]¶
Compute numerical derivatives of the gibbs free energy of a mineral under given conditions, and check these values against those provided analytically by the equation of state
- Parameters
- mmineral
The mineral for which the equation of state is to be checked for consistency
- Pfloat
The pressure at which to check consistency
- Tfloat
The temperature at which to check consistency
- tolfloat
The fractional tolerance for each of the checks
- verboseboolean
Decide whether to print information about each check
- including_shear_propertiesboolean
Decide whether to check shear information, which is pointless for liquids and equations of state without shear modulus parameterizations
- Returns
- consistency: boolean
If all checks pass, returns True
- burnman.tools.eos.check_anisotropic_eos_consistency(m, P=1000000000.0, T=2000.0, tol=0.0001, verbose=False)[source]¶
Compute numerical derivatives of the gibbs free energy of a mineral under given conditions, and check these values against those provided analytically by the equation of state
- Parameters
- mmineral
The mineral for which the equation of state is to be checked for consistency
- Pfloat
The pressure at which to check consistency
- Tfloat
The temperature at which to check consistency
- tolfloat
The fractional tolerance for each of the checks
- verboseboolean
Decide whether to print information about each check
- Returns
- consistency: boolean
If all checks pass, returns True