This page was generated from
docs/Examples/EOS_calculations/Example5f_H2O_CO2_EOS.ipynb.
Interactive online version:
.
H2O-CO2 EOS calculations
DiadFit contains the EOS of Duan and Zhang (2006)
This code was converted from the C code from Yoshimura et al. (2023) - we thank them for making such a clear description of the mistakes in the original Duan and Zhang paper! https://doi.org/10.2465/jmps.221224a
There are two main types of calculations you might want to do with a mixed EOS.
Example 1) If you work on fluid inclusions, its nice to assume they are a pure CO\(_2\) fluid. But in reality, you probably have some H\(_2\)O in the fluid phase. You can calculate a pressure for a given XH2O value, assuming the H has been lost .
Example 2) If you are an experimentalist, you probably know P and T, and you have some idea of XH2O. You might want to calculate molar volumes, densities, compressabilities, activities and fugacities. This is also useful for people doing EOS calculations for melt inclusion vapour bubbles etc.
The two base functions in DiadFit are:
pf.calculate_entrapment_P_XH2O - Calculates pressure for a known density, temperature and XH2O
pf.calc_prop_knownP_EOS_DZ2006 - Calculates properties (density, fugacity, compressability, molar volumes) for a H2O-CO2 mixed fluid for a known pressure, temperature and XH2O
[1]:
# If you havent already - install DiadFit.
#!pip install DiadFit --upgrade
If you havent installed CoolProp, you will need that package too If you use conda, in your command line run: conda install -c conda-forge coolprop
[2]:
import DiadFit as pf
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import math
pf.__version__
[2]:
'0.0.91'
1. Correcting fluid inclusions for H2O
The theory for this is layed out in Hansteen and Klugel (2008). In all but the most anhydrous, deep systems, the exsolved fluid phase will contain some H2O. It is generally assumed that this is lost from the fluid inclusion, so you are left with CO2 you have measured by Raman spectroscopy
First, you have to calculate the bulk density of the original CO2+H2O mix. Then, using this density, you can use the mixed EOS to calculate pressure. Here, we compare this to pure CO2 EOS.
You can make 2 different end member assumptions - either all H is lost or the H exists as a film. The math is layed out in this notebook incase of interest!
1a: Correcting a single FI composition
Say you measured a CO2 density of 0.3, and you think XH2O=0.1 (perhaps based on melt inclusion data - you can get XH2O from VESical)
Lets compare the calculations
[3]:
Diff_EOS=pf.calculate_entrapment_P_XH2O(XH2O=0, CO2_dens_gcm3=0.3, T_K=1200+273.15, T_K_ambient=37+273.15 )
Diff_EOS
[3]:
| P_kbar_pureCO2_SW96 | P_kbar_pureCO2_SP94 | P_kbar_pureCO2_DZ06 | P_kbar_mixCO2_DZ06_Hloss | P_kbar_mixCO2_DZ06_no_Hloss | P Mix_Hloss/P Pure DZ06 | P Mix_no_Hloss/P Pure DZ06 | rho_mix_calc_Hloss | rho_mix_calc_noHloss | CO2_dens_gcm3 | T_K | XH2O | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 1.055698 | 1.070557 | 1.052613 | 1.052613 | 1.052613 | 1.0 | 1.0 | 0.3 | 0.3 | 0.3 | 1473.15 | 0 |
We can see that the mixed fluid pressure is ~13% higher than doing the pure calculation if H is lost, and 12% higher if H2O is in a film
1b - Visualizing changes across a wider space
You dont have to feed a single value into this function, you could do an entire spreadsheet of values.
Here, lets do calculations for different densities for a range of XH2O values
[5]:
XH2O=np.linspace(0.0001, 0.2, 50) # 50 evenly spaced XH2O values
df_01=pf.calculate_entrapment_P_XH2O(XH2O=XH2O, CO2_dens_gcm3=0.1, T_K=1200+273.15)
df_02=pf.calculate_entrapment_P_XH2O(XH2O=XH2O, CO2_dens_gcm3=0.2, T_K=1200+273.15)
df_03=pf.calculate_entrapment_P_XH2O(XH2O=XH2O, CO2_dens_gcm3=0.3, T_K=1200+273.15)
df_05=pf.calculate_entrapment_P_XH2O(XH2O=XH2O, CO2_dens_gcm3=0.5, T_K=1200+273.15)
df_07=pf.calculate_entrapment_P_XH2O(XH2O=XH2O, CO2_dens_gcm3=0.7, T_K=1200+273.15)
df_09=pf.calculate_entrapment_P_XH2O(XH2O=XH2O, CO2_dens_gcm3=0.9, T_K=1200+273.15)
df_11=pf.calculate_entrapment_P_XH2O(XH2O=XH2O, CO2_dens_gcm3=1.1, T_K=1200+273.15)
# Create a dictionary to hold your density values and corresponding dataframe columns directly
densities = {
0.1: df_01['P Mix_Hloss/P Pure DZ06'],
0.2: df_02['P Mix_Hloss/P Pure DZ06'],
0.3: df_03['P Mix_Hloss/P Pure DZ06'],
0.5: df_05['P Mix_Hloss/P Pure DZ06'],
0.7: df_07['P Mix_Hloss/P Pure DZ06'],
0.9: df_09['P Mix_Hloss/P Pure DZ06'],
1.1: df_11['P Mix_Hloss/P Pure DZ06'],
}
of course, you can write this in a loop if you want!
[6]:
# Define CO2 densities to loop over
CO2_densities = [0.1, 0.2, 0.3, 0.5, 0.7, 0.9, 1.1]
# Initialize a dictionary to hold the dataframes
dfs = {}
# Loop over CO2 densities
for density in CO2_densities:
# Calculate entrapment_P_XH2O for each CO2 density
df = pf.calculate_entrapment_P_XH2O(XH2O=XH2O, CO2_dens_gcm3=density, T_K=1200 + 273.15)
# Store the dataframe in the dictionary with a key indicating the CO2 density
densities[density] = df['P Mix_Hloss/P Pure DZ06']
Lets make a nice plot showing this
[7]:
import matplotlib.pyplot as plt
import matplotlib.colors as mcolors
# Create figure and axes
fig, ax1 = plt.subplots(1, 1, figsize=(5, 4))
# Choose a colormap
cmap = plt.cm.plasma_r
# Create a Normalize object for your densities
norm = mcolors.Normalize(vmin=min(densities.keys()), vmax=max(densities.keys()))
# Create a ScalarMappable object with the normalization and colormap
sm = plt.cm.ScalarMappable(cmap=cmap, norm=norm)
sm.set_array([])
# Plot each line on ax1 using the density as the key to retrieve the corresponding dataframe column
for density, data in densities.items():
print(density)
ax1.plot(XH2O, data, label=f'{density} g/cm$^3$', color=cmap(norm(density)))
ax1.legend()
ax1.set_xlabel('X$_{H_{2}O}^{molar}$')
ax1.set_ylabel('P Mix EOS/P Pure EOS')
#fig.colorbar(sm, ax=ax1, label='Density (g/cm$^3$)')
ax1.set_xticks([0, 0.05, 0.1, 0.15, 0.2])
ax1.set_xlim([0, 0.2])
ax1.set_ylim([1, 1.5])
plt.show()
fig.savefig('H2OCO2.png', dpi=200, bbox_inches='tight')
0.1
0.2
0.3
0.5
0.7
0.9
1.1
1c: Calculations on an entire spreadsheet
Lets say to start with you want to apply a correction for a constant XH2O to all fluid inclusion densities
[8]:
# Load in the excel spreadsheet
# Download it here - https://github.com/PennyWieser/DiadFit/blob/main/docs/Examples/EOS_calculations/FI_densities.xlsx
df=pd.read_excel('FI_densities.xlsx')
df.head()
[8]:
| Sample | Density_g_cm3 | Temp in C | |
|---|---|---|---|
| 0 | FI1 | 0.1 | 1100 |
| 1 | FI2 | 0.4 | 1200 |
| 2 | FI3 | 0.7 | 1150 |
| 3 | FI4 | 1.1 | 1200 |
| 4 | FI5 | 1.5 | 1300 |
[9]:
# Lets perform a calculation for 10 mol% H2O (A mol fraction of 0.1)
df_10perc=pf.calculate_entrapment_P_XH2O(XH2O=0.1, CO2_dens_gcm3=df['Density_g_cm3'], T_K=df['Temp in C'], T_K_ambient=37+273.15)
df_10perc
[9]:
| P_kbar_pureCO2_SW96 | P_kbar_pureCO2_SP94 | P_kbar_pureCO2_DZ06 | P_kbar_mixCO2_DZ06_Hloss | P_kbar_mixCO2_DZ06_no_Hloss | P Mix_Hloss/P Pure DZ06 | P Mix_no_Hloss/P Pure DZ06 | rho_mix_calc_Hloss | rho_mix_calc_noHloss | CO2_dens_gcm3 | T_K | XH2O | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 0.218732 | 0.217395 | 0.218748 | 0.242265 | 0.241101 | 1.107509 | 1.102188 | 0.104545 | 0.104069 | 0.1 | 1100 | 0.1 |
| 1 | 1.230950 | 1.249353 | 1.227782 | 1.388440 | 1.352069 | 1.130852 | 1.101229 | 0.418182 | 0.410668 | 0.4 | 1200 | 0.1 |
| 2 | 2.965314 | 3.003832 | 2.981447 | 3.464521 | 3.253963 | 1.162027 | 1.091404 | 0.731818 | 0.709115 | 0.7 | 1150 | 0.1 |
| 3 | 8.841732 | 8.826510 | 9.107578 | 11.291880 | 9.695014 | 1.239833 | 1.064500 | 1.150000 | 1.095127 | 1.1 | 1200 | 0.1 |
| 4 | 23.436505 | 24.924885 | 29.177486 | 37.972866 | 29.435091 | 1.301444 | 1.008829 | 1.568182 | 1.470608 | 1.5 | 1300 | 0.1 |
But this isnt realistic - we know that XH2O will vary as a function of pressure, getting higher for the shallower pressure inclusions.
if you have melt inclusions for your system, you could use saturation pressures calculated in VESical, which returns XH2O, to determine how XH2O varies as a function of pressure
Here, we load in melt inclusion data from the 2018 eruption of Kilauea, where XH2O was calculated using MagmaSat
[10]:
# Download excel file here! https://github.com/PennyWieser/DiadFit/blob/main/docs/Examples/EOS_calculations/FI_densities.xlsx
df_MI=pd.read_excel('FI_densities.xlsx', sheet_name='MI_data')
df_MI.head()
[10]:
| Sat P (Bars) | XCO2 (molar) | |
|---|---|---|
| 0 | 310 | 0.959107 |
| 1 | 360 | 0.963423 |
| 2 | 430 | 0.968782 |
| 3 | 970 | 0.986911 |
| 4 | 600 | 0.978391 |
Lets see how XH2O varies with pressure
[11]:
x=df_MI['Sat P (Bars)']/10
y=1-df_MI['XCO2 (molar)']
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(10,5))
ax1.plot(x, y, '.r', label='MI')
# we only see FI with pressures of 20-100 MPa ish,
filt=(x>10)& (x<100)
Pf = np.poly1d(np.polyfit(x[filt], y[filt],
3))
Px = np.linspace(10, 100, 101)
ax2.plot(x[filt], y[filt], '.r')
Py = Pf(Px)
ax2.plot(Px, Py, '-r')
ax2.set_ylim([0, 0.2])
ax2.set_xlim([10, 105])
ax1.set_ylabel('XH2O')
ax1.set_xlabel('MI P (MPa)')
ax2.set_ylabel('XH2O')
ax2.set_xlabel('MI P (MPa)')
ax1.legend()
[11]:
<matplotlib.legend.Legend at 0x2b9566c4310>
Now we need to loop - so initially use the polynomial to guess a XH2O for each MI, then re-evaluate after we get a new pressure for the FI accounting for H2O
[12]:
# Load in FI data of DeVitre and Wieser (2024 - GPL) - https://github.com/PennyWieser/DiadFit/blob/main/docs/Examples/EOS_calculations/FI_densities.xlsx
df_FI=pd.read_excel('FI_densities.xlsx', sheet_name='FI_data')
df_FI.head()
[12]:
| Density_g_cm3 | |
|---|---|
| 0 | 0.176605 |
| 1 | 0.180578 |
| 2 | 0.180858 |
| 3 | 0.147317 |
| 4 | 0.102414 |
[14]:
# We start the iteration with a calculation at XH2O=0
P_initial=pf.calculate_entrapment_P_XH2O(XH2O=0, CO2_dens_gcm3=df_FI['Density_g_cm3'], T_K=1200+273.15, T_K_ambient=37+273.15)
P=P_initial['P_kbar_mixCO2_DZ06_Hloss']*100
[15]:
# Loop 10 times
for _ in range(10):
XH2O_calc = Pf(P)
print(XH2O_calc)
Pressure_H_loss = pf.calculate_entrapment_P_XH2O(XH2O=XH2O_calc, CO2_dens_gcm3=df_FI['Density_g_cm3'],
T_K=1200 + 273.15, T_K_ambient=37 + 273.15).P_kbar_mixCO2_DZ06_Hloss
P = Pressure_H_loss * 100 # Update P for the next iteration
# Optional: Print values to track changes in each iteration
[0.02616126 0.02562257 0.02558634 0.031926 0.04871835 0.06122687
0.02328426 0.03806389 0.03871544 0.0275518 0.02906754 0.02620721
0.02759615 0.026134 0.02620046 0.025854 0.02618501 0.02637775
0.02890327 0.02612303 0.02944233 0.04910864 0.04902438 0.04347711
0.04346895 0.04301599 0.05048917 0.06524242 0.08810179 0.08907062
0.05684403 0.06110404 0.00952658 0.03477215 0.05515149 0.06786117
0.06472027 0.06472427 0.06215855 0.07779917 0.06735852 0.06625022
0.06652548 0.05867082 0.07294691 0.07104699 0.07350334 0.02173245
0.05160996 0.04294449 0.05336954 0.05315722 0.06499151 0.03809115
0.03076985 0.03732762 0.02468089 0.03068704 0.04238007 0.03711157
0.03314194 0.01761598 0.02266898 0.04853351 0.04925092 0.05023459
0.06987516 0.05483234 0.05707121 0.05062823 0.02252232 0.03709871
0.07250375 0.03585812 0.08084905 0.08065717 0.03585929 0.07264176
0.02499965 0.05609238 0.06596978 0.06112307 0.04234836 0.03814723
0.07123552 0.07877681 0.07340111 0.02225983 0.03828862 0.06231994
0.05520107 0.03846753 0.06337383 0.07552583 0.08978508 0.06171236
0.0328908 0.03684592 0.05771167 0.04135907 0.04889922 0.05924628
0.05839519 0.05342139 0.06106911 0.03691097 0.04167174 0.04217851
0.0236916 0.03990502 0.03841082 0.04062817 0.03704548 0.03532414
0.03878014 0.04318756 0.04230108 0.0417218 0.04535187 0.03569798
0.04024244 0.07309899 0.05973546]
Iteration 1: XH2O_calc = 0.02616125591287602, Pressure_H_loss = 0.5718466673784245
[0.02559573 0.02510724 0.02507439 0.03085334 0.04640834 0.05822758
0.02297478 0.03650191 0.03710402 0.02685868 0.02823934 0.02563742
0.02689902 0.025571 0.0256313 0.02531706 0.02561728 0.02579214
0.02808951 0.02556105 0.02858135 0.04677381 0.04669489 0.04151864
0.04151106 0.04109001 0.04806814 0.06207352 0.08456448 0.08554329
0.05406001 0.05811037 0.00847897 0.03346714 0.0524585 0.06459681
0.06157186 0.0615757 0.05911749 0.07429162 0.06411152 0.06304313
0.06330827 0.05579341 0.06953362 0.06768353 0.0700768 0.02144687
0.04912082 0.04102358 0.0507769 0.05057684 0.06183239 0.03652709
0.02979459 0.03582207 0.02425323 0.02971882 0.04049931 0.03562269
0.03196878 0.01683141 0.02239278 0.04623532 0.04690709 0.04782927
0.06654593 0.05215698 0.05427529 0.04819866 0.02225033 0.03561083
0.06910145 0.03446699 0.0773077 0.07711734 0.03446806 0.06923599
0.02454246 0.05334826 0.06277313 0.05812853 0.04046986 0.0365789
0.06786681 0.07525624 0.06997697 0.02199051 0.03670954 0.05927178
0.05250534 0.03687487 0.06028044 0.07205641 0.08626659 0.05869113
0.03173824 0.03537761 0.05488265 0.03955184 0.04657769 0.05634052
0.05553153 0.05082577 0.05807704 0.03543762 0.03984187 0.04031217
0.02335116 0.03820454 0.03682245 0.03887431 0.03556171 0.03397518
0.03716384 0.04124947 0.04042596 0.03988831 0.04326385 0.03431946
0.03851698 0.06968202 0.05680601]
Iteration 2: XH2O_calc = 0.02559572820733616, Pressure_H_loss = 0.5715116449327785
[0.02560767 0.02511735 0.02508438 0.03088906 0.04651978 0.05837955
0.02297882 0.03656614 0.0371713 0.02687578 0.02826257 0.02564952
0.02691629 0.02558285 0.02564337 0.02532795 0.0256293 0.02580485
0.02811206 0.02557286 0.02860616 0.04688682 0.04680756 0.04160776
0.04160014 0.04117707 0.04818654 0.06223356 0.08471682 0.08569333
0.05420023 0.05826205 0.00859845 0.03351594 0.05259343 0.06476065
0.061731 0.06173486 0.05927157 0.07445884 0.06427472 0.06320477
0.06347032 0.05593889 0.0697014 0.06785038 0.07024473 0.0214508
0.04924343 0.04111031 0.05090584 0.05070504 0.06199201 0.03659145
0.02982515 0.03588284 0.02426045 0.02974902 0.04058349 0.03568245
0.03201004 0.01686912 0.02239615 0.04634602 0.04702066 0.04794669
0.06671188 0.05229087 0.0544162 0.04831759 0.02225364 0.03567053
0.06926907 0.03452087 0.07747241 0.07728226 0.03452195 0.06940366
0.0245506 0.05348619 0.06293435 0.05828026 0.0405539 0.03664352
0.06803377 0.07542283 0.07014487 0.02199384 0.03677482 0.05942621
0.05264044 0.03694099 0.06043706 0.07222442 0.08641485 0.05884421
0.03177835 0.03543613 0.05502543 0.03963137 0.04668986 0.05648755
0.05567625 0.05095489 0.05822864 0.03549644 0.03992283 0.04039544
0.02335594 0.03827736 0.03688831 0.03895048 0.03562116 0.03402657
0.03723142 0.04133729 0.04050978 0.0399695 0.04336121 0.03437259
0.03859136 0.06984984 0.05695432]
Iteration 3: XH2O_calc = 0.02560766818913658, Pressure_H_loss = 0.5715187138578254
[0.02560742 0.02511715 0.02508418 0.03088787 0.04651441 0.05837187
0.02297877 0.03656349 0.03716849 0.02687536 0.02826192 0.02564926
0.02691586 0.0255826 0.02564312 0.02532773 0.02562905 0.02580457
0.02811143 0.02557261 0.02860544 0.04688135 0.04680211 0.0416037
0.04159609 0.04117313 0.04818076 0.06222549 0.08471028 0.08568697
0.05419318 0.05825438 0.00858488 0.03351412 0.05258668 0.06475244
0.06172297 0.06172683 0.05926377 0.07445089 0.06426654 0.06319664
0.06346218 0.05593154 0.06969317 0.06784212 0.07023652 0.02145075
0.0492374 0.04110639 0.05089943 0.05069868 0.06198396 0.0365888
0.02982419 0.03588039 0.02426033 0.02974808 0.04057972 0.03568005
0.03200859 0.01686732 0.02239611 0.04634069 0.04701516 0.04794097
0.06670363 0.05228418 0.05440911 0.04831177 0.0222536 0.03566814
0.06926083 0.03451878 0.07746477 0.0772746 0.03451987 0.06939543
0.02455046 0.05347926 0.06292624 0.05827259 0.04055014 0.03664086
0.06802552 0.07541497 0.07013666 0.02199379 0.03677212 0.0594184
0.05263368 0.03693824 0.06042915 0.0722163 0.08640862 0.05883647
0.03177695 0.0354338 0.05501824 0.03962787 0.04668444 0.05648012
0.05566895 0.05094847 0.05822097 0.03549409 0.03991925 0.04039173
0.02335587 0.03827424 0.03688558 0.03894717 0.03561878 0.03402461
0.03722859 0.04133331 0.04050604 0.03996591 0.04335667 0.03437054
0.03858816 0.06984162 0.05694682]
Iteration 4: XH2O_calc = 0.025607415970367434, Pressure_H_loss = 0.5715185649873413
[0.02560742 0.02511715 0.02508419 0.03088791 0.04651467 0.05837226
0.02297877 0.0365636 0.03716861 0.02687537 0.02826194 0.02564927
0.02691587 0.0255826 0.02564312 0.02532773 0.02562905 0.02580458
0.02811145 0.02557262 0.02860546 0.04688162 0.04680238 0.04160389
0.04159627 0.04117331 0.04818104 0.0622259 0.08471056 0.08568724
0.05419354 0.05825477 0.00858642 0.03351419 0.05258702 0.06475285
0.06172338 0.06172723 0.05926417 0.07445127 0.06426695 0.06319705
0.06346259 0.05593191 0.06969357 0.06784253 0.07023692 0.02145075
0.0492377 0.04110657 0.05089975 0.05069899 0.06198437 0.03658891
0.02982422 0.03588049 0.02426033 0.0297481 0.04057989 0.03568015
0.03200864 0.0168674 0.02239611 0.04634095 0.04701543 0.04794125
0.06670404 0.05228452 0.05440947 0.04831205 0.0222536 0.03566823
0.06926123 0.03451887 0.07746512 0.07727495 0.03451995 0.06939583
0.02455046 0.05347961 0.06292664 0.05827297 0.04055031 0.03664097
0.06802592 0.07541534 0.07013706 0.0219938 0.03677223 0.0594188
0.05263402 0.03693836 0.06042955 0.07221669 0.08640888 0.05883686
0.031777 0.03543389 0.0550186 0.03962803 0.0466847 0.05648049
0.05566931 0.05094879 0.05822136 0.03549419 0.03991941 0.04039189
0.02335588 0.03827437 0.03688569 0.03894731 0.03561887 0.03402468
0.03722871 0.04133349 0.0405062 0.03996607 0.04335688 0.03437062
0.03858829 0.06984202 0.0569472 ]
Iteration 5: XH2O_calc = 0.025607421281928694, Pressure_H_loss = 0.5715185680676197
[0.02560742 0.02511715 0.02508419 0.03088791 0.04651466 0.05837224
0.02297877 0.0365636 0.0371686 0.02687537 0.02826194 0.02564927
0.02691587 0.0255826 0.02564312 0.02532773 0.02562905 0.02580458
0.02811145 0.02557262 0.02860546 0.0468816 0.04680236 0.04160388
0.04159627 0.0411733 0.04818103 0.06222588 0.08471055 0.08568723
0.05419352 0.05825475 0.00858624 0.03351418 0.052587 0.06475283
0.06172336 0.06172721 0.05926415 0.07445125 0.06426693 0.06319703
0.06346256 0.05593189 0.06969355 0.06784251 0.0702369 0.02145075
0.04923768 0.04110656 0.05089973 0.05069898 0.06198435 0.0365889
0.02982422 0.03588049 0.02426033 0.0297481 0.04057988 0.03568014
0.03200864 0.0168674 0.02239611 0.04634094 0.04701542 0.04794123
0.06670402 0.0522845 0.05440945 0.04831204 0.0222536 0.03566823
0.06926121 0.03451886 0.07746511 0.07727494 0.03451994 0.06939581
0.02455046 0.05347959 0.06292662 0.05827295 0.0405503 0.03664097
0.0680259 0.07541532 0.07013704 0.0219938 0.03677222 0.05941878
0.052634 0.03693835 0.06042953 0.07221667 0.08640887 0.05883684
0.031777 0.03543389 0.05501858 0.03962802 0.04668469 0.05648048
0.0556693 0.05094878 0.05822134 0.03549418 0.0399194 0.04039189
0.02335588 0.03827437 0.03688569 0.03894731 0.03561887 0.03402468
0.0372287 0.04133348 0.0405062 0.03996606 0.04335687 0.03437062
0.03858829 0.069842 0.05694718]
Iteration 6: XH2O_calc = 0.025607421172027175, Pressure_H_loss = 0.5715185679542237
[0.02560742 0.02511715 0.02508419 0.03088791 0.04651466 0.05837224
0.02297877 0.0365636 0.0371686 0.02687537 0.02826194 0.02564927
0.02691587 0.0255826 0.02564312 0.02532773 0.02562905 0.02580458
0.02811145 0.02557262 0.02860546 0.0468816 0.04680236 0.04160388
0.04159627 0.0411733 0.04818103 0.06222588 0.08471055 0.08568723
0.05419352 0.05825475 0.00858626 0.03351418 0.052587 0.06475284
0.06172336 0.06172721 0.05926415 0.07445125 0.06426693 0.06319703
0.06346257 0.05593189 0.06969355 0.06784251 0.0702369 0.02145075
0.04923768 0.04110656 0.05089973 0.05069898 0.06198435 0.0365889
0.02982422 0.03588049 0.02426033 0.0297481 0.04057988 0.03568014
0.03200864 0.0168674 0.02239611 0.04634094 0.04701542 0.04794123
0.06670402 0.0522845 0.05440945 0.04831204 0.0222536 0.03566823
0.06926121 0.03451886 0.07746511 0.07727494 0.03451994 0.06939581
0.02455046 0.05347959 0.06292662 0.05827295 0.0405503 0.03664097
0.0680259 0.07541532 0.07013704 0.0219938 0.03677223 0.05941878
0.052634 0.03693835 0.06042953 0.07221668 0.08640887 0.05883684
0.031777 0.03543389 0.05501858 0.03962802 0.04668469 0.05648048
0.0556693 0.05094878 0.05822134 0.03549418 0.0399194 0.04039189
0.02335588 0.03827437 0.03688569 0.03894731 0.03561887 0.03402468
0.0372287 0.04133348 0.0405062 0.03996606 0.04335688 0.03437062
0.03858829 0.069842 0.05694718]
Iteration 7: XH2O_calc = 0.025607421176073022, Pressure_H_loss = 0.5715185684384301
[0.02560742 0.02511715 0.02508419 0.03088791 0.04651466 0.05837224
0.02297877 0.0365636 0.0371686 0.02687537 0.02826194 0.02564927
0.02691587 0.0255826 0.02564312 0.02532773 0.02562905 0.02580458
0.02811145 0.02557262 0.02860546 0.0468816 0.04680236 0.04160388
0.04159627 0.0411733 0.04818103 0.06222588 0.08471055 0.08568723
0.05419352 0.05825475 0.00858626 0.03351418 0.052587 0.06475284
0.06172336 0.06172721 0.05926415 0.07445125 0.06426693 0.06319703
0.06346257 0.0559319 0.06969355 0.06784251 0.0702369 0.02145075
0.04923768 0.04110656 0.05089973 0.05069898 0.06198435 0.0365889
0.02982422 0.03588049 0.02426033 0.0297481 0.04057988 0.03568014
0.03200864 0.0168674 0.02239611 0.04634094 0.04701542 0.04794123
0.06670402 0.0522845 0.05440945 0.04831204 0.0222536 0.03566823
0.06926121 0.03451886 0.07746511 0.07727494 0.03451994 0.06939581
0.02455046 0.05347959 0.06292662 0.05827295 0.0405503 0.03664097
0.06802591 0.07541532 0.07013704 0.0219938 0.03677223 0.05941878
0.052634 0.03693835 0.06042953 0.07221668 0.08640887 0.05883684
0.031777 0.03543389 0.05501858 0.03962802 0.04668469 0.05648048
0.0556693 0.05094878 0.05822134 0.03549418 0.0399194 0.04039189
0.02335588 0.03827437 0.03688569 0.03894731 0.03561887 0.03402468
0.0372287 0.04133348 0.0405062 0.03996606 0.04335688 0.03437062
0.03858829 0.069842 0.05694718]
Iteration 8: XH2O_calc = 0.025607421158796995, Pressure_H_loss = 0.5715185674605416
[0.02560742 0.02511715 0.02508419 0.03088791 0.04651466 0.05837224
0.02297877 0.0365636 0.0371686 0.02687537 0.02826194 0.02564927
0.02691587 0.0255826 0.02564312 0.02532773 0.02562905 0.02580458
0.02811145 0.02557262 0.02860546 0.0468816 0.04680236 0.04160388
0.04159627 0.0411733 0.04818103 0.06222588 0.08471055 0.08568723
0.05419352 0.05825475 0.00858626 0.03351418 0.052587 0.06475284
0.06172336 0.06172721 0.05926415 0.07445125 0.06426693 0.06319703
0.06346257 0.0559319 0.06969355 0.06784251 0.0702369 0.02145075
0.04923768 0.04110656 0.05089973 0.05069898 0.06198435 0.0365889
0.02982422 0.03588049 0.02426033 0.0297481 0.04057988 0.03568014
0.03200864 0.0168674 0.02239611 0.04634094 0.04701542 0.04794123
0.06670402 0.0522845 0.05440945 0.04831204 0.0222536 0.03566823
0.06926121 0.03451886 0.07746511 0.07727494 0.03451994 0.06939581
0.02455046 0.05347959 0.06292662 0.05827295 0.0405503 0.03664097
0.0680259 0.07541532 0.07013704 0.0219938 0.03677222 0.05941878
0.052634 0.03693835 0.06042953 0.07221668 0.08640887 0.05883684
0.031777 0.03543389 0.05501858 0.03962802 0.04668469 0.05648048
0.0556693 0.05094878 0.05822134 0.03549418 0.0399194 0.04039189
0.02335588 0.03827437 0.03688569 0.03894731 0.03561887 0.03402468
0.0372287 0.04133348 0.0405062 0.03996606 0.04335688 0.03437062
0.03858829 0.069842 0.05694718]
Iteration 9: XH2O_calc = 0.025607421193687155, Pressure_H_loss = 0.5715185674965405
[0.02560742 0.02511715 0.02508419 0.03088791 0.04651466 0.05837224
0.02297877 0.0365636 0.0371686 0.02687537 0.02826194 0.02564927
0.02691587 0.0255826 0.02564312 0.02532773 0.02562905 0.02580458
0.02811145 0.02557262 0.02860546 0.0468816 0.04680236 0.04160388
0.04159627 0.0411733 0.04818103 0.06222588 0.08471055 0.08568723
0.05419352 0.05825475 0.00858626 0.03351418 0.052587 0.06475284
0.06172336 0.06172721 0.05926415 0.07445125 0.06426693 0.06319703
0.06346257 0.0559319 0.06969355 0.06784251 0.0702369 0.02145075
0.04923768 0.04110656 0.05089973 0.05069898 0.06198435 0.0365889
0.02982422 0.03588049 0.02426033 0.0297481 0.04057988 0.03568014
0.03200864 0.0168674 0.02239611 0.04634094 0.04701542 0.04794123
0.06670402 0.0522845 0.05440945 0.04831204 0.0222536 0.03566823
0.06926121 0.03451886 0.07746511 0.07727494 0.03451994 0.06939581
0.02455046 0.05347959 0.06292662 0.05827296 0.0405503 0.03664097
0.0680259 0.07541532 0.07013704 0.0219938 0.03677222 0.05941878
0.052634 0.03693835 0.06042953 0.07221668 0.08640887 0.05883684
0.031777 0.03543389 0.05501858 0.03962802 0.04668469 0.05648048
0.0556693 0.05094878 0.05822134 0.03549418 0.0399194 0.04039189
0.02335588 0.03827437 0.03688569 0.03894731 0.03561887 0.03402468
0.0372287 0.04133348 0.0405062 0.03996606 0.04335688 0.03437062
0.03858829 0.069842 0.05694718]
Iteration 10: XH2O_calc = 0.025607421192402738, Pressure_H_loss = 0.5715185684552788
[16]:
## Lets compare results as a function of pressure
fig, (ax2) = plt.subplots(1, 1, figsize=(10,5))
ax2.plot( XH2O_calc, P_initial['P_kbar_pureCO2_DZ06']*100/P, '.r')
ax2.set_ylabel('P (Pure CO2/Mixed EOS)')
ax2.set_xlabel('XH2O')
#ax1.plot([0, 100], [0, 100], '-k')
[16]:
Text(0.5, 0, 'XH2O')
2. Calculating properties of mixed fluid.
Say you are an experimentalist, or working on another problem where you know P, T, XH2O, but want to calculate the properties of the mixed H2O-CO2 fluid (density, compressability, fugacity, molar volumes etc)
Download an example spreadsheet of the known conditions (e.g., experimental conditions here) https://github.com/PennyWieser/DiadFit/blob/main/docs/Examples/EOS_calculations/Exp_Cond.xlsx
[ ]:
## Lets load in a spreadsheet of experimental conditions - Get from here: https://github.com/PennyWieser/DiadFit/blob/main/docs/Examples/EOS_calculations/Exp_Cond.xlsx
df_exp=pd.read_excel('Exp_Cond.xlsx')
df_exp.head()
| Experiment | P_kbar | XH2O | T_K | |
|---|---|---|---|---|
| 0 | test1 | 1.0 | 0.1 | 1200 |
| 1 | test2 | 2.0 | 0.2 | 1400 |
| 2 | test3 | 4.0 | 0.0 | 800 |
| 3 | test4 | 1.0 | 0.0 | 1500 |
| 4 | test5 | 0.5 | 1.0 | 1200 |
[ ]:
df_calc=pf.calc_prop_knownP_EOS_DZ2006(P_kbar=df_exp['P_kbar'], T_K=df_exp['T_K'], XH2O=df_exp['XH2O'])
df_calc.head()
| P_kbar | T_K | XH2O | XCO2 | Molar Volume (cm3/mol) | Density (g/cm3) | Compressability_factor | fugacity_H2O (kbar) | fugacity_CO2 (kbar) | activity_H2O | activity_CO2 | |
|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 1.0 | 1200 | 0.1 | 0.9 | 125.333500 | 0.330319 | 1.256179 | 0.092945 | 1.163275 | 0.100954 | 0.900294 |
| 1 | 2.0 | 1400 | 0.2 | 0.8 | 86.992117 | 0.446017 | 1.494677 | 0.424996 | 2.668323 | 0.213429 | 0.804635 |
| 2 | 4.0 | 800 | 0.0 | 1.0 | 45.153235 | 0.974460 | 2.715342 | 0.000000 | 17.198213 | 0.000000 | 1.000000 |
| 3 | 1.0 | 1500 | 0.0 | 1.0 | 155.085905 | 0.283714 | 1.243502 | 0.000000 | 1.264244 | 0.000000 | 1.000000 |
| 4 | 0.5 | 1200 | 1.0 | 0.0 | 191.178832 | 0.094153 | 0.958063 | 0.477024 | 0.000000 | 1.000000 | 0.000000 |
[ ]:
# You can also do the calc for a single condition, say a melt inclusion where you know pressure, temp and XH2O from a solubility model
df_calc2=pf.calc_prop_knownP_EOS_DZ2006(P_kbar=3, T_K=1400, XH2O=0.2)
df_calc2
| P_kbar | T_K | XH2O | XCO2 | Molar Volume (cm3/mol) | Density (g/cm3) | Compressability_factor | fugacity_H2O (kbar) | fugacity_CO2 (kbar) | activity_H2O | activity_CO2 | |
|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 3 | 1400 | 0.2 | 0.8 | 67.007214 | 0.579042 | 1.726952 | 0.695761 | 5.283101 | 0.223567 | 0.806015 |