Edit on GitHub

Phase Diagrams

It is possible to employ phase diagrams as lookup tables, which can be used to compute density, melt fraction or seismic properties, for example.

Perple_X

A popular way to compute phase diagrams (pseudosections) for a given (fixed) chemical composition as a function of pressure and temperature is by using Perple_X. Once a stable assemblage is computed, through Gibbs energy minimisation, it allows you to output a large number of properties such as density, melt fraction or seismic velocities.

Within GeoParams, we can import this diagram (provided it is formatted in the same manner as LaMEM expects Perple_X input to be).

GeoParams.PerpleX_LaMEM_Diagram — Function
PD_Data = PerpleX_LaMEM_Diagram(fname::String; CharDim = nothing)

Reads a precomputed phase diagram in the LaMEM/Perple_X format (which is a phase diagram computed using Perple_X, but formatted in a manner that is readable using LaMEM). The data is stored in the PhaseDiagram_LookupTable structure.

If the CharDim object is specified, the values of all diagrams will be non-dimensionalized.

Example

julia> PD_Data = PerpleX_LaMEM_Diagram("./test_data/Peridotite.in")
Perple_X/LaMEM Phase Diagram Lookup Table: 
                      File    :   ./test_data/Peridotite.in
                      T       :   293.0 - 1573.000039
                      P       :   1.0e7 - 2.9999999944e9
                      fields  :   :meltRho, :meltRho, :meltFrac, :rockRho, :Rho, :rockVp
                                  :rockVs, :rockVpVs, :meltVp, :meltVs, :meltVpVs
                                  :Vp, :Vs, :VpVs

Once imported, the properties on the diagram can be interpolated in a simple manner:

julia> PD_Data.Rho(1500,1e7)
3042.836820256982

This also works for vectors or arrays:

julia> T = [1500 1800; 1233 1300]
julia> P = [1e8 1e9; 1e7 1e7]
julia> rho = PD_Data.Rho.(T,P)

(Note the dot . in front of the bracket while evaluating arrays).

The fields that are available depend on what is listed in the diagram file. The units of the fields are automatically evaluated, and employed to non-dimensionalize the parameters if CharDim is specified.

Algorithm

Internally, we employ linear interpolation, as provided by the Interpolations.jl package. Values outside the range of the diagram are set to the boundary of the diagram. The interpolation object is directly encoded in the PhaseDiagram_LookupTable` object.

source