IN PROGRESS: Calculator: Solid Phase Properties

This notebook calculates the activity of endmembers given the composition of a solid phase. It replaces the functionality previously implemented by the MELTS Supplemental Calculator (https://melts.ofm-research.org/CalcForms/)

This calculator is a work in progress – please let us know if you find any errors!

Open this code in an executable MyBinder instance (MyBinder links may be slow to load– please be patient!):

Information

The following solution phases are implemented in this calculator:

• ‘Ol’ (olivine)

• ‘Fsp’ (feldspar)

• ‘Opx’ (orthopyroxene)

• ‘Cpx’ (clinopyroxene)

• ‘SplS’ (spinel solid solution)

• ‘Gt’ (garnet)

This notebook is experimental and has not yet been extensively tested! Proceed with caution!

User Input

1. Define phase composition in wt% oxides, temperature, and pressure below. Components not included in the phase model will be ignored.

2. Once you have input your data, click “Run” –> “Run All Cells” at the top of the notebook. This will run the entire notebook and produce an output table at the bottom of the file.

phase = 'SplS'

oxide_wt_dict = {
        'SiO2':     0,
        'TiO2':     0.118399,
        'Al2O3':    31.4741,
        'Fe2O3':    3.623535369,
        'Cr2O3':    33.6659,
        'FeO':      15.77978837,
        'MnO':      0.248043,
        'MgO':      13.0798,
        'NiO':      0.140058,
        'CoO':      0.000,
        'CaO':      0.019049,
        'Na2O':     0.0039,
        'K2O':      0.000,
        'P2O5':     0.000,
 }

TK = 794.12 + 273.15
Pbar = 6000

# Advanced options
database='MELTS_v1_2'

Calculations

The user does not need to modify this code– it will run automatically when “Run All Cells” is triggered.

###-----------------------------------------------------------------------------------------
### 0. Import modules
###-----------------------------------------------------------------------------------------

import numpy as np
import pandas as pd
from thermoengine import model, rockychem

###-----------------------------------------------------------------------------------------
### 1. Create phase instance
###-----------------------------------------------------------------------------------------

phs_obj = model.Database(database_name="MELTS_v1_2").get_phase(phase)
phs_endmems = rockychem.convert_phase_comp(phase, database=database, wt_oxides=oxide_wt_dict, to='mol_endmems')

if phase == 'SplS':
    # Custom method for calculating endmembers; more accurate for magnetite activity than default endmember calculation
    mol_oxides = rockychem.convert_bulk_comp(wt_oxides=oxide_wt_dict, to='mol_oxides')

    mol_cat = pd.Series([mol_oxides['TiO2'],
                        2*mol_oxides['Al2O3'],
                        2*mol_oxides['Cr2O3'],
                        2*mol_oxides['Fe2O3'],
                        mol_oxides['FeO'],
                        mol_oxides['MgO']],
                        ['Ti', 'Al', 'Cr', 'Fe3', 'Fe2', 'Mg'])

    cat_O = pd.Series([2, 3/2, 3/2, 3/2, 1, 1],
        ['Ti', 'Al', 'Cr', 'Fe3', 'Fe2', 'Mg'])
    mol_O = sum(mol_cat*cat_O)
    norm_cat = 4*mol_cat/mol_O

    magnesioaluminate = norm_cat['Mg']
    chromite = norm_cat['Cr']/2
    ulvospinel = norm_cat['Ti']
    magnetite = norm_cat['Fe3']/2
    hercynite = norm_cat['Fe2'] - (chromite + 2*ulvospinel + magnetite)

    phs_endmems = np.array([hercynite, magnesioaluminate, chromite, ulvospinel, magnetite])

###-----------------------------------------------------------------------------------------
### 2. Calculate Gibbs Free Energies of Formation (Composition-Independent)
###-----------------------------------------------------------------------------------------

end_array = np.identity(len(phs_endmems))
G_array = np.zeros_like(phs_endmems)
for i in range(len(phs_endmems)):
    G = phs_obj.gibbs_energy(TK,Pbar,mol=end_array[i,:])
    G_array[i] = G

###-----------------------------------------------------------------------------------------
### 3. Calculate Activities (Composition-Dependent)
###-----------------------------------------------------------------------------------------

activities = phs_obj.activity(TK, Pbar, mol=phs_endmems)

###-----------------------------------------------------------------------------------------
### 4. Combine and Display Data
###-----------------------------------------------------------------------------------------

data = {'Endmember': phs_obj.endmember_names,
        'Formula': phs_obj.endmember_formulas,
        'Gibbs Energy (kJ)': G_array/1000,
        'Activity': activities}
endmember_df = pd.DataFrame(data)
pd.set_option('display.max_columns', None)
pd.set_option('display.max_rows', None)
endmember_df

	Endmember	Formula	Gibbs Energy (kJ)	Activity
0	Hercynite	FeAl2O4	-2136.592698	0.203389
1	Magnesioaluminate	MgAl2O4	-2452.327826	0.480026
2	Chromite	FeCr2O4	-1664.092435	0.346179
3	Ulvospinel	Fe2TiO4	-1755.727916	0.017019
4	Magnetite	Fe3O4	-1358.466940	0.020706