Note
Go to the end to download the full example code.
MagmaForge: Fractional Decompression#
Demonstrates how to run a fractional (isentropic) decompression routine using MagmaForge.
Open this code in an executable MyBinder instance (MyBinder links may be slow to load– please be patient!):
Initialization#
Import necessary packages:
from thermoengine import magmaforge
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
Define a bulk composition:
DMM_comp_adj = {
'SiO2': 44.7354,
'TiO2': 0.130074,
'Al2O3': 3.98226,
'Fe2O3': 0.275456,
'Cr2O3': 0.570323,
'FeO': 8.01254,
'MgO': 38.752,
'NiO': 0.240136,
'CaO': 3.1718,
'Na2O': 0.130074
}
P_GPa = 1.45
Find the solidus#
Currently, ThermoEngine cannot equilibrate subsolidus, as a different equilibration routine is required when the liquid is not present. We plan to get subsolidus equilibration implemented soon.
sys = magmaforge.System(comp=DMM_comp_adj,
T_liquidus=True,
P_GPa=P_GPa,
database='MELTS_pMELTS',
)
sys.crystallize(method='equil',
T_step=5,
fix_fO2=False,
calc_args={'debug':0})
solidus_T = sys.T_C
print(f"Solidus found at {solidus_T :.1f}˚C")
Solidus found at 1396.1˚C
System and Calculations#
Define the system using magmaforge:
sys = magmaforge.System(comp=DMM_comp_adj,
T_C = solidus_T,
P_GPa=P_GPa,
database='MELTS_pMELTS'
)
Melt via isentropic decompression:
For fractional melting, the fraction retained is treated as a fixed porosity. By default, this value is set to 0.01.
Fractionation happens at the beginning of each step, prior to equilibration, so that the final state is an equilibrium state.
sys.melt_by_decompression(P_final_GPa=0.25,
method='frac',
P_step_GPa=0.01,
fix_fO2=False,
)
<thermoengine.magmaforge.system.System object at 0xffff9ae72fb0>
Output#
pd.set_option('display.max_columns', None)
pd.set_option('display.max_row', None)
phase_masses = sys.history.get_phase_mass_table(phases='present', index='P_GPa')
phase_masses['Aggregated Liquid'] = sys.history.get_fractionated_masses()
print(phase_masses)
Olivine Clinopyroxene Orthopyroxene Spinel Liquid \
1.45 54.761415 17.398883 23.437538 3.640274 0.761581
1.44 54.755699 17.292687 23.464876 3.633656 0.852773
1.43 54.750861 17.183296 23.492285 3.626705 0.913349
1.42 54.746608 17.071835 23.519775 3.619541 0.917464
1.41 54.742957 16.958262 23.547337 3.612157 0.921635
1.40 54.739930 16.842514 23.574968 3.604547 0.925870
1.39 54.737549 16.724533 23.602662 3.596704 0.930170
1.38 54.735833 16.604263 23.630408 3.588621 0.934533
1.37 54.734807 16.481643 23.658198 3.580290 0.938964
1.36 54.734495 16.356618 23.686015 3.571704 0.943461
1.35 54.734922 16.229134 23.713843 3.562854 0.948025
1.34 54.736113 16.099117 23.741669 3.553732 0.952667
1.33 54.738097 15.966526 23.769463 3.544329 0.957373
1.32 54.740902 15.831300 23.797201 3.534638 0.962154
1.31 54.744557 15.693381 23.824853 3.524648 0.967012
1.30 54.749094 15.552719 23.852384 3.514350 0.971943
1.29 54.754545 15.409239 23.879765 3.503734 0.976963
1.28 54.760945 15.262941 23.906928 3.492791 0.982035
1.27 54.768328 15.113724 23.933848 3.481509 0.987215
1.26 54.776732 14.961552 23.960466 3.469878 0.992471
1.25 54.786194 14.806402 23.986712 3.457887 0.997796
1.24 54.796757 14.648188 24.012539 3.445523 1.003233
1.23 54.808460 14.486902 24.037861 3.432776 1.008735
1.22 54.821349 14.322498 24.062603 3.419633 1.014331
1.21 54.835468 14.154935 24.086684 3.406081 1.020018
1.20 54.850867 13.984174 24.110014 3.392108 1.025800
1.19 54.867592 13.810201 24.132491 3.377701 1.031663
1.18 54.885696 13.632980 24.154013 3.362845 1.037625
1.17 54.905234 13.452461 24.174476 3.347526 1.043697
1.16 54.926261 13.268651 24.193755 3.331730 1.049847
1.15 54.948834 13.081517 24.211728 3.315444 1.056101
1.14 54.973014 12.891043 24.228264 3.298651 1.062451
1.13 54.998864 12.697211 24.243227 3.281336 1.068903
1.12 55.026447 12.500010 24.256474 3.263485 1.075452
1.11 55.055832 12.299428 24.267854 3.245081 1.082100
1.10 55.087088 12.095451 24.277216 3.226108 1.088852
1.09 55.120285 11.888087 24.284395 3.206552 1.095691
1.08 55.155499 11.677314 24.289230 3.186394 1.102642
1.07 55.192804 11.463145 24.291549 3.165621 1.109676
1.06 55.232279 11.245569 24.291181 3.144216 1.116813
1.05 55.274006 11.024589 24.287948 3.122163 1.124042
1.04 55.318062 10.800226 24.281671 3.099447 1.131344
1.03 55.364533 10.572473 24.272167 3.076054 1.138743
1.02 55.413505 10.341341 24.259255 3.051968 1.146220
1.01 55.465061 10.106853 24.242749 3.027176 1.153761
1.00 55.519286 9.869032 24.222466 3.001667 1.161362
0.99 55.576268 9.627896 24.198223 2.975428 1.169025
0.98 55.636087 9.383495 24.169841 2.948450 1.176711
0.97 55.698833 9.135844 24.137142 2.920724 1.184450
0.96 55.764582 8.885012 24.099954 2.892243 1.192184
0.95 55.833414 8.631050 24.058110 2.863004 1.199921
0.94 55.905403 8.374032 24.011454 2.833005 1.207632
0.93 55.980618 8.114035 23.959835 2.802246 1.215305
0.92 56.059121 7.851162 23.903117 2.770732 1.222903
0.91 56.140968 7.585517 23.841173 2.738471 1.230415
0.90 56.226205 7.317232 23.773896 2.705473 1.237801
0.89 56.314869 7.046449 23.701192 2.671753 1.245042
0.88 56.406982 6.773339 23.622992 2.637333 1.252089
0.87 56.502558 6.498074 23.539241 2.602234 1.258933
0.86 56.601593 6.220861 23.449911 2.566486 1.265524
0.85 56.704070 5.941921 23.355002 2.530123 1.271828
0.84 56.809956 5.661492 23.254536 2.493183 1.277817
0.83 56.919197 5.379834 23.148565 2.455707 1.283452
0.82 57.031726 5.097224 23.037174 2.417745 1.288699
0.81 57.147459 4.813943 22.920465 2.379345 1.293541
0.80 57.266285 4.530308 22.798589 2.340566 1.297922
0.79 57.388087 4.246625 22.671709 2.301465 1.301846
0.78 57.512722 3.963220 22.540028 2.262105 1.305272
0.77 57.640035 3.680423 22.403769 2.222551 1.308190
0.76 57.769856 3.398560 22.263183 2.182868 1.310589
0.75 57.901997 3.117968 22.118549 2.143127 1.312450
0.74 58.036262 2.838970 21.970160 2.103395 1.313782
0.73 58.172441 2.561890 21.818333 2.063742 1.314571
0.72 58.310315 2.287042 21.663402 2.024237 1.314827
0.71 58.449657 2.014727 21.505709 1.984948 1.314561
0.70 58.590236 1.745239 21.345614 1.945944 1.313774
0.69 58.731813 1.478854 21.183480 1.907288 1.312486
0.68 58.874147 1.215834 21.019678 1.869046 1.310714
0.67 59.016995 0.956428 20.854584 1.831278 1.308470
0.66 59.160115 0.700865 20.688571 1.794044 1.305779
0.65 59.303265 0.449360 20.522014 1.757399 1.302659
0.64 59.446203 0.202110 20.355285 1.721397 1.299135
0.63 59.576600 0.000000 20.186206 1.688763 1.266484
0.62 59.646275 0.000000 20.007297 1.669902 1.121912
0.61 59.716575 0.000000 19.827859 1.651255 1.121424
0.60 59.787530 0.000000 19.647808 1.632813 1.120979
0.59 59.859168 0.000000 19.467068 1.614568 1.120573
0.58 59.931528 0.000000 19.285543 1.596512 1.120219
0.57 60.004620 0.000000 19.103207 1.578641 1.119867
0.56 60.078491 0.000000 18.919936 1.560945 1.119586
0.55 60.153163 0.000000 18.735680 1.543420 1.119325
0.54 60.228663 0.000000 18.550372 1.526060 1.119094
0.53 60.305021 0.000000 18.363938 1.508859 1.118900
0.52 60.382259 0.000000 18.176329 1.491815 1.118725
0.51 60.460397 0.000000 17.987496 1.474925 1.118568
0.50 60.539472 0.000000 17.797350 1.458182 1.118459
0.49 60.619483 0.000000 17.605900 1.441589 1.118330
0.48 60.700459 0.000000 17.413079 1.425143 1.118233
0.47 60.782411 0.000000 17.218862 1.408843 1.118143
0.46 60.865350 0.000000 17.023228 1.392689 1.118053
0.45 60.949281 0.000000 16.826166 1.376683 1.117962
0.44 61.034213 0.000000 16.627659 1.360826 1.117874
0.43 61.120134 0.000000 16.427737 1.345123 1.117756
0.42 61.207045 0.000000 16.226402 1.329575 1.117632
0.41 61.294931 0.000000 16.023687 1.314187 1.117480
0.40 61.383775 0.000000 15.819636 1.298965 1.117296
0.39 61.473549 0.000000 15.614310 1.283916 1.117069
0.38 61.564232 0.000000 15.407759 1.269044 1.116811
0.37 61.655779 0.000000 15.200085 1.254359 1.116487
0.36 61.748145 0.000000 14.991386 1.239869 1.116105
0.35 61.841291 0.000000 14.781747 1.225581 1.115673
0.34 61.935149 0.000000 14.571314 1.211505 1.115155
0.33 62.029662 0.000000 14.360217 1.197651 1.114565
0.32 62.124760 0.000000 14.148597 1.184027 1.113897
0.31 62.220367 0.000000 13.936617 1.170643 1.113140
0.30 62.316412 0.000000 13.724427 1.157508 1.112305
0.29 62.412807 0.000000 13.512210 1.144632 1.111374
0.28 62.509467 0.000000 13.300139 1.132023 1.110354
0.27 62.606305 0.000000 13.088391 1.119688 1.109246
0.26 62.703241 0.000000 12.877128 1.107635 1.108061
0.25 62.800179 0.000000 12.666540 1.095870 1.106783
Aggregated Liquid
1.45 0.000000
1.44 0.000000
1.43 0.033196
1.42 0.124468
1.41 0.217343
1.40 0.311861
1.39 0.408074
1.38 0.506034
1.37 0.605791
1.36 0.707400
1.35 0.810916
1.34 0.916396
1.33 1.023906
1.32 1.133500
1.31 1.245245
1.30 1.359206
1.29 1.475450
1.28 1.594057
1.27 1.715074
1.26 1.838600
1.25 1.964707
1.24 2.093460
1.23 2.224965
1.22 2.359286
1.21 2.496514
1.20 2.636738
1.19 2.780054
1.18 2.926543
1.17 3.076308
1.16 3.229460
1.15 3.386081
1.14 3.546281
1.13 3.710164
1.12 3.877839
1.11 4.049411
1.10 4.224993
1.09 4.404699
1.08 4.588630
1.07 4.776914
1.06 4.969651
1.05 5.166963
1.04 5.368961
1.03 5.575742
1.02 5.787425
1.01 6.004114
1.00 6.225902
0.99 6.452877
0.98 6.685133
0.97 6.922726
0.96 7.165744
0.95 7.414220
0.94 7.668195
0.93 7.927683
0.92 8.192688
0.91 8.463180
0.90 8.739117
0.89 9.020421
0.88 9.306992
0.87 9.598688
0.86 9.895354
0.85 10.196786
0.84 10.502749
0.83 10.812977
0.82 11.127168
0.81 11.444983
0.80 11.766067
0.79 12.090006
0.78 12.416393
0.77 12.744773
0.76 13.074685
0.75 13.405652
0.74 13.737177
0.73 14.068770
0.72 14.399926
0.71 14.730147
0.70 15.058945
0.69 15.385831
0.68 15.710335
0.67 16.032000
0.66 16.350383
0.65 16.665060
0.64 16.975630
0.63 17.281708
0.62 17.554376
0.61 17.682650
0.60 17.810633
0.59 17.938385
0.58 18.065961
0.57 18.193429
0.56 18.320806
0.55 18.448178
0.54 18.575577
0.53 18.703048
0.52 18.830640
0.51 18.958382
0.50 19.086305
0.49 19.214467
0.48 19.342855
0.47 19.471512
0.46 19.600450
0.45 19.729678
0.44 19.859200
0.43 19.989022
0.42 20.119120
0.41 20.249489
0.40 20.380102
0.39 20.510930
0.38 20.641929
0.37 20.773066
0.36 20.904272
0.35 21.035486
0.34 21.166654
0.33 21.297684
0.32 21.428499
0.31 21.559012
0.30 21.689127
0.29 21.818758
0.28 21.947799
0.27 22.076152
0.26 22.203719
0.25 22.330412
Plot System Evolution#
Automatic plotting functions for decompression melting coming soon! The dotted lines show the instantaneous fractional melt compositions.
oxides_to_plot = ['MgO', 'SiO2', 'Al2O3', 'FeO', 'Fe2O3', 'CaO']
pressures = sys.history.get_pressures('GPa')
temperatures = sys.history.get_temperatures('C')
phase_masses = sys.history.get_phase_mass_table(phases='present', index='P_GPa')
mass_agg_melt = sys.history.get_fractionated_masses('aggregated')
phase_masses['Agg Liquid'] = mass_agg_melt
agg_liq_comps = sys.history.get_fractionated_comps(unit='wt_oxides', method='aggregated', norm=True)
instant_liq_comps = sys.history.get_fractionated_comps(unit='wt_oxides', method="incremental", norm=True)
phase_colors = {'Spinel': 'black',
'Garnet': 'purple',
'Olivine': 'yellowgreen',
'Clinopyroxene' : 'green',
'Orthopyroxene': 'darkolivegreen',}
# Create plot
fig, ax = plt.subplots(1, 3, sharey='row', gridspec_kw={'wspace': 0.1, 'width_ratios': [1,1,2]}, figsize=(8,4))
# Plot phase masses
for ph in phase_masses.columns:
if (ph != 'Liquid') and (ph != 'Agg Liquid'):
ax[0].plot(phase_masses[ph], pressures,
color=phase_colors[ph],label=ph)
ax[0].plot(phase_masses['Agg Liquid'], pressures,
lw=2, c='red', label='Agg. Melt')
# Plot temperature
ax[1].plot(temperatures, pressures, c='k', lw=2)
# Plot liquid composition
i=0
for ox in oxides_to_plot:
ax[2].plot(agg_liq_comps[ox], pressures, label=ox, c=f"C{i}")
ax[2].plot(instant_liq_comps[ox], pressures, c=f"C{i}", ls=':')
i += 1
# Add instantaneous melts to legend
ax[2].set_xlim(ax[2].get_xlim())
ax[2].set_ylim(ax[2].get_ylim())
ax[2].plot([100]*2, [100]*2, label=" ", lw=0)
ax[2].plot([100]*2, [100]*2, c='k', ls=':', label="Instantaneous\nMelts")
# Format axes
ax[0].invert_yaxis()
ax[0].legend(ncol=2, loc='upper left', bbox_to_anchor=(-0.4, 1.3))
ax[2].legend(ncol=3, loc='upper left', bbox_to_anchor=(-0.3, 1.3))
plt.subplots_adjust(top=0.8)
ax[0].set_ylabel('Pressure (GPa)')
ax[0].set_xlabel('Phase Mass')
ax[1].set_xlabel('Temperature (˚C)')
ax[2].set_xlabel('Aggregated Liquid Comp (wt%)')
plt.show()
Save Output#
All of the output can be saved to csv files for future use, for example, uncomment the following:
# phase_masses.to_csv('pMELTS_phase_masses.csv')
# agg_liq_comps.to_csv('pMELTS_aggregated_liquid_comps.csv')
Total running time of the script: (0 minutes 40.138 seconds)