This page was generated from docs/Examples/Cpx_Cpx_Liq_Thermobarometry/Cpx_Liq_Thermobarometry.ipynb. Interactive online version: Binder badge.

Python Notebook Download

Clinopyroxene-only and Clinopyroxene-Liquid Thermobarometry.

You need to install Thermobar once on your machine, if you haven’t done this yet, uncomment the line below (remove the #)

[1]:

#!pip install Thermobar

First, load the necessary python things

[2]:

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import Thermobar as pt

Now, load the data

[3]:

out=pt.import_excel('Cpx_Liq_Example.xlsx', sheet_name="Sheet1")
my_input=out['my_input']
Liqs=out['Liqs']
Cpxs=out['Cpxs']

Inspect the data to check it loaded properly

[4]:

Liqs.head()

[4]:
SiO2_Liq TiO2_Liq Al2O3_Liq FeOt_Liq MnO_Liq MgO_Liq CaO_Liq Na2O_Liq K2O_Liq Cr2O3_Liq P2O5_Liq H2O_Liq Fe3Fet_Liq NiO_Liq CoO_Liq CO2_Liq Sample_ID_Liq
0 51.1 0.93 17.5 8.91 0.18 6.09 11.50 3.53 0.17 0.0 0.15 3.8 0.0 0.0 0.0 0.0 0
1 51.5 1.19 19.2 8.70 0.19 4.98 10.00 3.72 0.42 0.0 0.14 6.2 0.0 0.0 0.0 0.0 1
2 59.1 0.54 19.1 5.22 0.19 3.25 7.45 4.00 0.88 0.0 0.31 6.2 0.0 0.0 0.0 0.0 2
3 52.5 0.98 19.2 8.04 0.20 4.99 9.64 4.15 0.21 0.0 0.14 6.2 0.0 0.0 0.0 0.0 3
4 56.2 0.34 20.4 5.88 0.20 2.58 7.18 6.02 1.02 0.0 0.23 6.2 0.0 0.0 0.0 0.0 4 
[5]:

Cpxs.head()

[5]:
SiO2_Cpx TiO2_Cpx Al2O3_Cpx FeOt_Cpx MnO_Cpx MgO_Cpx CaO_Cpx Na2O_Cpx K2O_Cpx Cr2O3_Cpx Sample_ID_Cpx
0 51.5 0.50 3.70 5.18 0.09 15.8 22.8 0.24 0.0 0.66 0
1 50.3 0.73 4.12 5.83 0.00 15.0 22.7 0.24 0.0 0.28 1
2 47.3 1.75 7.85 6.51 0.14 13.1 22.5 0.25 0.0 0.22 2
3 51.1 0.63 4.41 5.66 0.13 15.6 22.6 0.23 0.0 0.27 3
4 51.0 0.56 4.14 7.33 0.20 14.4 22.4 0.31 0.0 0.09 4

Getting help

  • At any point, you can do help(pt.function) to get some more information

  • For example, here we get information on inputs/outputs for Cpx-Liq thermometry, including the equation options

[6]:

help(pt.calculate_cpx_liq_temp)

Help on function calculate_cpx_liq_temp in module Thermobar.clinopyroxene_thermobarometry:

calculate_cpx_liq_temp(*, equationT, cpx_comps=None, liq_comps=None, meltmatch=None, P=None, eq_tests=False, H2O_Liq=None, Fe3Fet_Liq=None, sigma=1, Kd_Err=0.03)
    Clinopyroxene-Liquid thermometry, calculates temperature in Kelvin
    (and equilibrium tests as an option)

    Parameters
    -------
    cpx_comps: pandas.DataFrame
        Clinopyroxene compositions with column headings SiO2_Cpx, MgO_Cpx etc.

    liq_comps: pandas.DataFrame
        Liquid compositions with column headings SiO2_Liq, MgO_Liq etc.

    Or:

    meltmatch: pandas.DataFrame
        Combined dataframe of cpx-Liquid compositions
        Used for calculate_cpx_liq_press_temp_matching function.

    EquationT: str
        Choice of equation:
        Cpx-Liquid
        |  T_Put1996_eqT1  (P-indep, H2O-indep)
        |  T_Mas2013_eqTalk1  (P-indep, H2O-indep, alk adaption of T1)
        |  T_Brug2019  (P-indep, H2O-indep)
        |  T_Put1996_eqT2 (P-dep, H2O-indep)
        |  T_Mas2013_eqTalk2  (P-dep, H2O-indep, alk adaption of T2)
        |  T_Put1999  (P-dep, H2O-indep)
        |  T_Put2003  (P-dep, H2O-indep)
        |  T_Put1999  (P-dep, H2O-indep)
        |  T_Put2008_eq33  (P-dep, H2O-dep)
        |  T_Mas2013_eqalk33  (P-dep, H2O-dep, alk adaption of eq33)
        |  T_Mas2013_Palk2012 (P-indep, H2O_dep)
        |  T_Petrelli2020_Cpx_Liq (gives voting)
        |  T_Jorgenson2022_Cpx_Liq (gives voting)
        |  T_Petrelli2020_Cpx_Liq_onnx (gives consistent result every time)


    P: float, int, pandas.Series, str  ("Solve")
        Pressure in kbar
        Only needed for P-sensitive thermometers.
        If enter P="Solve", returns a partial function
        Else, enter an integer, float, or panda series


    eq_tests: bool
        If False, just returns pressure (default) as a panda series
        If True, returns pressure, Values of Eq tests,
        Values of Eq tests (Kd, EnFs, DiHd, CaTs, CrCaTs),
        as well as user-entered cpx and liq comps and components.


    Returns
    -------
    If eq_tests=False
        pandas.Series: Temperature in Kelvin
    If eq_tests=True
        pandas.DataFrame: Temperature in Kelvin +
        Eq Tests + cpx+liq comps + components

Example 1 - Calculating Temperature

1a - Temperature for a known pressure and water content

  • Here, we calculate temperature using the H2O content given in the H2O_Liq column in the user-entered spreadsheet (the default), and P=5 kbar

  • There are a number of equations (see help above), but here we use T_Put2008_eq33 for temperature

[7]:

Temp_T33=pt.calculate_cpx_liq_temp(cpx_comps=Cpxs, liq_comps=Liqs,
                        equationT="T_Put2008_eq33", P=5)-273.15 # Default Kelvin
Temp_T33

C:\Users\penny\anaconda3\lib\site-packages\pandas\core\indexing.py:2120: FutureWarning: In a future version, the Index constructor will not infer numeric dtypes when passed object-dtype sequences (matching Series behavior)
  new_ix = Index(new_ix)

[7]:

0     1091.562867
1     1031.605285
2      999.303442
3     1032.536881
4      979.726983
5     1076.821062
6     1096.927022
7     1097.227100
8     1113.361468
9     1084.298534
10    1091.178760
11    1071.691587
12    1089.996302
13    1028.359107
14    1028.359107
15    1028.359107
16    1093.611526
17            NaN
18    1087.546061
19    1104.700975
dtype: float64

1b - Temperature, overwriting the spreadsheet water content in the function itself

  • Here, we are reseting water to 0 wt%. You can see the temperatures are much higher

[6]:

Temp_T33_0H2O=pt.calculate_cpx_liq_temp(cpx_comps=Cpxs, liq_comps=Liqs,
                                        equationT="T_Put2008_eq33",
                                        P=5, H2O_Liq=0)-273.15
Temp_T33_0H2O

[6]:

0    1142.969951
1    1109.923744
2    1073.681511
3    1110.970579
4    1051.769286
dtype: float64

1c- Lets use the thermometer of Brugman and Till, 2019.

  • This returns a number of warnings, because the authors recomend a compositional calibration range, which our entered cpx compositions lie outside of.

[7]:

Temp_TBrug=pt.calculate_cpx_liq_temp(cpx_comps=Cpxs, liq_comps=Liqs,
                                    equationT="T_Brug2019", P=5)-273.15
Temp_TBrug

Youve selected a P-independent function

C:\Users\penny\AppData\Local\Temp\ipykernel_4796\1899851194.py:1: UserWarning: Some inputted CPX compositions have Cpx Mg#>0.65;.
  Temp_TBrug=pt.calculate_cpx_liq_temp(cpx_comps=Cpxs, liq_comps=Liqs,
C:\Users\penny\AppData\Local\Temp\ipykernel_4796\1899851194.py:1: UserWarning: Some inputted CPX compositions have Al2O3>7 wt%;.
  Temp_TBrug=pt.calculate_cpx_liq_temp(cpx_comps=Cpxs, liq_comps=Liqs,
C:\Users\penny\AppData\Local\Temp\ipykernel_4796\1899851194.py:1: UserWarning: Some inputted Liq compositions have  SiO2<70 wt%;
  Temp_TBrug=pt.calculate_cpx_liq_temp(cpx_comps=Cpxs, liq_comps=Liqs,
c:\users\penny\onedrive - oregon state university\postdoc\pymme\mybarometers\thermobar_outer\src\Thermobar\clinopyroxene_thermobarometry.py:1165: UserWarning: which is outside the recomended calibration range of Brugman and Till (2019)
  w.warn("which is outside the recomended calibration range of Brugman and Till (2019)")

[7]:

0    1797.689417
1    1632.682392
2    1202.363348
3    1558.356437
4    1028.947694
dtype: float64

1d - We can also specify eq_tests=True to get a full dataframe back with all the components, as well as a number of equilibrium test values

  • You could then extract just the temps using Temp_T33_0H2O_EqTests[‘T_K_calc’] or any other column you want the same way

[9]:

Temp_T33_0H2O_EqTests=pt.calculate_cpx_liq_temp(cpx_comps=Cpxs, liq_comps=Liqs,
                                                equationT="T_Put2008_eq33",
                                                P=5, H2O_Liq=0, eq_tests=True)
Temp_T33_0H2O_EqTests

Using Fe3FeT from input file to calculate Kd Fe-Mg

[9]:
P_kbar_calc T_K_calc Eq Tests Neave2017? Delta_Kd_Put2008 Delta_Kd_Mas2013 Delta_EnFs_Mollo13 Delta_EnFs_Put1999 Delta_CaTs_Put1999 Delta_DiHd_Mollo13 Delta_DiHd_Put1999 ... Delta_EnFs_I_M_Mollo13 CaTs_Pred_Put1999 Delta_CaTs_I_M_Put1999 CrCaTs_Pred_Put1999 Delta_CrCaTs_I_M_Put1999 CaTi_Pred_Put1999 Delta_CaTi_I_M_Put1999 Jd_Pred_Put1999 Delta_Jd_Put1999 Delta_Jd_I_M_Put1999
0 5 1416.119951 False 0.042816 0.094027 0.010858 0.021183 0.016957 0.089424 0.001586 ... -0.010858 0.013418 -0.016957 0.0 0.009562 0.043320 0.002307 0.016184 0.000872 0.000872
1 5 1383.073744 False 0.036792 0.090778 0.007575 0.021866 0.023796 0.104285 0.024207 ... -0.007575 0.013021 -0.023796 0.0 0.004122 0.056932 0.011994 0.017429 0.000099 0.000099
2 5 1346.831511 False 0.058659 0.172858 0.043395 0.013382 0.075254 0.006413 0.103915 ... 0.043395 0.016630 -0.075254 0.0 0.003245 0.029418 0.042249 0.018331 0.000246 0.000246
3 5 1384.120579 False 0.034331 0.101144 0.012488 0.028698 0.031526 0.084045 0.024697 ... -0.012488 0.014301 -0.031526 0.0 0.003909 0.050962 0.009764 0.019332 0.002997 0.002997
4 5 1324.919286 False 0.022156 0.097447 0.021106 0.049019 0.032121 0.071950 0.035481 ... -0.021106 0.010990 -0.032121 0.0 0.001315 0.050713 0.014761 0.027953 0.005738 0.005738

5 rows × 130 columns

Example 2 - Calculating pressure for a known temperature

2a - Pressure at fixed temperature (T=1300 K), and pressures from Neave and Putirka (2017)

[11]:

P_FixedTNeave=pt.calculate_cpx_liq_press(cpx_comps=Cpxs, liq_comps=Liqs,
                                       equationP="P_Neave2017", T=1300)
P_FixedTNeave

[11]:

0    0.634602
1    1.655874
2    1.146083
3    1.028349
4    3.854147
dtype: float64

2b - Equation 30 from Putirka (2008), overwriting input water, return equilibrium tests

  • Here we change equation P, overwrite the H2O content in the function, and ask for equilibrium tests.

  • We are selecting equation 30 from Putirka (2008) this time, T=1300 K

[12]:

Temp_P30_0H2O=pt.calculate_cpx_liq_press(cpx_comps=Cpxs, liq_comps=Liqs,
                                         equationP="P_Put2008_eq30",
                                         T=1300, H2O_Liq=0, eq_tests=True)
Temp_P30_0H2O

Using Fe3FeT from input file to calculate Kd Fe-Mg

[12]:
P_kbar_calc T_K_calc Eq Tests Neave2017? Delta_Kd_Put2008 Delta_Kd_Mas2013 Delta_EnFs_Mollo13 Delta_EnFs_Put1999 Delta_CaTs_Put1999 Delta_DiHd_Mollo13 Delta_DiHd_Put1999 ... Delta_EnFs_I_M_Mollo13 CaTs_Pred_Put1999 Delta_CaTs_I_M_Put1999 CrCaTs_Pred_Put1999 Delta_CrCaTs_I_M_Put1999 CaTi_Pred_Put1999 Delta_CaTi_I_M_Put1999 Jd_Pred_Put1999 Delta_Jd_Put1999 Delta_Jd_I_M_Put1999
0 -2.090909 1300 False 0.015390 0.116829 0.002696 0.003912 0.019173 0.089851 0.350659 ... -0.002696 0.011202 -0.019173 0.0 0.009562 0.047180 0.006167 0.016165 0.000891 0.000891
1 -2.604072 1300 True 0.016996 0.108121 0.006605 0.005373 0.025431 0.049622 0.286210 ... 0.006605 0.011385 -0.025431 0.0 0.004122 0.060657 0.015718 0.017407 0.000077 0.000077
2 -2.145841 1300 False 0.069927 0.183575 0.074576 0.029940 0.076260 0.125978 0.249574 ... 0.074576 0.015625 -0.076260 0.0 0.003245 0.030292 0.041375 0.018308 0.000223 0.000223
3 -2.334004 1300 False 0.014291 0.117657 0.003134 0.008701 0.033245 0.066274 0.281808 ... 0.003134 0.012582 -0.033245 0.0 0.003909 0.054136 0.012939 0.019308 0.002973 0.002973
4 -0.177987 1300 False 0.016126 0.102892 0.007476 0.045035 0.032562 0.004383 0.113023 ... -0.007476 0.010549 -0.032562 0.0 0.001315 0.051697 0.015744 0.027927 0.005712 0.005712

5 rows × 130 columns

2c - As above, but setting a fixed Fe3Fet_Liq ratio

  • Can overwrite the Fe3Fet in the input spreadsheet to a different value, affects calculations of delta Kd as this uses just Fe2+ in the melt

  • Note, it is debated whether Kd Fe-Mg should be calculated with all Fe (to do that here, specify Fe3Fet_Liq=0, Putirka), or using just Fe2+ (e.g., Neave and Putirka, 2017)

  • you can compare the delta Kd Put 2008 in this option from the answers above. You can see, by adding 30% Fe3+, you have become further from equilibrium

[13]:

Temp_P30_0H2O_30Fe=pt.calculate_cpx_liq_press(cpx_comps=Cpxs,
                    liq_comps=Liqs, equationP="P_Put2008_eq30",
                    T=1300, H2O_Liq=0, Fe3Fet_Liq=0.3, eq_tests=True)
Temp_P30_0H2O_30Fe

[13]:
P_kbar_calc T_K_calc Eq Tests Neave2017? Delta_Kd_Put2008 Delta_Kd_Mas2013 Delta_EnFs_Mollo13 Delta_EnFs_Put1999 Delta_CaTs_Put1999 Delta_DiHd_Mollo13 Delta_DiHd_Put1999 ... Delta_EnFs_I_M_Mollo13 CaTs_Pred_Put1999 Delta_CaTs_I_M_Put1999 CrCaTs_Pred_Put1999 Delta_CrCaTs_I_M_Put1999 CaTi_Pred_Put1999 Delta_CaTi_I_M_Put1999 Jd_Pred_Put1999 Delta_Jd_Put1999 Delta_Jd_I_M_Put1999
0 -2.090909 1300 False 0.080646 0.212865 0.002696 0.003912 0.019173 0.089851 0.350659 ... -0.002696 0.011202 -0.019173 0.0 0.009562 0.047180 0.006167 0.016165 0.000891 0.000891
1 -2.604072 1300 False 0.078351 0.203468 0.006605 0.005373 0.025431 0.049622 0.286210 ... 0.006605 0.011385 -0.025431 0.0 0.004122 0.060657 0.015718 0.017407 0.000077 0.000077
2 -2.145841 1300 False 0.202528 0.316176 0.074576 0.029940 0.076260 0.125978 0.249574 ... 0.074576 0.015625 -0.076260 0.0 0.003245 0.030292 0.041375 0.018308 0.000223 0.000223
3 -2.334004 1300 False 0.082216 0.214165 0.003134 0.008701 0.033245 0.066274 0.281808 ... 0.003134 0.012582 -0.033245 0.0 0.003909 0.054136 0.012939 0.019308 0.002973 0.002973
4 -0.177987 1300 False 0.079595 0.198613 0.007476 0.045035 0.032562 0.004383 0.113023 ... -0.007476 0.010549 -0.032562 0.0 0.001315 0.051697 0.015744 0.027927 0.005712 0.005712

5 rows × 130 columns

Example 3 - Iterating pressure and temperature

  • In reality, unlesa you are an experimentalist, you rarely know one of pressure or temperature

  • In Keith Putirka’s spreadsheets, you can link up columns to iterate P and T towards a solution, this can be done here using the function calculate_cpx_liq_press_temp…

3a - Iterating equation 30 from Putirka (2008) for P, and equation 33 from Putirka (2008) for T

  • Without specifying anything else, you get a dataframe with columns for calculated pressure and temperature

[14]:

PT_iter_30_31=pt.calculate_cpx_liq_press_temp(cpx_comps=Cpxs, liq_comps=Liqs,
                                              equationP="P_Put2008_eq30",
                                              equationT="T_Put2008_eq33")
PT_iter_30_31

[14]:
P_kbar_calc T_K_calc
0 2.530914 1352.408784
1 1.786845 1290.151507
2 1.171520 1255.933868
3 2.143416 1292.669093
4 2.763538 1243.469600

3b - Same as above, but with eq_tests=True

  • Get all equilibrium tests, and input compostions as a larger dataframe.

[15]:

PT_iter_30_31_EqTests=pt.calculate_cpx_liq_press_temp(cpx_comps=Cpxs, liq_comps=Liqs, equationP="P_Put2008_eq30",
                                                      equationT="T_Put2008_eq33", eq_tests=True)
PT_iter_30_31_EqTests

Using Fe3FeT from input file to calculate Kd Fe-Mg

[15]:
P_kbar_calc T_K_calc Eq Tests Neave2017? Delta_Kd_Put2008 Delta_Kd_Mas2013 Delta_EnFs_Mollo13 Delta_EnFs_Put1999 Delta_CaTs_Put1999 Delta_DiHd_Mollo13 Delta_DiHd_Put1999 ... Delta_EnFs_I_M_Mollo13 CaTs_Pred_Put1999 Delta_CaTs_I_M_Put1999 CrCaTs_Pred_Put1999 Delta_CrCaTs_I_M_Put1999 CaTi_Pred_Put1999 Delta_CaTi_I_M_Put1999 Jd_Pred_Put1999 Delta_Jd_Put1999 Delta_Jd_I_M_Put1999
0 2.530914 1352.408784 True 0.027981 0.106599 0.011304 0.009468 0.018174 0.014780 0.160622 ... -0.011304 0.012201 -0.018174 0.0 0.009562 0.045314 0.004301 0.016178 0.000878 0.000878
1 1.786845 1290.151507 True 0.014591 0.110154 0.009917 0.002832 0.025619 0.009195 0.327542 ... -0.009917 0.011198 -0.025619 0.0 0.004122 0.061147 0.016209 0.017420 0.000090 0.000090
2 1.171520 1255.933868 False 0.080784 0.193537 0.046615 0.050574 0.077199 0.134845 0.434655 ... 0.046615 0.014686 -0.077199 0.0 0.003245 0.031201 0.040466 0.018319 0.000234 0.000234
3 2.143416 1292.669093 False 0.012502 0.119081 0.015465 0.006444 0.033388 0.022622 0.311452 ... -0.015465 0.012439 -0.033388 0.0 0.003909 0.054443 0.013245 0.019323 0.002988 0.002988
4 2.763538 1243.469600 False 0.002154 0.115101 0.025085 0.033393 0.033542 0.029538 0.345537 ... -0.025085 0.009568 -0.033542 0.0 0.001315 0.054153 0.018201 0.027943 0.005728 0.005728

5 rows × 130 columns

Example 3 - Cpx-only Barometry

  • Very similar to above, just don’t need liq_comps input

3a -Pressure only, using equation 32b (at T=1300 K), and H2O=0

  • This equation requires H2O content in the liquid. If you don’t enter anything, it assumes H2O=0

  • else specify using H2O_Liq=….

  • it prints a warning telling you that by defualt, this is what the function is doing

[16]:

eq32b_noH=pt.calculate_cpx_only_press(cpx_comps=Cpxs, T=1300,
       equationP="P_Put2008_eq32b")
eq32b_noH

c:\users\penny\onedrive - oregon state university\postdoc\pymme\mybarometers\thermobar_outer\src\Thermobar\clinopyroxene_thermobarometry.py:2051: UserWarning: This Cpx-only barometer is sensitive to H2O content of the liquid.  By default, this function uses H2O=0 wt%, else you can enter a value of H2O_Liq in the function
  w.warn('This Cpx-only barometer is sensitive to H2O content of the liquid. '

[16]:

0   -1.567381
1   -1.202221
2   -0.317946
3   -0.996608
4    0.467437
dtype: float64

3b - Pressure only, using 5 wt% water

[17]:

eq32b_5H=pt.calculate_cpx_only_press(cpx_comps=Cpxs, T=1300,
       equationP="P_Put2008_eq32b", H2O_Liq=5)
eq32b_5H

[17]:

0    0.697619
1    1.062779
2    1.947054
3    1.268392
4    2.732437
dtype: float64

3c - Temperature-only using eq 32d at 5 kbar

[19]:

eq32d_5kbar=pt.calculate_cpx_only_temp(cpx_comps=Cpxs, equationT="T_Put2008_eq32d",
         P=5)
eq32d_5kbar

[19]:

0    1457.849197
1    1441.107847
2    1415.899813
3    1455.722443
4    1441.512524
dtype: float64

3d - Iterating P from 32b, and T from 32d, with H2O=5

[20]:

eq32b_32d_5H=pt.calculate_cpx_only_press_temp(cpx_comps=Cpxs, equationT="T_Put2008_eq32d",
       equationP="P_Put2008_eq32b", H2O_Liq=5)
eq32b_32d_5H

[20]:
P_kbar_calc T_K_calc SiO2_Cpx TiO2_Cpx Al2O3_Cpx FeOt_Cpx MnO_Cpx MgO_Cpx CaO_Cpx Na2O_Cpx K2O_Cpx Cr2O3_Cpx Sample_ID_Cpx
0 3.889323 1448.656494 51.5 0.50 3.70 5.18 0.09 15.8 22.8 0.24 0 0.66 0
1 3.721484 1430.647515 50.3 0.73 4.12 5.83 0.00 15.0 22.7 0.24 0 0.28 1
2 3.982623 1407.721620 47.3 1.75 7.85 6.51 0.14 13.1 22.5 0.25 0 0.22 2
3 4.567865 1452.151021 51.1 0.63 4.41 5.66 0.13 15.6 22.6 0.23 0 0.27 3
4 5.941303 1449.216073 51.0 0.56 4.14 7.33 0.20 14.4 22.4 0.31 0 0.09 4

Example 4 - Plotting a Cpx-Liq Rhodes diagram to assess Fe-Mg equilibrium using fixed Kd values

  • The function calculate_cpx_rhodes_diagram_lines calculates the lines needed for the plot in a number of ways

  • There is disagrement in the literature as to whether Kd Fe-Mg should be assessed using just Fe2+ in the melt, or FeT, so we show both scenarios here.

Step 1 - Calculate Mg# for liq and cpxs

  • A number of functions in thermobar let you do it, we use this one here because it returns Mg#s for both phases

[22]:

cpx_comps_Fe3=pt.calculate_clinopyroxene_liquid_components(liq_comps=Liqs,
                cpx_comps=Cpxs, Fe3Fet_Liq=0.2)

  • Mgnos of Cpx are stored in the column accesed by cpx_comps_Fe3[‘Mgno_Cpx’]

  • Mgnos of Liq are stored in the column accesed by cpx_comps_Fe3[‘Mgno_Liq_noFe3’], for no Fe3, or cpx_comps_Fe3[‘Mgno_Liq_Fe2’] using just Fe2+ (e.g., 20% of Fe is 3+ here)

Step 2 - Calculate equilibrium lines to show on the rhodes diagram

  • You tell the function the min and max glass Mg# you want to show, e.g., the xspan of your plot

  • It returns lines for Kd=0.28+-0.08 after Putirka (2008)

[23]:

# Want to calculate Mg# to show on diagram between say 0.4 and 0.7 for the glass
eq_lines_1=pt.calculate_cpx_rhodes_diagram_lines(Min_Mgno=0.4, Max_Mgno=0.7)
eq_lines_1

[23]:
Mg#_Liq Eq_Cpx_Mg# (Kd=0.28) Eq_Cpx_Mg# (Kd=0.2) Eq_Cpx_Mg# (Kd=0.36)
0 0.400000 0.704225 0.769231 0.649351
1 0.403030 0.706845 0.771462 0.652217
2 0.406061 0.709445 0.773672 0.655065
3 0.409091 0.712025 0.775862 0.657895
4 0.412121 0.714586 0.778032 0.660707
... ... ... ... ...
95 0.687879 0.887273 0.916801 0.859588
96 0.690909 0.888681 0.917874 0.861287
97 0.693939 0.890081 0.918941 0.862979
98 0.696970 0.891473 0.920000 0.864662
99 0.700000 0.892857 0.921053 0.866337

100 rows × 4 columns

Step 3 - Combine these on a plot

  • You might need to adust the x and y limits of this plot for your own data

[25]:

fig, (ax1) = plt.subplots(1, 1, figsize = (6,5))
ax1.plot(eq_lines_1['Mg#_Liq'], eq_lines_1['Eq_Cpx_Mg# (Kd=0.2)'], ':r', label="K$_d$=0.2")
ax1.plot(eq_lines_1['Mg#_Liq'], eq_lines_1['Eq_Cpx_Mg# (Kd=0.36)'], ':r', label="K$_d$=0.36")
ax1.plot(eq_lines_1['Mg#_Liq'], eq_lines_1['Eq_Cpx_Mg# (Kd=0.28)'], '-r', label="K$_d$=0.28")

ax1.plot(cpx_comps_Fe3['Mgno_Liq_noFe3'], cpx_comps_Fe3['Mgno_Cpx'], '*k', mfc='yellow', ms=8, label="No Fe3")
ax1.plot(cpx_comps_Fe3['Mgno_Liq_Fe2'], cpx_comps_Fe3['Mgno_Cpx'], 'dk', mfc='green', ms=8, label="20% Fe3")
ax1.legend()
ax1.set_xlabel('Mg# Glass')
ax1.set_ylabel('Mg# Cpx')
# adjust x and y limits
ax1.set_xlim([0.4, 0.7])
ax1.set_ylim([0.6, 0.95])

[25]:

(0.6, 0.95)
../../_images/Examples_Cpx_Cpx_Liq_Thermobarometry_Cpx_Liq_Thermobarometry_49_1.png

Example 5 - Rhodes diagram using equation 35 of Putirka to calculate Kd as a function of T.

  • Here, we plot the equilibrium fields on the Rhodes diagram using equation 35 of Putirka, which is T-sensitive

  • Must specify T in Kelvin. Then return column for default 0.28, as well as the results for Putirka eq 35.

[26]:

eq_lines_2=pt.calculate_cpx_rhodes_diagram_lines(Min_Mgno=0.4, Max_Mgno=0.7, T=1300)
eq_lines_2.head()

[26]:
Mg#_Liq Eq_Cpx_Mg# (Kd=0.28) Eq_Cpx_Mg# (Kd=0.2) Eq_Cpx_Mg# (Kd=0.36) Kd_Eq35_P2008 Eq_Cpx_Mg# (Kd from Eq 35 P2008) Eq_Cpx_Mg# (Eq 35 P2008)+0.08 Eq_Cpx_Mg# (Eq 35 P2008)-0.08
0 0.400000 0.704225 0.769231 0.649351 0.239475 0.735720 0.676036 0.806964
1 0.403030 0.706845 0.771462 0.652217 0.239475 0.738165 0.678791 0.808921
2 0.406061 0.709445 0.773672 0.655065 0.239475 0.740589 0.681528 0.810858
3 0.409091 0.712025 0.775862 0.657895 0.239475 0.742993 0.684246 0.812775
4 0.412121 0.714586 0.778032 0.660707 0.239475 0.745377 0.686945 0.814673 
[27]:

fig, (ax1) = plt.subplots(1, 1, figsize = (6,5))
ax1.plot(eq_lines_2['Mg#_Liq'], eq_lines_2['Eq_Cpx_Mg# (Kd from Eq 35 P2008)'], '-r', label="K$_d$=Put eq 35")
ax1.plot(eq_lines_2['Mg#_Liq'], eq_lines_2['Eq_Cpx_Mg# (Eq 35 P2008)+0.08'], ':r', label="K$_d$=Put eq 35 + 0.08")
ax1.plot(eq_lines_2['Mg#_Liq'], eq_lines_2['Eq_Cpx_Mg# (Eq 35 P2008)-0.08'], ':r', label="K$_d$=Put eq 35 - 0.08")
ax1.set_xlim([0.4, 0.7])
ax1.plot(cpx_comps_Fe3['Mgno_Liq_noFe3'], cpx_comps_Fe3['Mgno_Cpx'], '*k', mfc='yellow', ms=8, label="No Fe3")
ax1.plot(cpx_comps_Fe3['Mgno_Liq_Fe2'], cpx_comps_Fe3['Mgno_Cpx'], 'dk', mfc='green', ms=8, label="20% Fe3")
ax1.legend()
ax1.set_xlabel('Mg# Glass')
ax1.set_ylabel('Mg# Cpx')

# adjust x and y limits
ax1.set_xlim([0.4, 0.7])
ax1.set_ylim([0.6, 0.95])

[27]:

(0.6, 0.95)
../../_images/Examples_Cpx_Cpx_Liq_Thermobarometry_Cpx_Liq_Thermobarometry_52_1.png

Example 6 - You can also specify a minimum and maximum Kd value you wish to calculate Rhodes lines for (here 0.2, 0.3)

[28]:

eq_lines_3=pt.calculate_cpx_rhodes_diagram_lines(Min_Mgno=0.4, Max_Mgno=0.7, KdMin=0.2, KdMax=0.3)
eq_lines_3.head()

[28]:
Mg#_Liq Eq_Cpx_Mg# (Kd=0.28) Eq_Cpx_Mg# (Kd=0.2) Eq_Cpx_Mg# (Kd=0.36) Eq_Cpx_Mg# (KdMin=0.2) Eq_Cpx_Mg# (KdMax=0.3)
0 0.400000 0.704225 0.769231 0.649351 0.769231 0.689655
1 0.403030 0.706845 0.771462 0.652217 0.771462 0.692348
2 0.406061 0.709445 0.773672 0.655065 0.773672 0.695021
3 0.409091 0.712025 0.775862 0.657895 0.775862 0.697674
4 0.412121 0.714586 0.778032 0.660707 0.778032 0.700309 
[29]:

fig, (ax1) = plt.subplots(1, 1, figsize = (6,5))
ax1.plot(eq_lines_3['Mg#_Liq'], eq_lines_3['Eq_Cpx_Mg# (KdMin=0.2)'], ':r', label="K$_d$=Put eq 35 + 0.08")
ax1.plot(eq_lines_3['Mg#_Liq'], eq_lines_3['Eq_Cpx_Mg# (KdMax=0.3)'], ':r', label="K$_d$=Put eq 35 - 0.08")
ax1.set_xlim([0.4, 0.7])
ax1.plot(cpx_comps_Fe3['Mgno_Liq_noFe3'], cpx_comps_Fe3['Mgno_Cpx'], '*k', mfc='yellow', ms=8, label="No Fe3")
ax1.plot(cpx_comps_Fe3['Mgno_Liq_Fe2'], cpx_comps_Fe3['Mgno_Cpx'], 'dk', mfc='green', ms=8, label="20% Fe3")
ax1.legend()
ax1.set_xlabel('Mg# Glass')
ax1.set_ylabel('Mg# Cpx')

# adjust x and y limits
ax1.set_xlim([0.4, 0.7])
ax1.set_ylim([0.6, 0.95])

[29]:

(0.6, 0.95)
../../_images/Examples_Cpx_Cpx_Liq_Thermobarometry_Cpx_Liq_Thermobarometry_55_1.png

Example 7 - Can get all options by specifying a temp, and a min and max Kd

  • Can then plot them however you want.

[30]:

eq_lines_4=pt.calculate_cpx_rhodes_diagram_lines(Min_Mgno=0.4, Max_Mgno=0.7, T=1300, KdMin=0.2, KdMax=0.3)
eq_lines_4.head()

[30]:
Mg#_Liq Eq_Cpx_Mg# (Kd=0.28) Eq_Cpx_Mg# (Kd=0.2) Eq_Cpx_Mg# (Kd=0.36) Kd_Eq35_P2008 Eq_Cpx_Mg# (Kd from Eq 35 P2008) Eq_Cpx_Mg# (Eq 35 P2008)+0.08 Eq_Cpx_Mg# (Eq 35 P2008)-0.08 Eq_Cpx_Mg# (KdMin=0.2) Eq_Cpx_Mg# (KdMax=0.3)
0 0.400000 0.704225 0.769231 0.649351 0.239475 0.735720 0.676036 0.806964 0.769231 0.689655
1 0.403030 0.706845 0.771462 0.652217 0.239475 0.738165 0.678791 0.808921 0.771462 0.692348
2 0.406061 0.709445 0.773672 0.655065 0.239475 0.740589 0.681528 0.810858 0.773672 0.695021
3 0.409091 0.712025 0.775862 0.657895 0.239475 0.742993 0.684246 0.812775 0.775862 0.697674
4 0.412121 0.714586 0.778032 0.660707 0.239475 0.745377 0.686945 0.814673 0.778032 0.700309 
[ ]: