Equilibration

Overview

For a composite with fixed phases at a given pressure, temperature and composition, equilibrium is reached when the following relationships are satisfied:

\[0_i = R_{ij} \mu_j\]

where \(\mu_j\) are the chemical potentials of all of the endmembers in all of the phases, and \(R_{ij}\) is an independent set of balanced reactions between endmembers.

It is generally true that at a fixed composition, one can choose two equilibrium constraints (such as fixed temperature, pressure, entropy, volume, phase proportion or some composition constraint) and solve for the remaining unknowns. In BurnMan, this can be achieved using the burnman.equilibrate() function.

Implemented function

burnman.equilibrate(composition, assemblage, equality_constraints, free_compositional_vectors=[], tol=None, store_iterates=False, store_assemblage=True, max_iterations=100.0, verbose=False)[source]

A function that finds the thermodynamic equilibrium state of an assemblage subject to given equality constraints by solving a set of nonlinear equations related to the chemical potentials and other state variables of the system.

The user chooses an assemblage (e.g. olivine, garnet and orthopyroxene) and \(2+n_c\) equality constraints, where \(n_c\) is the number of bulk compositional degrees of freedom. The equilibrate function attempts to find the remaining unknowns that satisfy those constraints.

There are a number of equality constraints implemented in burnman. These are given as a list of lists. Each constraint should have the form: [<constraint type>, <constraint>], where <constraint type> is one of ‘P’, ‘T’, ‘S’, ‘V’, ‘X’, ‘phase_fraction’, or ‘phase_composition’. The format of the <constraint> object depends on the constraint type:

  • P: float or numpy.array of

    pressures [Pa]

  • T: float or numpy.array of

    temperatures [K]

  • S: float or numpy.array of

    entropies [J/K]

  • V: float or numpy.array of

    volumes [m:math:^3]

  • phase_fraction: tuple with the form (<phase> <fraction(s)>),

    where <phase> is one of the phase objects in the assemblage and <fraction(s)> is a float or numpy.array corresponding to the desired phase fractions.

  • phase_composition: tuple with the form (<phase> <constraint>),

    where <phase> is one of the phase objects in the assemblage and <constraint> has the form (site_names, n, d, v), where \((nx)/(dx) = v\), n and d are constant vectors of site coefficients, and v is a float or numpy.array. For example, a constraint of the form ([Mg_A, Fe_A], [1., 0.], [1., 1.], [0.5]) would correspond to equal amounts Mg and Fe on the A site (Mg_A / (Mg_A + Fe_A) = 0.5).

  • X: list of a numpy.array and a float or numpy.array,

    where the equality satisfies the linear equation arr[0] x = arr[1]. The first numpy.array is fixed, and the second is to be looped over by the equilibrate function. This is a generic compositional equality constraint.

Parameters:

  • composition (dict) – The bulk composition that the assemblage must satisfy.

  • assemblage (burnman.Composite or burnman.RelaxedComposite) – The assemblage to be equilibrated.

  • equality_constraints (list of list) – The list of equality constraints. See above for valid formats.

  • free_compositional_vectors (list of dict) – A list of dictionaries containing the compositional freedom of the solution. For example, if the list contains the vector {‘Mg’: 1., ‘Fe’: -1}, that implies that the bulk composition is equal to composition + \(a\) (n_Mg - n_Fe), where the value of \(a\) is to be determined by the solve. Vector given in atomic (molar) units of elements.

  • tol (float or numpy.array) – The tolerance for the nonlinear solver. This can be a single float, which is then applied to all parameters, or a numpy array with the same length as the number of parameters in the solve (2 + number of phases + number of free compositional vectors). The default is None, which sets the tolerances to 1 Pa for pressure, 1e-3 K for temperature, 1e-6 for phase amounts and 1e-6 for free compositional vectors. One of the termination criteria for the nonlinear solver is that the absolute change in each parameter is less than the corresponding tolerance.

  • store_iterates (bool) – Whether to store the parameter values for each iteration in each solution object.

  • store_assemblage (bool) – Whether to store a copy of the assemblage object in each solution object.

  • max_iterations (int) – The maximum number of iterations for the nonlinear solver.

  • verbose (bool) – Whether to print output updating the user on the status of equilibration.

Returns:

Solver solution object (or a list, or a 2D list of solution objects) created by burnman.optimize.nonlinear_solvers.damped_newton_solve(), and a namedtuple object created by burnman.tools.equilibration.get_equilibration_parameters(). See documentation of these functions for the return types. If store_assemblage is True, each solution object also has an attribute called assemblage, which contains a copy of the input assemblage with the equilibrated properties. So, for a 2D grid of solution objects, one could call either sols[0][1].x[0] or sols[0][1].assemblage.pressure to get the pressure.

Return type:

tuple