Skip to content

Commit

Permalink
Merge pull request #182 from thermotools/hydrate
Browse files Browse the repository at this point in the history
Initial hydrate support
  • Loading branch information
vegardjervell authored Dec 10, 2024
2 parents ca75d8e + 9e679b7 commit 42cff52
Show file tree
Hide file tree
Showing 38 changed files with 10,969 additions and 410 deletions.
2 changes: 1 addition & 1 deletion .github/prerelease_bodyFile.md
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@ For stable releases of `ThermoPack`, see releases named `vX.Y.Z`.
# Wheels for Python
To install this release, download and unzip the appropriate zip file, then install with `pip install thermopack -f wheel-<version>-<system>/`, where `wheel-<version>-<system>` is the directory created by unzipping the file.

For macOS, the wheels marked `macOS-12` are built for `x86_64` machines (using intel chips), while `macOS-latest` wheels are built for `arm64` machines (Apple Silicon, i.e. M1, M2, etc.).
Wheels for macOS are currently only built for `arm64` machines (Apple Silicon, i.e. M1, M2, etc.). Wheels for ThermoPack running on macOS with x86 chips (Intel), may be available for older versions of ThermoPack.

You can also use
```bash
Expand Down
14 changes: 7 additions & 7 deletions .github/workflows/cibuildwheel.yml
Original file line number Diff line number Diff line change
Expand Up @@ -18,13 +18,13 @@ jobs:
os: [macOS, ubuntu, windows]
version: [latest]
tp_version: [v2, v3]
include:
- version: 12
os: macOS
tp_version: v2
- version: 12
os: macOS
tp_version: v3
# include:
# - version: 12
# os: macOS
# tp_version: v2
# - version: 12
# os: macOS
# tp_version: v3

steps:
- uses: actions/checkout@v4
Expand Down
2 changes: 1 addition & 1 deletion .github/workflows/unittests.yml
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@ jobs:
strategy:
fail-fast: false
matrix:
os: [ubuntu-latest, macos-latest, macos-12]
os: [ubuntu-latest, macos-latest]
toolchain:
- {compiler: gcc, version: 13}
include:
Expand Down
10 changes: 8 additions & 2 deletions MSVStudio/thermopack.def
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
; File generated by thermopack/addon/pyUtils/exports/export_list.py
; Time stamp: 2024-10-04T15:56:57.529901
; Time stamp: 2024-12-09T14:49:32.643496
; @[email protected] : Declares the module parameters.

LIBRARY
Expand Down Expand Up @@ -98,6 +98,10 @@ EXPORTS
eostv_mp_enthalpy_tvp_
fundamental_measure_theory_mp_fmt_energy_density_
hardsphere_bmcsl_mp_calc_bmcsl_gij_fmt_
hydrate_mp_init_hydrate_model_
hydrate_mp_fugacity_water_in_hydrate_tpx_
hydrate_mp_fugacity_water_in_hydrate_tvn_
hydrate_curves_mp_map_hydrate_appearance_curve_
ideal_mp_set_standard_entropy_
ideal_mp_get_standard_entropy_
ideal_mp_set_enthalpy_of_formation_
Expand All @@ -117,6 +121,7 @@ EXPORTS
lj_splined_mp_ljs_wca_model_control_
lj_splined_mp_ljs_wca_set_pure_params_
lj_splined_mp_ljs_wca_get_pure_params_
multi_phase_envelope_tv_mp_multi_phase_envelope_plot_tv_
mut_solver_mp_solve_mu_t_
mut_solver_mp_solve_lnf_t_
mut_solver_mp_map_meta_isotherm_
Expand Down Expand Up @@ -196,6 +201,7 @@ EXPORTS
saturation_curve_mp_envelopeplot_
saturation_curve_mp_envelope_isentrope_cross_
saturation_curve_mp_pure_fluid_saturation_wrapper_
saturation_tv_mp_envelope_plot_tv_
saturation_point_locators_mp_locate_saturation_property_
saturation_point_locators_mp_property_index_from_string_
saturation_point_locators_mp_sat_points_based_on_prop_
Expand Down Expand Up @@ -229,4 +235,4 @@ EXPORTS
thermopack_var_mp_get_eos_identification_
tp_solver_mp_twophasetpflash_
sv_solver_mp_twophasesvflash_
uv_solver_mp_twophaseuvflash_
uv_solver_mp_twophaseuvflash_
50 changes: 46 additions & 4 deletions Makefile.code
Original file line number Diff line number Diff line change
Expand Up @@ -162,14 +162,16 @@ OBJE = \
$(ODIR)/puresaturation.o \
$(ODIR)/thermo_utils.o \
$(ODIR)/critical.o \
$(ODIR)/stability.o \
$(ODIR)/saturation.o \
$(ODIR)/solideos.o \
$(ODIR)/vls.o \
$(ODIR)/stability.o \
$(ODIR)/tuning.o \
$(ODIR)/tp_solver.o \
$(ODIR)/state_functions.o \
$(ODIR)/saturation_curve.o \
$(ODIR)/saturation_tv.o \
$(ODIR)/multi_phase_envelope_tv.o \
$(ODIR)/solid_saturation.o \
$(ODIR)/saturation_point_locators.o \
$(ODIR)/joule_thompson_inversion.o \
Expand All @@ -183,6 +185,8 @@ OBJE = \
$(ODIR)/binaryplot.o \
$(ODIR)/speed_of_sound.o \
$(ODIR)/saftvrmie_testing.o \
$(ODIR)/hydrate.o \
$(ODIR)/hydrate_curves.o \
$(ODIR)/consistency.o \
$(ODIR)/eoslibinit.o \
$(ODIR)/complexmodelinit.o \
Expand Down Expand Up @@ -503,6 +507,23 @@ $(ODIR)/hardsphere_wca.o: $(SRC)/hardsphere_wca.f90 \
$(SRC)/thermopack_constants.f90 $(SRC)/hardsphere_wca.f90 \
$(SRC)/saftvrmie_options.f90
@$(call make_object, $(SRC)/hardsphere_wca.f90)
$(ODIR)/hydrate.o: $(SRC)/hydrate.f90 \
$(SRC)/thermopack_constants.f90 $(SRC)/thermopack_var.f90 \
$(SRC)/eostv.f90 $(SRC)/saturation.f90 \
$(SRC)/stringmod.f90 $(SRC)/solideos.f90 \
$(SRC)/nonlinear_solvers.f90 $(SRC)/numconstants.f90 \
$(SRC)/quadratures.f90 $(SRC)/eos.f90
@$(call make_object, $(SRC)/hydrate.f90)
$(ODIR)/hydrate_curves.o: $(SRC)/hydrate_curves.f90 \
$(SRC)/hydrate.f90 $(SRC)/thermopack_constants.f90 \
$(SRC)/thermopack_var.f90 $(SRC)/eostv.f90 \
$(SRC)/saturation.f90 $(SRC)/multi_phase_envelope_tv.f90 \
$(SRC)/saturation_tv.f90 $(SRC)/nonlinear_solvers.f90 \
$(SRC)/stability.f90 $(SRC)/numconstants.f90 \
$(SRC)/eosdata.f90 $(SRC)/cubic_eos.f90 \
$(SRC)/utilities.f90 $(SRC)/eos.f90 \
$(SRC)/thermo_utils.f90
@$(call make_object, $(SRC)/hydrate_curves.f90)
$(ODIR)/hyperdual_mod.o: $(SRC)/hyperdual_mod.f90
@$(call make_object, $(SRC)/hyperdual_mod.f90)
$(ODIR)/hyperdual_utility.o: $(SRC)/hyperdual_utility.f90 \
Expand Down Expand Up @@ -562,6 +583,16 @@ $(ODIR)/mbwrdata.o: $(SRC)/mbwrdata.f90 \
$(ODIR)/mixdatadb.o: $(SRC)/mixdatadb.f90 \
$(SRC)/cubic_eos.f90 $(SRC)/assocschemeutils.f90
@$(call make_object, $(SRC)/mixdatadb.f90)
$(ODIR)/multi_phase_envelope_tv.o: $(SRC)/multi_phase_envelope_tv.f90 \
$(SRC)/eos.f90 $(SRC)/eostv.f90 \
$(SRC)/thermopack_constants.f90 $(SRC)/nonlinear_solvers.f90 \
$(SRC)/thermopack_var.f90 $(SRC)/numconstants.f90 \
$(SRC)/thermo_utils.f90 $(SRC)/saturation.f90 \
$(SRC)/saturation_tv.f90 $(SRC)/saturation_curve.f90 \
$(SRC)/cubic_eos.f90 $(SRC)/utilities.f90 \
$(SRC)/eosdata.f90 $(SRC)/stability.f90 \
$(SRC)/critical.f90
@$(call make_object, $(SRC)/multi_phase_envelope_tv.f90)
$(ODIR)/multiparameter_base.o: $(SRC)/multiparameter_base.f90 \
$(SRC)/numconstants.f90 $(SRC)/thermopack_constants.f90 \
$(SRC)/hyperdual_mod.f90
Expand Down Expand Up @@ -765,16 +796,17 @@ $(ODIR)/saturation.o: $(SRC)/saturation.f90 \
$(SRC)/eos.f90 $(SRC)/thermopack_constants.f90 \
$(SRC)/thermopack_var.f90 $(SRC)/nonlinear_solvers.f90 \
$(SRC)/numconstants.f90 $(SRC)/puresaturation.f90 \
$(SRC)/thermo_utils.f90 $(SRC)/eostv.f90
$(SRC)/thermo_utils.f90 $(SRC)/eostv.f90 \
$(SRC)/stability.f90
@$(call make_object, $(SRC)/saturation.f90)
$(ODIR)/saturation_curve.o: $(SRC)/saturation_curve.f90 \
$(SRC)/eos.f90 $(SRC)/thermopack_constants.f90 \
$(SRC)/thermopack_var.f90 $(SRC)/numconstants.f90 \
$(SRC)/puresaturation.f90 $(SRC)/saturation.f90 \
$(SRC)/nonlinear_solvers.f90 $(SRC)/critical.f90 \
$(SRC)/thermo_utils.f90 $(SRC)/eostv.f90 \
$(SRC)/eosdata.f90 $(SRC)/utilities.f90 \
$(SRC)/solideos.f90
$(SRC)/eosdata.f90 $(SRC)/solideos.f90 \
$(SRC)/stability.f90 $(SRC)/utilities.f90
@$(call make_object, $(SRC)/saturation_curve.f90)
$(ODIR)/saturation_point_locators.o: $(SRC)/saturation_point_locators.f90 \
$(SRC)/saturation.f90 $(SRC)/saturation_curve.f90 \
Expand All @@ -786,6 +818,16 @@ $(ODIR)/saturation_point_locators.o: $(SRC)/saturation_point_locators.f90 \
$(SRC)/thermo_utils.f90 $(SRC)/critical.f90 \
$(SRC)/solideos.f90 $(SRC)/vls.f90
@$(call make_object, $(SRC)/saturation_point_locators.f90)
$(ODIR)/saturation_tv.o: $(SRC)/saturation_tv.f90 \
$(SRC)/eostv.f90 $(SRC)/thermopack_constants.f90 \
$(SRC)/thermopack_var.f90 $(SRC)/nonlinear_solvers.f90 \
$(SRC)/numconstants.f90 $(SRC)/saturation.f90 \
$(SRC)/saturation_curve.f90 $(SRC)/puresaturation.f90 \
$(SRC)/cubic_eos.f90 $(SRC)/utilities.f90 \
$(SRC)/thermo_utils.f90 $(SRC)/eos.f90 \
$(SRC)/eosdata.f90 $(SRC)/critical.f90 \
$(SRC)/stability.f90
@$(call make_object, $(SRC)/saturation_tv.f90)
$(ODIR)/single_component.o: $(SRC)/single_component.f90 \
$(SRC)/eos_parameters.f90 $(SRC)/eosdata.f90 \
$(SRC)/mbwr_additional.f90 $(SRC)/compdata.f90 \
Expand Down
90 changes: 90 additions & 0 deletions addon/pyExamples/hydrate_curves.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,90 @@
#!/usr/bin/python
#Modify system path
import sys
sys.path.insert(0,'../pycThermopack/')
from thermopack.cubic import cubic
import numpy as np
import matplotlib.pyplot as plt
import itertools as it
import csv

cb = cubic("CO2,H2O","SRK")
cb.init_hydrate()

Z_H2O = np.array([30.0, 100.0, 250.0])*1e-6
initial_pressure = 1.0e3
minimum_temperature = 100.0
maximum_pressure = 150.0e5
p_scaling = 1.0e-6
colors = [ "black", "blue", "red", "green", "m"]
linestyles = [ "-", "--", ":", "-.", "."]

z = np.zeros((2))
for i in range(len(Z_H2O)):
z[0] = (1-Z_H2O[i])
z[1] = Z_H2O[i]
fluid, water = cb.get_multi_phase_envelope_tv(initial_pressure, z,
minimum_temperature,
maximum_pressure)
t_vals, p_vals = (fluid.t, fluid.p)
tw_vals, pw_vals = (water.t, water.p)
t_hyd_vals, p_hyd_vals = cb.get_hydrate_apperance_curve(initial_pressure, z,
minimum_temperature,
maximum_pressure)
plt.plot(t_vals, p_vals*p_scaling, linestyle="-", color=colors[i])
#plt.plot(tw_vals, pw_vals*p_scaling, linestyle=":", color=colors[i],
# label="H2O: {} ppm".format(round(z[-1]*1e6)))
plt.plot(t_hyd_vals, p_hyd_vals*p_scaling, linestyle="--", color=colors[i],
label="Hyd. H2O: {} ppm".format(round(z[-1]*1e6)))

# z = np.ones((2))*0.5
# t_hyd_vals, p_hyd_vals = cb.get_hydrate_apperance_curve(initial_pressure, z,
# minimum_temperature,
# maximum_pressure)
# plt.plot(t_hyd_vals, p_hyd_vals*p_scaling, linestyle="--", color="orange",
# label="Hyd. excess H2O")

plt.title("SRK: CO2 and variable H2O")
leg = plt.legend(loc="best", numpoints=1)
leg.get_frame().set_linewidth(0.0)
plt.ylabel(r"$P$ (MPa)")
plt.xlabel(r"$T$ (K)")
plt.tight_layout()
plt.savefig("hydrate_co2.pdf")

plt.figure()
cb.init("CO2,N2,H2O","SRK","Classic","Classic")
cb.init_hydrate()
Z_CO2 = np.array([0.85,0.15])
Z_H2O = np.array([30.0, 100.0, 250.0])*1e-6
initial_pressure = 5.0e3

z = np.zeros((3))
for i in range(len(Z_H2O)):
z[0:2] = (1-Z_H2O[i])*Z_CO2
z[2] = Z_H2O[i]
fluid, water = cb.get_multi_phase_envelope_tv(initial_pressure, z,
minimum_temperature,
maximum_pressure)
t_vals, p_vals = (fluid.t, fluid.p)
tw_vals, pw_vals = (water.t, water.p)
t_hyd_vals, p_hyd_vals = cb.get_hydrate_apperance_curve(initial_pressure, z,
minimum_temperature,
maximum_pressure)

plt.plot(t_vals, p_vals*p_scaling, linestyle="-", color=colors[i])
# plt.plot(tw_vals, pw_vals*p_scaling, linestyle=":", color=colors[i],
# label="H2O: {} ppm".format(round(z[-1]*1e6)))
plt.plot(t_hyd_vals, p_hyd_vals*p_scaling, linestyle="--", color=colors[i],
label="Hyd. H2O: {} ppm".format(round(z[-1]*1e6)))


plt.title("SRK: CO2 (0.85), N2 (0.15) and variable H2O")
leg = plt.legend(loc="best", numpoints=1)
leg.get_frame().set_linewidth(0.0)
plt.ylabel(r"$P$ (MPa)")
plt.xlabel(r"$T$ (K)")
plt.tight_layout()
plt.savefig("hydrate_co2_n2.pdf")
plt.show()
plt.clf()
88 changes: 88 additions & 0 deletions addon/pyExamples/multi_phase_curves.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,88 @@
#!/usr/bin/python
#Modify system path
import sys
sys.path.insert(0,'../pycThermopack/')
from thermopack.cubic import cubic
import numpy as np
import matplotlib.pyplot as plt
import itertools as it
import csv

cb = cubic("CO2,N2,H2O","SRK")
Z_CO2 = np.array([0.85,0.15])
Z_H2O = np.array([0.05,1e-2,1e-3,1e-4])
Z_H2O = np.array([5.0, 30.0, 100.0, 1000.0])*1e-6
initial_pressure = 5.0e3
minimum_temperature = 140.0
maximum_pressure = 150.0e5
p_scaling = 1.0e-6
colors = [ "black", "blue", "red", "green"]
linestyles = [ "-", "--", ":", "-."]

z = np.zeros((3))
for i in range(len(Z_H2O)):
z[0:2] = (1-Z_H2O[i])*Z_CO2
z[2] = Z_H2O[i]
fluid, water = cb.get_multi_phase_envelope_tv(initial_pressure, z,
minimum_temperature,
maximum_pressure)
t_vals, p_vals = (fluid.t, fluid.p)
tw_vals, pw_vals = (water.t, water.p)
plt.plot(t_vals, p_vals*p_scaling, linestyle="-", color=colors[i])
plt.plot(tw_vals, pw_vals*p_scaling, linestyle="--", color=colors[i],
label="H2O: {} ppm".format(round(z[-1]*1e6)))

cb.init("CO2,N2","SRK","Classic","Classic")
Z_CO2 = np.array([0.85,0.15])
T, P, v = cb.get_envelope_twophase_tv(initial_pressure, Z_CO2,
maximum_pressure=maximum_pressure,
minimum_temperature=minimum_temperature)
plt.plot(T, P*p_scaling, linestyle="-", color="cyan", label="No water")

plt.title("SRK: CO2 (0.85), N2 (0.15) and variable H2O")
leg = plt.legend(loc="best", numpoints=1)
leg.get_frame().set_linewidth(0.0)
plt.ylabel(r"$P$ (MPa)")
plt.xlabel(r"$T$ (K)")

plt.figure()
cb.init("CO2,H2O","SRK","Classic","Classic")
initial_pressure = 1.0e3
minimum_temperature = 100.0
maximum_pressure = 150.0e5
p_scaling = 1.0e-6
colors = [ "black", "blue", "red", "green"]
linestyles = [ "-", "--", ":", "-."]

arrays = []
z = np.zeros((2))
for i in range(len(Z_H2O)):
z[0] = (1-Z_H2O[i])
z[1] = Z_H2O[i]
fluid, water = cb.get_multi_phase_envelope_tv(initial_pressure, z,
minimum_temperature,
maximum_pressure)
t_vals, p_vals = (fluid.t, fluid.p)
tw_vals, pw_vals = (water.t, water.p)
arrays += [t_vals, p_vals, tw_vals, pw_vals]
plt.plot(t_vals, p_vals*p_scaling, linestyle="-", color=colors[i])
plt.plot(tw_vals, pw_vals*p_scaling, linestyle="--", color=colors[i],
label="H2O: {} ppm".format(round(z[-1]*1e6)))

# with open('co2_h2o.csv', 'w') as f:
# csv.writer(f, delimiter='\t').writerows(it.zip_longest(*arrays, fillvalue =np.NaN))

cb.init("CO2","SRK","Classic","Classic")
Z_CO2 = np.array([1.0])
T, P, v = cb.get_envelope_twophase_tv(initial_pressure, Z_CO2,
maximum_pressure=maximum_pressure,
minimum_temperature=minimum_temperature)
plt.plot(T, P*p_scaling, linestyle="-", color="cyan", label="No water")

plt.title("SRK: CO2 and variable H2O")
leg = plt.legend(loc="best", numpoints=1)
leg.get_frame().set_linewidth(0.0)
plt.ylabel(r"$P$ (MPa)")
plt.xlabel(r"$T$ (K)")

plt.show()
4 changes: 4 additions & 0 deletions addon/pyExamples/properties.py
Original file line number Diff line number Diff line change
Expand Up @@ -46,7 +46,11 @@
print("Compressibility factor: {}".format(zFac))

sos = srk.speed_of_sound(T, P, x=z, y=z, z=z, betaV=0.0, betaL=1.0, phase=srk.LIQPH)
sos_z = srk.speed_of_sound_tv(T, v, n=z)
sos_n = srk.speed_of_sound_tv(T, v*2, n=z*2)
print("Saturated liquid speed of sound: {} (m/s)".format(sos))
print("Saturated liquid speed of sound: {} (m/s)".format(sos_z))
print("Saturated liquid speed of sound: {} (m/s)".format(sos_n))

# Properties evaluated in temperature volume and mol numbers
p, dp = srk.pressure_tv(T, v, z, dpdv=True)
Expand Down
Loading

0 comments on commit 42cff52

Please sign in to comment.