Perple_X Seismic Velocity Calculations

 

 


 

Contents

Introduction

Run-time options

Thermodynamic data files

Computational details

A Worked Example

Running Build: problem setup

Running Vertex: phase diagram section computation

Running Pssect: phase diagram section depiction

Running Werami: extraction of seismic velocities

Graphical Representation

Running Pstable

Matlab Plotting

Figures:

Figure 1: Perple_X program structure

Figure 2: Metabasalt phase relations

Figure 3: Metabasalt P-wave velocity


 

Introduction

 

This document details an application in which seismic velocities are recovered from a phase equilibrium model for a metabasalt as a function of pressure and temperature, i.e., an isochemical phase diagram section. The procedure to recover seismic velocities from a phase diagram section is no different than that to recover any other thermodynamic property with Perple_X; thus this document provides a generic template for the construction of phase diagram sections and recovery of thermodynamic properties. Readers who are not interested in the specifics of seismic wave-speed calculations should skip to the worked example section.

 

Recognition of the thermodynamic constraints on seismic wave-speeds has led some geophysicists to use free energy minimization to compute mineral modes and compositions as a function of bulk-rock composition, temperature, and pressure. These modes and compositions are then used to estimate seismic wave-speeds from semi-empirical models. This approach may be justified by the uncertainties involved, but it is thermodynamically inconsistent in that it amounts to using one thermodynamic database to establish the amounts of the minerals and a different thermodynamic database to compute the wave-speeds. For geophysicists who remain nonetheless enamored with empirical estimates of elastic modulii, Perple_X can be used to make calculations in this manner. However, I personally advocate using the same thermodynamic data for the computation of both the stable mineralogy and seismic velocities. The approach actually used in Perple_X is determined by the structure of the data base and options specified by the user as discussed below.

 

Run-time options

 

Perple_X programs read run-time options from a file with the default name perplex_option.dat, options of particular relevance to seismic velocity calculations are (see perplex_options for more details):

 

Thermodynamic data files

 

Typically, thermodynamic data bases for phase equilibrium calculations are designed for the computation of isostatic phase equilibria (i.e., for systems without differential stresses), however seismic wave propagation involves shear stresses. Therefore to compute seismic velocities with Perple_X conventional thermodynamic data bases must be expanded to define a shear modulus that characterizes the thermodynamic response to shear stress. Any Perple_X thermodynamic data file can be modified to include empirical functions for the adiabatic shear and bulk moduli (see thermodynamic_data_files). However, even if a data file includes empirical functions for the bulk modulus, the isostatic equation of state will be used to compute the bulk modulus unless the explicit_bulk_modulus option is set to T (true).

 

Data files that have been prepared for seismic velocity calculations assume either an isostatic or a non-isostatic equation of state. The data files hp02ver.dat and cr_hp02ver.dat are currently the only files in the first category, both being variations on the Holland & Powell (1998) thermodynamic data base. This data base is a conventional petrologic data base, its main virtue is that it can be used to predict the evolution of fluid- and melt-bearing silicate rocks. The Holland & Powell data base is designed primarily  for low pressures (<4 GPa), but performs reasonably well to pressures of ~12 GPa. It cannot be used at higher pressures because it lacks data for the phases that form at depths of the mantle transition zone. The shear moduli in these data files are ad-hoc and should be used with caution.

 

The second category of data files (sfo05ver.dat, stx08ver.dat, stx11ver.dat) are based on the Mie-Gruneisen formulations of Lars Stixrude and coworkers (Stixrude & Bukowinskii 1993; Stixrude & Lithgow-Bertelloni 2005a,b). These data files were developed specifically for seismic velocity calculations in mantle lithologies and can be used for the entire range of physical conditions realized in the Earth's mantle. The limitations of these data files are that they do not include the data necessary to model fluids or melts and that they assume a simplified chemical model (Na2O-CaO-MgO-FeO-Al2O3-SiO2).

 

For additional information, in particular sources, for the aforementioned data files see perplex_citation.

 

Computational details

 

As customary in seismic velocity computations, it is assumed that all phases are isotropic and therefore completely characterized by two independent elastic moduli. The density of any phase is computed as

 

ρ = N/diff(G,P),

 

where G is the molar Gibbs energy and N is the molar formula weight. For isostatic formulations, if the run-time option explicit_bulk_modulus = T, the adiabatic bulk modulus of any phase is computed as described in perplex_thermodynamic_data_files. Otherwise the adiabatic bulk modulus is computed as

 

KS = -diff(G,P)[diff(G,P,P)+diff(G,P,T)2/diff(G,T,T)].

 

The adiabatic shear moduli (μS) of individual phases are computed as described in perplex_thermodynamic_data_files, depending on whether the equation of state is isostatic or non-isostatic.

 

Seismic compressional and shear wave velocities through rock are determined by: the corresponding velocities within the constituent minerals; the amounts of the minerals; and rock texture. For a rock in thermodynamic equilibrium, the first two factors are determined entirely by thermodynamic properties and are computed automatically by Perple_X. The third factor, texture, determines the relation of the bulk rock seismic velocities to those of the constituent minerals. As texture is not a thermodynamic property, the averaging schemes used to compute bulk rock velocities are equivocal. Two popular averaging schemes, Hashin-Shtrikman (HS) and Voigt-Reuss-Hill (VRH) are supported in Perple_X (bounds, vrh/hs_weighting). Once the aggregate moduli have been computed, the compressional and shear wave-speeds are

 

vp = [(KS + 4μS/3)/ρ]1/2, vs = [μS/ρ]1/2

 

No provision is made within Perple_X to account for irreversible effects (dissipation).


 

A Worked Example

 

The example provided here (originally from Connolly & Kerrick 2002) illustrates the computation of seismic velocities through an altered metabasalt at P-T conditions relevant to subduction zones at depths of 90-240 km. The type of calculation in this example generates an isochemical phase diagram section. Such sections are commonly, albeit incorrectly, referred to as pseudosections (Appendix of Connolly & Petrini 2002).

 

The calculation and interpretation of a phase diagram section consists of five steps (Figure 1): the problem definition; the phase diagram section computation and depiction; computational assessment; digital interpretation; and graphical interpretation. Each step involves different programs that are described in the following sections:

 

 

Figure 1. Perple_X program structure for gridded minimization calculations (i.e., for phase diagram section and phase fractionation calculations). The file names given as examples are not prescribed. VERTEX, PSSECT and WERAMI request the name of the problem definition file generated by BUILD (e.g., my_project). PSTABLE and Perple_X MATLAB scripts request a file named by adding an integer suffix to this name (e.g., my_project_1), where the integer indicates the number of calculations made with WERAMI. All files generated by VERTEX are named by combining the root name of the problem definition file with different file type suffixes. For typical calculations, the files required by subsequent programs, in addition to the solution, thermodynamic and option files, are: my_project.dat (problem definition file), my_project.plt (plot file), my_project.blk (assemblage file), my_project.arf (auto-refine data file), and my_project.tof (auto-refine flag file).


 


 

Running BUILD

 

The following dialog details the prompts from the program BUILD (normal font) and user responses (bold font). Explanatory comments are written in red. Additional explanation on some prompts is available via the (help) links.

 

c\jamie\perp666> build

 

Perple_X version 6.6.5.9, compiled 5/8/2011.

 

NO is the default (blank) answer to all Y/N prompts (help)

 

Enter a name for this project (the name will be used as the root for all output file names) [default = my_project] (help)

in17a

 

Project names should not include blanks or "." characters, but they may include directory information. The project name can be up to 100 characters long, but because Perple_X output files are named using project name plus various suffixes (e.g., .dat, .plt, .prt, .arf, .tof, .tab) it is unwise to specify names that are longer than 92 characters.

 

The problem definition file will be named: in17a.dat

 

Enter thermodynamic data file name, left justified, [default = hp02ver.dat]: (help)

hp02ver.dat

 

The file hp02ver.dat is a version of the Holland and Powell (1998, revised 2002) data base modified (Connolly 2005) to include shear modulii relevant to the high pressure mineralogy of metabasalts. The hp02ver.dat file is suitable for calculations to pressures approaching 12 GPa. See above for a discussion of thermodynamic data bases for seismic wave-speed calculations. Sources and brief descriptions of commonly used files are at thermodynamic_data_files.

 

Enter the computational option file name, left justified, [default = perplex_option.dat]:
See: www.perplex.ethz.ch/perplex_options.html
perplex_option.dat
 
Reading computational options from: perplex_option.dat
 

Computational options used by Perple_X are controlled via the computational option file. This file is normally named perplex_option.dat. Provision is made here to specify an alternative name because it sometimes desirable to have customized option files for specific purposes (e.g., high resolution or solvus calculations).

 

The first major block of prompts from BUILD determine the manner in which the chemical components defined within the thermodynamic data file are to be used. VERTEX permits components to be treated in four different ways: saturated phase components are components whose chemical potentials are determined by the (assumed) stability of a phase, most commonly a fluid, containing these components; saturated components are components whose chemical potentials are determined by the (assumed) stability of a pure phase(s) and or solutions containing only the saturated-phase and saturated components; mobile components are components whose chemical potentials are defined explicitly; and thermodynamic components are components whose chemical potentials are the dependent variables of a phase diagram calculation. Phase diagram calculations require the specification of at least one thermodynamic component.

 

The current data base components are:

 NA2O  MGO   AL2O3 SIO2  K2O   CAO   TIO2  MNO   FEO   O2    H2O   CO2

Transform them (Y/N)? (help)

n

Calculations with a saturated FLUID (Y/N)? (help)

n

 

Here the system is not assumed to be saturated with respect to a fluid phase, rather fluid stability and composition is a function of the physical conditions. This closed system model requires that the fluid components are specified as thermodynamic components (below). There are problems with such a model, but the closed system model for water-CO2 bearing systems is certainly less arbitrary (and probably more realistic) than the petrological artifice of assuming that the system is saturated with respect to a fluid with a fixed composition. More realistic models for devolatilization require path-dependent models, such models can be implemented in Perple_X in phase fractionation calculation mode (and for phase fractionation with mass transfer see Connolly 2005). The main problem with the closed system model is that it implies that fluid generated by devolatilization remains within the system. This problem can be overcome by instructing WERAMI to compute the seismic properties of the solid aggregate (a limiting case).

 

Calculations with saturated components (Y/N)? (help)

n

Use chemical potentials, activities or fugacities as independent variables (Y/N)? (help)

n

Select thermodynamic components from the set: (help)

 NA2O  MGO   AL2O3 SIO2  K2O   CAO   TIO2  MNO   FEO   O2    H2O   CO2

Enter names, left justified, 1 per line, [cr] to finish:

SIO2

TIO2

AL2O3

FEO

MGO

CAO

NA2O

K2O

H2O

CO2

 

**warning ver016** you are going to treat a saturated (fluid) phase component as a thermodynamic component, this may not be what you want to do. (help)

 

Select fluid equation of state: (help)

 

   0 - X(CO2) Modified Redlich-Kwong (MRK/DeSantis/Holloway)

   1 - X(CO2) Kerrick & Jacobs 1981 (HSMRK)

   2 - X(CO2) Hybrid MRK/HSMRK

   3 - X(CO2) Saxena & Fei 1987 pseudo-virial expansion

   4 - Bottinga &Richet 1981 (CO2 RK)

   5 - X(CO2) Holland &Powell 1991, 1998 (CORK)

   6 - X(CO2) Hybrid Haar et al 1979/CORK (TRKMRK)

   7 - f(O2/CO2)-f(S2) Graphite buffered COHS MRK fluid

   8 - f(O2/CO2)-f(S2) Graphite buffered COHS hybrid-EoS fluid

   9 - Max X(H2O) GCOH fluid Cesare & Connolly 1993

  10 - X(O) GCOH-fluid hybrid-EoS Connolly & Cesare 1993

  11 - X(O) GCOH-fluid MRK Connolly & Cesare 1993

  12 - X(O)-f(S2) GCOHS-fluid hybrid-EoS Connolly & Cesare 1993

  13 - X(H2) H2-H2O hybrid-EoS

  14 - X(CO2) Pitzer & Sterner 1994; Holland & Powell mixing 2003

  15 - X(H2) low T H2-H2O hybrid-EoS

  16 - X(O) H-O HSMRK/MRK hybrid-EoS

  17 - X(O) H-O-S HSMRK/MRK hybrid-EoS

  18 - X(CO2) Delany/HSMRK/MRK hybrid-EoS, for P >10 kb

  19 - X(O)-X(S) COHS hybrid-EoS Connolly & Cesare 1993

  20 - X(O)-X(C) COHS hybrid-EoS Connolly & Cesare 1993

  21 - X(CO2) Halbach &Chatterjee 1982, P >10 kb, hybrid-Eos

  22 - X(CO2) DHCORK, hybrid-Eos

  23 - Toop-Samis Silicate Melt

  24 - f(O2/CO2)-N/C Graphite saturated COHN MRK fluid
  25 - H2O-CO2-NaCl Aranovich et al. 2010

 

5

 

The second major block of prompts defines the explicit computational variables.

 

Specify computational mode: (help)

1 - Unconstrained minimization
2 - Constrained minimization on a 2d grid [default]
3 - Constrained minimization on a 1d grid
4 - Output pseudocompound data
5 - Phase fractionation calculations

Use unconstrained minimization for Schreinemakers projections or phase diagrams with > 2 independent variables. Use constrained minimization for phase diagrams or phase diagram sections with < 3 independent variables.

2

 

Here we are concerned with a normal P-T section; however, in gridded minimization it is possible to construct sections as a function of any arbitrarily defined composition (variable X(C1) in this prompt) or along a geothermal gradient, for examples of such calculations refer to the T-X, P-X and X-X section and phase abundances along an arbitary path tutorials.

 

The data base has P(bars) and T(K) as default independent potentials. Make one dependent on the other, e.g., as along a geothermal gradient (y/n)? (help)

n

 

Select x-axis variable: (help)

 

     1 - P(bars)

     2 - T(K)

     3 - Composition X(C1)* (user defined)

 

*X(C1) can not be selected as the y-axis variable

2

 

 

Enter minimum and maximum values, respectively, for: T(K)

673 1473

Enter minimum and maximum values, respectively, for: P(bars)

30000 80000

 

The following blurb provides information on the resolution with which phase relations are mapped as a function of the section variables. This resolution is controlled via the grid_parameters keywords specified in the perplex_option.dat file. When the option auto_refine = auto (the default) calculations involve two iterations: an initial low resolution exploratory stage establishes the range of stable phase compositions; and a second high resolution auto-refine stage computes the final results.

 

For gridded minimization, grid resolution is determined by the number of levels (grid_levels) and the resolution at the lowest level in the X- and Y-directions (x_nodes and y_nodes) these parameters are currently set for the exploratory and autorefine cycles as follows:
 

stage        grid_levels  xnodes  ynodes    effective resolution
exploratory       1          20      20        20 x  20 nodes
auto-refine       4          60      60       473 x 473 nodes

 

To change these options edit or create the file perplex_option.dat
See: www.perplex.ethz.ch/perplex_options.html
 

Specify component amounts by weight (Y/N)?

y

Enter weight amounts of the components: (help)

 SIO2  TIO2  AL2O3 FEO   MGO   CAO   NA2O  K2O   H2O   CO2

for the bulk composition of interest:

45 1.1 15.26 9.84 6.54 12.66 2.03 0.55 2.63 2.9
 

Output a print file (Y/N)? (help)

n

 

The print file is usually not useful for gridded minimization, unless the user wants an explicit list of the phases considered in the calculation.

 

The final block of prompts from BUILD define the phases for the calculation. 

 

Exclude endmember phases (Y/N)? (help)

y

 

For inexperienced users it is best to begin by including all possible endmember phases (i.e., stoichiometric or pure phases and stoichiometric chemical potential buffers). This allows the user to identify flaws in her perception of what the stable phases should be and/or problems in the thermodynamic data.

 

The endmember phases are identified by abbreviated names, in general these abbreviations are defined in the header section of the thermodynamic data file, in a few cases separate lists of endmember abbreviations are available (perplex_thermodynamic_data_files).

 

Here, experience leads me to make six exclusions: Clinozoisite (cz) is excluded to suppress the zoisite-clinozoisite transition, which is primarily a distraction in terms of interpreting phase equilibria. The water endmember of the Holland & Powell melt model (h2oL) is excluded because the endmember is (incorrectly) more stable than water at high pressure and low temperature (this difficulty can also be circumvented by specifying the T_melt keyword in perplex_option.dat). The "aluminum-free chlorite" endmember (afchl) of the Baker, Holland & Powell chlorite model is excluded because the use of this endmember requires many more pseudocompounds to resolve chlorite compositions yet it has almost no effect on the computed phase relations. oilm is a special ilmenite endmember required for the Ilm(WPH) solution model, it should be excluded from any calculations that do not use the model. tiGL is the liquid TiO2 endmember for the PMELTs model, its stability field seems to excessive relative to Holland & Powell's data for Ti-bearing solids.

 

Do you want to be prompted for phases (Y/N)?

n

Enter names, left justified, 1 per line, press <enter> to finish:

cz

afchl

h2oL

oilm

tiGL

 

Include solution phases (Y/N)? (help)

y

Enter the solution model file name [default = solution_model.dat]: (help)

solution_model.dat
 
...diagnostics omitted... 
 

BUILD reads the solution model file and determines which solutions are consistent with the specified problem. The diagnostics (omitted here) can usually be ignored. They are useful if you discover that a solution you wish to include is not present in the list of solutions BUILD offers.

 

Select phases from the following list, enter 1 per line, press <enter> to finish: (help)

  Bio(TCC)    Opx(HP)     CrOpx(HP)   Gt(MPF)     oAmph(DP)   Omph(GHP)
  cAmph(DP)   Chl(HP)     O(SG)       melt(HP)    pMELTS(G)   Pheng(HP)
  Sapp(HP)    Sapp(KWP)   Osm(HP)     GaHcSp      F           F(salt)
  T           Scap        St(HP)      Ctd(HP)     Carp        hCrd
  Sud(Livi)   Sud         Cumm        Anth        Gl          Tr
  GlTrTsPg    Amph(DHP)   Amph(DPW)   feldspar    feldspar_B  Pl(h)
  Kf          San         San(TH)     MaPa        Mica(CF)    Mica(CHA1)
  Mica(CHA)   O(HP)       Cpx(l)      Cpx(h)      Mont        Do(HP)
  M(HP)       Do(AE)      Cc(AE)      Sp(JR)      Sp(GS)      Sp(HP)
  IlGkPy      Neph(FB)    GrPyAlSp(B  Gt(GCT)     Gt(HP)      Gt(EWHP)
  Gt(WPH)     A-phase     Chum        Atg         B           P
  OrFsp(C1)   AbFsp(C1)   Pl(I1,HP)   Fsp(C1)     oCcM(HP)    O(stx)
  Sp(stx)     Gt(stx)     Opx(stx)    Cpx(stx)    Cpx(stx7)   o-Amph
  Omph(HP)    CrGt        Cpx(HP)     CrSp        Ca-Amph(D)  Na-Amph(D)
  O(stx7)     Sp(stx7)    Sp(WPC)     Cpx(m)      Ol(m)       Pl(stx8)
  Sp(stx8)    O(stx8)     Opx(stx8)   Cpx(stx8)   Gt(stx8)    Pu(M)
  Act(M)      Stlp(M)     Mica(M)     Carp(M)     Sud(M)      GlTrTsMr
  Carb(M)     TiBio(HP)

For details on these models see: www.perplex.ethz.ch/perplex_solution_model_glossary.html or read the commentary in the solution model file.

 

F

Chl(HP)

Mica(CF)

San

Pl(h)

Gt(HP)

Do(HP)

M(HP)

Opx(HP)

Omph(GHP)

Bio(TCC)

T

 

Because fluid is not specified as a saturated phase it is necessary to specify explicitly that a fluid solution may be present (i.e., the solution model F).

 

Enter calculation title:

Staudigel closed system model

 

BUILD creates a problem definition file (in17a.dat), named by adding the suffix ".dat" to the project name, that summarizes the users input. To avoid the tedious procedure of running BUILD for minor modifications, the problem definition file can be edited in a text-editor to modify physicochemical conditions, solution model choices or other parameters.

Top


 

Running VERTEX

 

If the auto_refine option is on calculations are done in two stages: an initial low quality exploratory stage establishes the range of phase compositions; and a final high resolution auto-refine stage. These stages are executed in a single calculation, as done here, only if the value of the auto_refine keyword is auto. Before running VERTEX verify that the auto_refine keyword is set to auto in your copy of perplex_option.dat.

 

c\jamie\perp666> vertex

 

Perple_X version 6.6.5.9, compiled 5/14/2011.

 

Enter the project name (the name assigned in BUILD) [default = my_project]:

in17a

 

VERTEX begins by asking for a project name, i.e., the name specified by the first BUILD prompt.

 

Reading problem definition from file: in17a.dat
Reading thermodynamic data from file: hp02ver.dat
Reading solution models from file: solution_model.dat
Writing print output to file: not requested
Writing plot output to file: in17a.plt
Writing bulk composition plot output to file: in17a.blk

Reading computational options from: perplex_option.dat

Writing pseudocompound glossary to file: not requested
Writing auto refine summary to file: not requested
Writing computational option summary to file: not requested

 

Once the project name has been entered VERTEX lists the input/output files to be used/created for the calculation. Several potential output files are optional user information file, these include: the print file (requested in the problem definition file), a pseudocompound glossary that summarizes the static discretization points used for a calculation (requested with the pseudocompound_file keyword in perplex_option.dat), the auto refine summary, which lists the ranges of the stable solution phases (requested with the auto_refine_file keyword in perplex_option.dat), the computational option summary, which lists the options in perplex_option.dat that are relevant to the calculation (requested with the option_list_files keyword in perplex_option.dat).

 

The file summary is followed by an echo of the computational options for the particular calculation as specified in perplex_option.dat:

 

Perple_X computational option settings for VERTEX:

    Keyword:               Value:     Permitted values [default]:

  Auto-refine options:

    auto_refine            aut        off [manual] auto
    auto_refine_factor_I    2.0       >=1 [3]
    auto_refine_slop       0.010      [0.01]

  Free energy minimization options:

    zero_mode              0.1E-05    0->1 [1e-6]; < 0 => off
    zero_bulk              0.1E-05    0->1 [1e-6]; < 0 => off
    iteration limit         4         >0 [4]; iteration keyword value 1
    refinement factor       3         >2*r [12]; iteration keyword value 2
    refinement points       4         3->7 [4]; iteration keyword value 3
    reach_factor           0.9        >0.5 [0.9]

  2D grid options:

    x_nodes                40 / 60    [20/40], >0, <2048; effective x-resolution   40 / 473 nodes
    y_nodes                40 / 60    [20/40], >0, <2048; effective y-resolution   40 / 473 nodes
    grid_levels             1 / 4     [1/4], >0, <10
    linear_model             on       off [on]

  Solution subdivision options:

    initial_resolution     0.07       0->1 [0.1], 0 => off
    stretch_factor         0.016      >0 [0.0164]
    subdivision_override   off        [off] lin str
    hard_limits            off        [off] on

  Thermodynamic options:

    solvus_tolerance       aut        [aut] or 0->1; aut = automatic, 0 => p=c pseudocompounds, 1 => homogenize
    speciation_tolerance   0.1E-02    0->1 [1e-3]; order-disorder speciation precision
    T_stop (K)                0.0     [0]
    T_melt (K)              873.0     [873]
    order_check            on         off [on]
    approx_alpha           T          [T] F
    Anderson-Gruneisen     F          [T] F
    site_check             T          [T] F

  Input/Output options:

    dependent_potentials   on         off [on]
    logarithmic_p          F          [F] T
    bad_number           NaN          [0.0]

  Information file output options:

    option_list_files      F          [F] T; echo computational options
    pseudocompound_file    F          [F] T; echo static pseudocompound compositions
    auto_refine_file       F          [F] T; echo auto-refine compositions

  Worst case (Cartesian) compositional resolution (mol %):

    Exploratory stage:   0.098
    Auto-refine stage:   0.049

To change these options see: www.perplex.ethz.ch/perplex_options.html
 

All "make definitions" consistent with the problem definition file are output next. Make definitions define thermodynamic entities (e.g., buffers, endmembers or reactions) as a linear combination of independent entries in the thermodynamic data file. The make definitions are in the header section of the thermodynamic data file and can be modified, deleted or created for your own purposes, see thermodynamic data files for more information.

 

Summary of valid make definitions:

          MGO   AL2O3 SIO2  K2O   TIO2  FEO   H2O
tbi       1.00  0.50  1.00  0.50  1.00  0.00  1.00
tbi1      2.00  0.50  3.00  0.50  1.00  0.00  0.00
fbr       0.00  0.00  0.00  0.00  0.00  1.00  1.00
fchum     0.00  0.00  4.00  0.00  0.00  9.00  1.00
fphA      0.00  0.00  2.00  0.00  0.00  7.00  3.00
fper      0.00  0.00  0.00  0.00  0.00  1.00  0.00
fatg      0.00  0.00 34.00  0.00  0.00 48.00 31.00
sil8L     0.00  1.60  1.60  0.00  0.00  0.00  0.00
fo8L      4.00  0.00  2.00  0.00  0.00  0.00  0.00
fa8L      0.00  0.00  2.00  0.00  0.00  4.00  0.00
q8L       0.00  0.00  4.00  0.00  0.00  0.00  0.00
 

The remaining section of the initial console output (abridged below) consists of information and warnings on the solution models to be used in the calculation. These messages inform the user as to which solution models in the solution model file are compatible with the problem definition, i.e., that the solution end-members are possible phases for the calculation. If only a subset of the solution model end-members are present, the messages indicate how the model will be reformulated in terms of this subset. If you do not get the solution models you requested, or if you get different solution compositions than you expect, these messages may be helpful. The diagnostics also indicate the number of static pseudocompounds used for each solution, this number is useful for determining which solution models are costly in terms of computer memory.

**warning ver114** the following endmembers are missing for Bio(TCC) 
  ffbi_i    fbi

**warning ver050** reformulating reciprocal solution: Bio(TCC) because of missing endmembers.(reformulation can be controlled explicitly by excluding additional endmembers).

     338 pseudocompounds generated for: Bio(TCC) 

**warning ver114** the following endmembers are missing for Chl(HP)
  mame_i    mafchl_i  mnchl     fafchl_i  afchl

**warning ver050** reformulating reciprocal solution: Chl(HP)    because of missing endmembers. (reformulation can be controlled explicitly by excluding additional endmembers).

     118 pseudocompounds generated for: Chl(HP)
      63 pseudocompounds generated for: Omph(HP)
       9 pseudocompounds generated for: F
     118 pseudocompounds generated for: T

**warning ver114** the following endmembers are missing for Amph(DPW)
  mfets     ffets_i

**warning ver050** reformulating reciprocal solution: Amph(DPW)  because of missing endmembers. (reformulation can be controlled explicitly by excluding additional endmembers).

    3141 pseudocompounds generated for: Amph(DPW)
       9 pseudocompounds generated for: Pl(h)
       9 pseudocompounds generated for: San
    4288 pseudocompounds generated for: Mica(CH)
       9 pseudocompounds generated for: Do(HP)

**warning ver114** the following endmembers are missing for M(HP)
  rhc
       9 pseudocompounds generated for: M(HP)

**warning ver114** the following endmembers are missing for Gt(HP)
  spss
      63 pseudocompounds generated for: Gt(HP)
     118 pseudocompounds generated for: Opx(HP)

Total number of pseudocompounds:    8292

** Starting exploratory computational stage **

The preceding line terminates the initial section of the console output. Often this section is followed by warning messages such as:

 

**warning ver369** failed to converge at T= 673.00 K, P= 72105.3 bar Using Murnaghan EoS, probably for Ghiorso et al. MELTS/PMELTS endmember data. Volume estimated using 3rd order Taylor series.

 

that indicate unstable numerical behavior of specialized equations of state, most commonly these result because the calibration is being used beyond the range of physical conditions for which it was intended. In general, these warnings have no significant consequences, but if you have doubts it is possible to eliminate the offending entity (e.g., here endmembers of the pMELTS(G) silicate melt model are causing problems at subsolidus temperatures and pressures far beyond the 4 GPa upper limit for the pMELTS calibration).

 

During the computation VERTEX outputs progress information:

 

  5.0% done with low level grid.
 10.0% done with low level grid.
 15.0% done with low level grid.


...output abridged...


 85.0% done with low level grid.
 90.0% done with low level grid.
 95.0% done with low level grid.
100.0% done with low level grid.

 

Beginning grid refinement stage.

 

After completion of the exploratory stage, VERTEX indicates which solution phases were not stable and the compositional ranges of the stable solution phases, this information, which is saved in in17a.arf and echoed upon request in more readable form in in17a_auto_refine.txt, is used to increase for the subsequent auto-refine stage.

 

WARNING: The following solutions were input, but are not stable:

     Bio(TCC) 
     Chl(HP)
     Pl(h)
     Opx(HP)

Endmember compositional ranges for model: Omph(GHP)

     Endmember   Minimum   Maximum
     di          0.40153   0.82889
     jd          0.07778   0.52822

Endmember compositional ranges for model: F

     Endmember   Minimum   Maximum
     CO2         0.00000   0.33032

Site fraction ranges for multisite model: T

  Site 1
     Species   Minimum   Maximum
        1      0.91556   0.94222

  Site 2
     Species   Minimum   Maximum
        1      0.99444   0.99889

 

..output abridged...

 

At this point if the auto_refine keyword were set to manual, VERTEX would terminate, allowing the user to examine intermediate results, but requiring that VERTEX be run a second time to obtain the final high resolution calculation. However this example has been run with auto_refine set to auto, thus the VERTEX automatically begins the auto-refine stage of the calculation:

 

Reading auto-refinement data from file: in17a.arf

 

At the beginning of auto-refine stage VERTEX eliminates all solution phases specified in the problem definition file (in17a.dat) that were not stable in the exploratory calculation and then restricts the compositional ranges for the remaining solutions to those found in the exploratory calculation (plus an increment [auto_refine_slop] to allow for inaccuracy), additionally the compositional resolution of the calculation in these restricted compositional resolution is increased by the auto_refine_factor.

 

Eliminating solution model: Chl(HP)    in auto-refinement.

Eliminating solution model: Pl(h)      in auto-refinement.

Eliminating solution model: Opx(HP)    in auto-refinement.

Eliminating solution model: Bio(TCC)   in auto-refinement.

     148 pseudocompounds generated for: Omph(GHP)
      11 pseudocompounds generated for: F
       2 pseudocompounds generated for: T
      31 pseudocompounds generated for: San
  164205 pseudocompounds generated for: Mica(CH)
       4 pseudocompounds generated for: Do(HP)
      14 pseudocompounds generated for: M(HP)
     117 pseudocompounds generated for: Gt(HP)

Total number of pseudocompounds:  164532

** Starting auto_refine computational stage **

During auto-refinement VERTEX may write diagnostics such as:

 

WARNING: composition of solution Do(HP) has reached an internal limit (0.759) on site 1 for species 1.

If this warning occurs during auto-refinement, the problem can be circumvented by increasing the auto_refine_slop specified in perplex_option.dat.

 

For the remainder of this calculation the internal limit has been relaxed to: 0.726

 

that commonly indicates that the composition of a phase has exceeded the range discovered in a previous exploratory stage. When this occurs VERTEX relaxes the range as indicated, which usually prevents any significant consequences.

 

..output abridged...

 

 98.3% done with low level grid.

100.0% done with low level grid.

 

In contrast to the exploratory stage, the grid_parameter keyword specifies a multi-level (i.e., 4-level) grid for the auto-refine stage, thus after the calculation of the first (or "low level") grid, VERTEX makes three refinements of the grid, each of these refinements takes roughly twice the time of the previous refinement.

 

Beginning grid refinement stage.

   649 grid cells to be refined at grid level 2

       ...working (   501 minimizations done)

       ...working (  1003 minimizations done)

       refinement at level 2 involved   1376 minimizations

  4976 minimizations required of the theoretical limit of  14161

  1202 grid cells to be refined at grid level 3

       ...working (   128 minimizations done)

       ...working (   629 minimizations done)

       ...working (  1131 minimizations done)

       ...working (  1634 minimizations done)
 

..output abridged...

 

VERTEX terminates as in the exploratory cycle with a summary of the phase compositions and writes a new auto-refine data file.

 

Top


Running PSSECT

PSSECT generates PostScript graphical output and can be used once the phase diagram section has been calculated with VERTEX. The default plot from PSSECT shows the computed phase relations (Figure 2) and is generated by the following dialog.

 

c\jamie\perp666> pssect

 

Perple_X version 6.6.5.9, compiled 5/14/2011.

Enter the project name (the name assigned in BUILD) [default = my_project]:

in17a

 

Perple_X plotting programs (PSSECT, PSTABLE, PSVDRAW, PSPTS) read options that control the graphical output (image size, orientation, etc) from a file named perplex_plot_option.dat, these options are summarized after the project name is read. PSSECT also reads the computational option file perplex_option.dat, but only to determine the parameters used by VERTEX for mapping by gridded minimization.

 

Reading computational options from: perplex_option.dat

 
Perple_X plot options are currently set as:

    Keyword:               Value:     Permitted values [default]:
    axis_label_scale       1.20       [1.2] (rel)
    bounding_box :
                               0      [0] x-min (pts)
                               0      [0] y-min (pts)
                             800      [800] x-length (pts)
                             800      [800] y-length (pts)
    field_fill             T          [T] F
    field_label            T          [T] F
    field_label_scale      0.75       [0.72] (rel)
    font                   Helvetica
    grid                   F          [F] T
    half_ticks             T          [T] F
    line_width             1.00       0-99 [1.] (pts)
    picture_transformation :
                            0.180     [0.18] x-scale (rel)
                            0.180     [0.18] y-scale (rel)
                             130.     [0.18] x-translation (pts)
                             220.     [0.18] y-translation (pts)
                             0.00     [0.0]  rotation (deg)
    plot_aspect_ratio      1.000      [1.0] x_axis_length/y_axis_length
    replicate_label        0.250      0->1 [0.025]
    splines                T          [T] F
    tenth_ticks            F          [F] T
    text_scale             1.000      [1.] (rel)

To change these options edit or create the plot option file
See: www.perplex.ethz.ch/perplex_plot_options.html

 

The graphical output file (in17a.ps) from PSSECT is (in fact) encapsulated PostScript.

 

PostScript will be written to file: in17a.ps

 

Answer yes to the following prompt to: specify rational numerical labels for the axes, restrict the range of the independent variables, or to show equilibria involving specific phases or phase assemblages.

 

Modify the default plot (y/n)?

n

 

The following diagnostic is issued because, if a phase assemblage is stable in more than one field of a section, then the fields may not be labeled. Whether multiple fields for the same assemblage are individually labeled is controlled by the replicate_label keyword in perplex_plot_options.dat.

 

There are 2 fields for: Omph(GHP) F Mica(CF) Do(HP) Gt(HP) coe ru arag kcym

 

Once PSSECT finishes, the output file (plot17a.ps) can be viewed with any graphics program (e.g., Adobe Illustrator, CorelDraw, Ghostview) or printer capable of interpreting encapsulated PostScript (Figure 2).

 

 
 
 
 
 
 
 
 

Figure 2 Computed metabasalt phase relations (in17a.ps), this figure differs from Figure 2a of Connolly & Kerrick (2002) because the published figure was calculated with different solution models. To simplify the plot the solution model abbreviations have been shortened (e.g., Omph(GHP) has been changed to Omph) and some field labels have been eliminated.
 

Top


Running WERAMI

 

WERAMI permits the user to extract information from a phase diagram section at a point (operational mode 1), throughout a region (mode 2), or along a line or curve (modes 3 and 4). The latter two options generate output that can be plotted with PSTAB, PyWerami or general-purpose programs such as Matlab. The dialog and output for each of these modes is detailed in the pseudosection tutorial. Here the second mode is used to extract P-wave velocities for the metabasalt over the entire range of pressure-temperature conditions represented by the phase diagram section.

 

c\jamie\perp666> werami

Perple_X version 6.6.5.9, compiled 5/14/2011.

Enter the project name (the name assigned in BUILD) [default = my_project]:

in17a

 

Although the phase equilibria analyzed by WERAMI are computed by VERTEX, seismic wave-speeds are only computed in WERAMI. Consequently, the computational options that influence wave-speeds can be changed without necessitating that the calculation be repeated with VERTEX. For example, a weakness of the thermodynamic formulation (Holland & Powell 1998) used for the present calculation is that it sometimes results in non-physical expansivities. This weakness is of no major consequence for phase equilibria computations, but can be problematic for the calculation of seismic or physical properties at elevated pressure. A simple remedy (Helffrich & Connolly 2009) for this problem is to change the value of the Anderson-Gruneisen keyword in perplex_option.dat from F to T after running VERTEX.

 

Reading computational options from: perplex_option.dat
Writing computational option summary to file: not requested

Perple_X computational option settings for WERAMI:

    Keyword:               Value:     Permitted values [default]:

  Input/Output options:

    spreadsheet            T          [F] T
    logarithmic_p          F          [F] T
    bad_number           NaN          [0.0]
    compositions           mol        wt  [mol]
    proportions            vol        wt  [vol] mol
    interpolation          on         off [on ]
    extrapolation          off        on  [off]
    melt_is_fluid          F          [F] T

  Information file output options:

    option_list_files      F          [F] T; echo computational options

  Thermodynamic options:

    approx_alpha           T          [T] F
    Anderson-Gruneisen     F          [T] F

  Seismic velocity options:

    bounds                 VRH        HS  [VRH]
    vrh/hs_weighting       0.5        0->1 [0.5]
    explicit_bulk_modulus  T          [F] T
    poisson_ratio          on         off [on ] all; Poisson ratio = 0.35

To change these options see: www.perplex.ethz.ch/perplex_options.html


Select operational mode:
    1 - properties at specified conditions
    2 - properties on a 2d grid
    3 - properties along a 1d path
    4 - as in 3, but input from file
    0 - EXIT

2

 

Select properties [enter 0 to finish]:
    1 - Specific Enthalpy (J/m3)
    2 - Density (kg/m3)
    3 - Specific heat capacity (J/K/m3)
    4 - Expansivity (1/K, for volume)
    5 - Compressibility (1/bar, for volume)
    6 - Weight (%) of a component
    7 - Mode (Vol, Mol, or Wt proportion) of a phase
    8 - Composition (Mol or Wt) of a solution
    9 - Gruneisen thermal ratio
   10 - Adiabatic bulk modulus (bar)
   11 - Adiabatic shear modulus (bar)
   12 - Sound velocity (km/s)
   13 - P-wave velocity (Vp, km/s)
   14 - S-wave velocity (Vs, km/s)
   15 - Vp/Vs
   16 - Specific entropy (J/K/m3)
   17 - Entropy (J/K/kg)
   18 - Enthalpy (J/kg)
   19 - Heat Capacity (J/K/kg)
   20 - Specific mass of a phase (kg/m3-system)
   21 - Poisson ratio
   22 - Molar Volume (J/bar)
   23 - Dependent potentials (J/mol, bar, K)
   24 - Assemblage Index
   25 - Modes of all phases
   26 - Sound velocity T derivative (km/s/K)
   27 - P-wave velocity T derivative (km/s/K)
   28 - S-wave velocity T derivative (km/s/K)
   29 - Adiabatic bulk modulus T derivative (bar/K)
   30 - Shear modulus T derivative (bar/K)
   31 - Sound velocity P derivative (km/s/bar)
   32 - P-wave velocity P derivative (km/s/bar)
   33 - S-wave velocity P derivative (km/s/bar)
   34 - Adiabatic bulk modulus P derivative (unitless)
   35 - Shear modulus P derivative (unitless)
   37 - Absolute amount (Vol, Mol, or Wt) of a phase
   38 - Multiple property output (PHEMGP format)

13

Calculate individual phase properties (y/n)?

n

Include fluid in computation of aggregate (or modal) properties (y/n)?

n

 

Select properties [enter 0 to finish]:
    1 - Specific Enthalpy (J/m3)
    2 - Density (kg/m3)
    3 - Specific heat capacity (J/K/m3)
    4 - Expansivity (1/K, for volume)
    5 - Compressibility (1/bar, for volume)
    6 - Weight (%) of a component
    7 - Mode (Vol, Mol, or Wt proportion) of a phase
    8 - Composition (Mol or Wt) of a solution
    9 - Gruneisen thermal ratio
   10 - Adiabatic bulk modulus (bar)
   11 - Adiabatic shear modulus (bar)
   12 - Sound velocity (km/s)
   13 - P-wave velocity (Vp, km/s)
   14 - S-wave velocity (Vs, km/s)
   15 - Vp/Vs
   16 - Specific entropy (J/K/m3)
   17 - Entropy (J/K/kg)
   18 - Enthalpy (J/kg)
   19 - Heat Capacity (J/K/kg)
   20 - Specific mass of a phase (kg/m3-system)
   21 - Poisson ratio
   22 - Molar Volume (J/bar)
   23 - Dependent potentials (J/mol, bar, K)
   24 - Assemblage Index
   25 - Modes of all phases
   26 - Sound velocity T derivative (km/s/K)
   27 - P-wave velocity T derivative (km/s/K)
   28 - S-wave velocity T derivative (km/s/K)
   29 - Adiabatic bulk modulus T derivative (bar/K)
   30 - Shear modulus T derivative (bar/K)
   31 - Sound velocity P derivative (km/s/bar)
   32 - P-wave velocity P derivative (km/s/bar)
   33 - S-wave velocity P derivative (km/s/bar)
   34 - Adiabatic bulk modulus P derivative (unitless)
   35 - Shear modulus P derivative (unitless)
   37 - Absolute amount (Vol, Mol, or Wt) of a phase
   38 - Multiple property output (PHEMGP format)

0

Change default variable range (y/n)?

n

Enter number of nodes in the T(K) and P(bar) directions: (help)

100 100

 

**warning ver637** Immiscibility occurs in one or more phases
interpolation will be turned off at all affected nodes.
To overide this feature at the risk of computing inconsistent properties
set solvus_tolerance = 1 and rerun VERTEX

Data ranges excluding values equal to bad_number (NaN) specified in perplex_option.dat:

        vp,km/s
min    8.032687
max    8.939233

Output has been written to the 2d tab format file: in17a_1.tab

2d tab format files can be processed with:

     PSTABLE - a Perple_X plotting program
     PERPLE_X_PLOT - a MATLAB plotting script
     PYWERAMI - petrol.natur.cuni.cz/~ondro/pywerami:home
     spread-sheet programs, e.g., EXCEL

for details on tab format refer to:
     perplex.ethz.ch/perplex/faq/perple_x_tab_file_format.txt

Select operational mode:
    1 - properties at specified conditions
    2 - properties on a 2d grid
    3 - properties along a 1d path
    4 - as in 3, but input from file
    0 - EXIT

 

0

 

As indicated in the above output, the seismic wave-speed data has been written to a file named in17a_1.tab, the structure of this file is described at perple_x_tab_file_format.

 

Top


Graphical Representation

 

Grid data generated by WERAMI can be plotted with external programs such as PYWERAMI or MatLab, alternatively (the relatively crude) Perple_X program PSTABLE can generate an encapsulated PostScript plot. PSTABLE and MatLab are demonstrated here.

 

PSTABLE

 

c\jamie\perp666> pstable

Perple_X version 6.6.5.9, compiled 5/14/2011.
 

Reading computational options from: perplex_option.dat

 
Perple_X plot options are currently set as:

    Keyword:               Value:     Permitted values [default]:
    axis_label_scale       1.20       [1.2] (rel)
    bounding_box :
                               0      [0] x-min (pts)
                               0      [0] y-min (pts)
                             800      [800] x-length (pts)
                             800      [800] y-length (pts)
    field_fill             T          [T] F
    field_label            T          [T] F
    field_label_scale      0.75       [0.72] (rel)
    font                   Helvetica
    grid                   F          [F] T
    half_ticks             T          [T] F
    line_width             1.00       0-99 [1.] (pts)
    picture_transformation :
                            0.180     [0.18] x-scale (rel)
                            0.180     [0.18] y-scale (rel)
                             130.     [0.18] x-translation (pts)
                             220.     [0.18] y-translation (pts)
                             0.00     [0.0]  rotation (deg)
    plot_aspect_ratio      1.000      [1.0] x_axis_length/y_axis_length
    replicate_label        0.250      0->1 [0.025]
    splines                T          [T] F
    tenth_ticks            F          [F] T
    text_scale             1.000      [1.] (rel)

To change these options edit or create the plot option file
See: www.perplex.ethz.ch/perplex_plot_options.html

 

Enter the tab file name [without the .tab suffix]:

in17a_1

 

Unlike VERTEX, WERAMI and PSSECT, PSTABLE does not request the project name, but rather any Perple_X tab format file; the name of a tab format file consists of a root and the suffix .tab (file type indicator). WERAMI generates tab file names by taking the project root and appending a counter to distinguish different calculations.

 

Plot the ratio of two dependent variables (Y/N)?
n

Select the dependent variable to be contoured:
1 - T(K)
2 - P(bar)
3 - vp,km/s

3

 

Although the WERAMI calculation illustrated previously requested only P-wave velocities as a function of the independent variables (pressure and temperature) specified in BUILD, the tab file in17a_1 includes the independent variables as pseudo-dependent variables hence the foregoing prompts. Whether WERAMI outputs the independent variables of a calculation as pseudo-dependent variables is determined by the spreadsheet keyword in perplex_option.dat.

 

Contour the ratio of values in separate tab files (y/n)?

If you answer yes the data from the file just read will define the
numerator of the ratio and you will be prompted next for a file
containing the data for the denominator.
n

PostScript will be written to file: in17a_1.ps
 

The Encapsulated PostScript output file (in17a_1.ps) from PSTABLE is named by adding the suffix .ps to the root of input file.

 

Modify the default plot (y/n)?
y

 

To illustrate some of the options available in PSTABLE, the default plot is modified here.

 

Contour log10 of the dependent variable (y/n)?
n
Reset plot limits (y/n)?
n
Contoured variable range:   8.03269    ->   8.93923
Range excluding zero values:   8.03269    ->   8.93923

Modify default contour interval (y/n)?
y
Enter min, max and interval for contours:
8.05 8.95 0.05
Echo contour data to file contor.dat (Y/N)?
n
Modify default axes numbering (y/n)?
y
Enter the starting value and interval for major tick marks on
the X-axis ( current values are:  673.      160.    )
Enter the new values:
673 100
Enter the starting value and interval for major tick marks on
the Y-axis ( current values are: 0.300E+05 0.100E+05)
Enter the new values:
3e4 1e4

 

The output (Figure 3) from PSTABLE can be printed on any PostScript printer or viewed on a PostScript viewer such as Ghostview. Alternatively, the PostScript plot may be imported into graphical editors such as Adobe Illustrator or CorelDraw. PSCONTOR does not label the plotted contours, but it does use a heavy solid line for the minimum contour and a heavy broken line to indicate the maximum contour value. If this information is not sufficient to determine the contour values, WERAMI can be used in mode 1 to query the value of the property of interest at any point within the phase diagram section.

 

MatLab Plotting

 

A MatLab script (perple_x_plot in Perple_X_6.6.7_MatLab_scripts.zip) is provided with Perple_X to allow users to view tab file data (see Perple_X tab file format for format details) quickly within Matlab. To make a MatLab plot, the user must start MatLab, change directory to the directory containing the perple_x_plot.m script, and type perple_x_plot (without the .m suffix), the rest is self-explanatory. The script provides several plot style options, two of which are shown in Figure 3.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Figure 3: Three plot variations of the computed metabasalt P-wave velocities (km/s). The uppermost plot is a 3D surface plot generated in MATLAB by the script perple_x_plot. The middle plot is a manually labeled contour plot generated in MATLAB by the script perple_x_plot. The lowermost plot is the plot generated by PSTABLE.
 

Back to Perple_X Home Page