Database Interoperability
So far, we have looked at the possibility of creating and manipulating any species, whether they exist or not. If we wanted to create a H₂O⁺⁴ molecule, it would not be a problem. However, you will admit that it is a little strange...
This is why ChemistryLab relies on existing databases, in particular Cemdata18 and PSI-Nagra-12-07. Cemdata18 is a chemical thermodynamic database for hydrated Portland cements and alkali-activated materials. PSI-Nagra is a Chemical Thermodynamic Database. The formalism adopted for these databases is that of Thermofun which is a universal open-source client that delivers thermodynamic properties of substances and reactions at the temperature and pressure of interest. The information is stored in json files.
Naming a bundled database
Several databases ship with ChemistryLab, in its data/ directory. Name them with datapath, which returns an absolute path and therefore does not depend on the working directory:
using ChemistryLab
path = datapath("cemdata18-thermofun.json")
isfile(path)trueThe value is an absolute path, so it is not shown here: it depends on where the package is installed, and printing it would put this machine's directories into the page. What it points at does not depend on the machine:
relpath(path, pkgdir(ChemistryLab))"data/cemdata18-thermofun.json"readdir(datapath())12-element Vector{String}:
"CEMDATA18-31-03-2022-phaseVol.dat"
"aq17-thermofun.json"
"cemdata18-merged.json"
"cemdata18-thermofun.json"
"cemdata18-zeolites.json"
"experimental"
"pitzer-reardon1990.toml"
"psinagra-12-07-thermofun.json"
"slop98-inorganic-thermofun.json"
"slop98-organic-thermofun.json"
"solid_solutions.toml"
"zeolites"This matters in practice: a script written as build_species("data/cemdata18-thermofun.json") only runs when the working directory happens to be the repository root, which is not the case in an editor whose REPL started elsewhere, nor in a documentation build. Written with datapath, the same line runs from anywhere.
Resolution of a bundled name
The readers are tolerant, so older code keeps working: a path is tried against the working directory first, and only then against the bundled data/ directory. A call that already resolves therefore keeps resolving to exactly the same file, and a local database always takes precedence over a bundled one of the same name. The fallback can only succeed for a name that is one of the bundled files, so a mistyped path to a file of your own still fails, with an error listing what is available.
The zeolite extension of CEMDATA18
cemdata18-zeolites.json is one of the bundled files, and unlike the others it is generated by this repository rather than vendored from upstream. It is CEMDATA18 unchanged, with 28 zeolites appended.
The reason it exists is a phase-list error that looks like a solver error. A pozzolanic or an alkali-activated binder at high alkalinity precipitates zeolites; without them in the species list the alkalis have nowhere to go but the pore solution, and the calculated pH comes out too high. CEMDATA18 carries five zeolites, and the families a blended cement actually needs — clinoptilolite, heulandite, mordenite, phillipsite, analcime, stilbite and the gismondine/faujasite/LTA series, in both their Na and their K forms — are not among them.
The data are transcribed, number by number, from two open-access papers by the laboratory that produced CEMDATA18 itself: the Na series from (Ma and Lothenbach, 2020) and the K series from (Ma and Lothenbach, 2021).
using ChemistryLab
base = build_species(datapath("cemdata18-thermofun.json"); verbose = false)
ext = build_species(datapath("cemdata18-zeolites.json"); verbose = false)
length(base), length(ext)(228, 256)The extension is a strict superset: every symbol of the base file is present with the same values, so a script that loads it instead of the base gets the same answer unless it names one of the new phases.
added = setdiff(Set(symbol.(ext)), Set(symbol.(base)))
println(length(added), " phases added")
for s in sort(collect(added))
print(s, " ")
end28 phases added
ANA-Na CAN-NO3-Na CHA-K CHA-Na CLI-K FAU-X-K FAU-X-Na FAU-Y-K FAU-Y-Na GIS-LSP-K GIS-LSP-Na GIS-P1-K HEU-K LEU-K LTA-4A-Na LTA-K LTA-Na MOR-K MOR-Na NAT-K NAT-Na NAT-tetra-K PHI-K PHI-Na PHI-NaK SOD-Cl-Na SOD-OH-Na STI-KNothing is overwritten, so both values stay available
The papers revise ΔfG⁰ for phases CEMDATA18 already carries — natrolite by about 20 kJ/mol, from new solubility measurements. The new entries therefore take distinct symbols: NAT-Na sits beside natrolite, and the caller chooses. Declaring both in one system would count the same substance twice.
What makes the merge defensible, and what would make it wrong
Merging two thermodynamic datasets is only legitimate if they share a reference state, and the usual failure is silent: an offset of a few kJ/mol on Na+ moves every dissolution equilibrium by an order of magnitude without any solver complaining.
Both papers publish log Ksp and ΔfG⁰ for each phase, referred to the CEMDATA18 primary species. data/zeolites/regenerate.jl recomputes one from the other through CEMDATA18's own aqueous Gibbs energies and refuses to write the file if any phase misses by more than 0.05 log units. All 28 agree to within 0.026. It refuses on two further grounds: a symbol that would overwrite a CEMDATA18 entry, and a dissolution reaction that does not balance in elements and charge when re-derived from the formula string.
data/zeolites/README.md records the whole provenance, including the three other candidate datasets that were examined and rejected on measurement.
Loading species from a database
The simplest way to load species from a ThermoFun-compatible JSON file is build_species, which reads the file and directly returns a Vector{Species} with compiled thermodynamic functions:
using ChemistryLab
all_species = build_species(datapath("cemdata18-merged.json"))Each species already carries its molar mass and temperature-dependent thermodynamic functions (Cp⁰, ΔₐH⁰, ΔₐS⁰, ΔₐG⁰, logK⁰) as SymbolicFuncs and NumericFuncs.
Low-level access
If you need the raw DataFrames (e.g. to inspect metadata or filter on database columns), the lower-level function read_thermofun_database is still available and returns three DataFrames (df_elements, df_substances, df_reactions):
df_elements, df_substances, df_reactions = read_thermofun_database(datapath("cemdata18-merged.json"))build_species(df_substances) can then be called on the filtered DataFrame.
Filtering species with speciation
In practice, only a small subset of the database is relevant to a given problem. speciation filters a species list to those whose atomic composition is a subset of the atoms found in a set of seed species:
# Keep only species that can form from the calcite / water system
species_calcite = speciation(all_species, split("Cal H2O@");
aggregate_state=[AS_AQUEOUS],
exclude_species=split("H2@ O2@ CH4@"))
dict_species_calcite = Dict(symbol(s) => s for s in species_calcite)Dict{String, Species{Int64}} with 11 entries:
"Ca(HCO3)+" => Ca(HCO3)+ {CaHCO3+} [Ca(HCO3)+ ◆ Ca(HCO₃)⁺]
"CaOH+" => CaOH+ {CaOH+} [Ca(OH)+ ◆ Ca(OH)⁺]
"HCO3-" => HCO3- {HCO3-} [HCO3- ◆ HCO₃⁻]
"H+" => H+ {H+} [H+ ◆ H⁺]
"CO2@" => CO2@ {CO2 aq} [CO2@ ◆ CO₂@]
"Ca(CO3)@" => Ca(CO3)@ {CaCO3 aq} [CaCO3@ ◆ CaCO₃@]
"CO3-2" => CO3-2 {CO3-2} [CO3-2 ◆ CO₃²⁻]
"Cal" => Cal {Calcite} [CaCO3 ◆ CaCO₃]
"Ca+2" => Ca+2 {Ca+2} [Ca+2 ◆ Ca²⁺]
"OH-" => OH- {OH-} [OH- ◆ OH⁻]
"H2O@" => H2O@ {H2O l} [H2O@ ◆ H₂O@]The aggregate_state keyword restricts results to aqueous species (the seed species themselves — Cal and H2O@ — are always kept through include_species internally). For example, the properties of Ca(HCO₃)⁺ can then be read as:
dict_species_calcite["Ca(HCO3)+"]Species{Int64}
name: CaHCO3+
symbol: Ca(HCO3)+
formula: Ca(HCO3)+ ◆ Ca(HCO₃)⁺
atoms: Ca => 1, H => 1, C => 1, O => 3
charge: 1
aggregate_state: AS_AQUEOUS
class: SC_AQSOLUTE
properties: M = 0.10109399996506409 kg mol⁻¹
Tref = 298.15 K
Pref = 100000.0 m⁻¹ kg s⁻²
Cp⁰ = NumericFunc [m² kg s⁻² K⁻¹ mol⁻¹] ◆ vars=(T, P) ◆ T=298.15 K, P=100000.0 m⁻¹ kg s⁻²
ΔₐH⁰ = NumericFunc [m² kg s⁻² mol⁻¹] ◆ vars=(T, P) ◆ T=298.15 K, P=100000.0 m⁻¹ kg s⁻²
S⁰ = NumericFunc [m² kg s⁻² K⁻¹ mol⁻¹] ◆ vars=(T, P) ◆ T=298.15 K, P=100000.0 m⁻¹ kg s⁻²
ΔₐG⁰ = NumericFunc [m² kg s⁻² mol⁻¹] ◆ vars=(T, P) ◆ T=298.15 K, P=100000.0 m⁻¹ kg s⁻²
V⁰ = NumericFunc [m³ mol⁻¹] ◆ vars=(T, P) ◆ T=298.15 K, P=100000.0 m⁻¹ kg s⁻²
Cp⁰_Tref = 233.6999206543 m² kg s⁻² K⁻¹ mol⁻¹
ΔₐH⁰_Tref = -1.231942e6 m² kg s⁻² mol⁻¹
S⁰_Tref = 66.944000244141 m² kg s⁻² K⁻¹ mol⁻¹
ΔₐG⁰_Tref = -1.146041e6 m² kg s⁻² mol⁻¹
V⁰_Tref = 1.3329811096191e-5 m³ mol⁻¹speciation signatures
speciation accepts seed arguments in three forms:
| Seed argument | Description |
|---|---|
Vector{Symbol} | Explicit list of atom symbols |
Vector{<:AbstractSpecies} | Species objects — their union of atoms defines the space |
Vector{<:AbstractString} | Species symbol strings — looked up in species_list |
Common keyword arguments:
| Keyword | Default | Description |
|---|---|---|
aggregate_state | all states | restrict to [AS_AQUEOUS], [AS_CRYSTAL], etc. |
class | all classes | restrict to [SC_AQSOLUTE], etc. |
exclude_species | [] | species (or symbols) to always exclude |
include_species | [] | species to always include regardless of composition |
Primary species extraction
It is also possible to retrieve primary species from the Cemdata18 database. Primary species are a minimal subset such that every other species can be expressed as their linear combination.
df_primaries = extract_primary_species(datapath("CEMDATA18-31-03-2022-phaseVol.dat"))
show(df_primaries, allcols=true, allrows=true)Building solid solution phases from a TOML file
Solid solution phases (e.g. C-S-H gel, AFm) group several end-member species with a mixing model. Because database species carry class = SC_COMPONENT rather than class = SC_SSENDMEMBER, they cannot be passed directly to SolidSolutionPhase.
build_solid_solutions automates the full pipeline:
Reads phase definitions from a TOML file.
Looks up each end-member in a pre-built species dictionary.
Requalifies end-members to
SC_SSENDMEMBERviawith_class.Constructs and returns a
Vector{SolidSolutionPhase}.
TOML format
Each [[solid_solution]] entry specifies a phase name, the list of end-member symbols (as they appear in the database), and the mixing model:
# Ideal solid solution (any number of end-members)
[[solid_solution]]
name = "CSHQ"
end_members = ["CSHQ-TobD", "CSHQ-TobH", "CSHQ-JenH", "CSHQ-JenD",
"KSiOH", "NaSiOH"]
model = "ideal"
source = "Lothenbach2015"
# Binary Redlich-Kister (exactly 2 end-members, parameters in J/mol)
[[solid_solution]]
name = "AFm"
end_members = ["monosulphate12", "monocarbonate"]
model = "redlich_kister"
a0 = 3000.0
a1 = 500.0
a2 = 0.0
source = "Lothenbach2019"The file data/solid_solutions.toml shipped with ChemistryLab.jl contains pre-calibrated entries for the main cemdata18 solid solutions (CSHQ, AFm, Hydrogarnet, Ettringite_ss, Hydrotalcite).
Usage
using ChemistryLab, DynamicQuantities
substances = build_species(datapath("cemdata18-thermofun.json"))
dict = Dict(symbol(s) => s for s in substances)
# Build all solid solution phases defined in the TOML
ss_phases = build_solid_solutions(datapath("solid_solutions.toml"), dict)
# Use them directly in ChemicalSystem
cs = ChemicalSystem(species, CEMDATA_PRIMARIES; solid_solutions = ss_phases)Phases whose end-members are not found in dict are skipped with a warning (pass skip_missing = false to raise an error instead).
Manual requalification
If you only need one or two end-members, you can requalify them individually with with_class instead of going through the TOML file:
em = with_class(dict["CSHQ-TobD"], SC_SSENDMEMBER)