VESIcal Part I: An open-source thermodynamic model engine for mixed volatile (H2O-CO2) solubility in silicate melts

16 Thermodynamics has been fundamental to the interpretation of geologic data 17 and modeling of geologic systems for decades. However, more recent advance18 ments in computational capabilities and a marked increase in researchers’ acces19 sibility to computing tools has outpaced the functionality and extensibility of 20 currently available modeling tools. Here we present VESIcal (Volatile Equilibria 21 and Saturation Identification calculator): the first comprehensive modeling tool 22 for H2O, CO2, and mixed (H2O-CO2) solubility in silicate melts that: a) allows 23 users access to seven of the most popular models, plus easy inter-comparison be24 tween models; b) provides universal functionality for all models (e.g., functions 25 for calculating saturation pressures, degassing paths, etc.); c) can process large 26 datasets (1,000’s of samples) automatically; d) can output computed data into an 27 excel spreadsheet for simple post-modeling analysis; e) integrates advanced plot28 ting capabilities directly within the tool; and f) provides all of these within the 29 framework of a python library, making the tool extensible by the user and allow30 ing any of the model functions to be incorporated into any other code capable of 31 calling python. The tool is presented within this manuscript, which is a Jupyter 32 notebook containing worked examples accessible to python users with a range of 33 skill levels. The basic functions of VESIcal can also be accessed via a web app 34 (https://vesical.anvil.app). The VESIcal python library is open-source and 35 available for download at https://github.com/kaylai/VESIcal. 36 Plain Language Summary 37 Geologists use numerical models to understand and predict how volcanoes be38 have during storage (pre-eruption), eruption, and the composition and amount of 39 volcanic gas released into the atmosphere of Earth and other planets. Most mod40 els are made by performing experiments on a limited dataset and creating a model 41 that applies to that dataset. Some models combine lots of these individual models 42 to make a generalized model that can apply to lots of di↵erent volcanoes. Many of 43 these di↵erent models exist, and they all have specific uses, limitations, and pitfalls. 44 Here we present the first tool, VESIcal, which acts as a simple interface to seven 45 of the most commonly used models. VESIcal is written in python, so users can use 46 VESIcal as an application or include it in their own models. VESIcal is the first tool 47 that allows geologists to easily model thousands of data points automatically and 48 provides a simple platform to compare results from di↵erent models in a way never 49 before possible. 50


Iacono-Marziano et al., 2012
H2O -CO2 95 10,500 (mostly <5000) 1 1100-1400 (preferably 1200-1300) 2 Predominantly mafic compositions: subalkaline and alkaline basalts-andesites 1 Range of calibration dataset, as authors do not specifically state a calibration range. We note that the vast majority of experiments were conducted at <5000 bar. 2 Authors state that most experiments were conducted between 1200-1300 C (whole range 1100-1400 C_ Shishkina Shishkina  Note, all calculations performed at 1200 C (the experimental temperature). Authors suggest results generally applicable between 1000-1400 C mentary Jupyter Notebooks S1-S7), which allow users to plot their data amongst  To facilitate both ease of use and flexibility in model application, we have 324 structured the code such that users can follow one of two computational paths: a 325 batch processing path and a single sample path (more advanced options). The level and extensibility (Fig. 4).

374
The documentation for any function can be viewed in a jupyter notebook by typ-375 ing the function followed by a question mark and executing the cell (e.g., "v.

397
In each example below, a generic "method structure" is given along with def-  True  False  True  True  smooth isobars  True  smooth isopleths  True  fractionate vapor 0.0 init vapor 0.0 SS = Single-sample. Batch = batch processing. Color of cells corresponds to the type of argument: green=required; orange=optional; gray=argument not used. Values in cells indicate the unit or type of data to input for required arguments or the default value in the case of optional arguments.

Initialize packages 498
For any code using the VESIcal library, the library must be imported for use.

499
Here we import VESIcal as v. Any time we wish to use a function from VESIcal, 500 that function must be preceded by 'v.' (e.g., v.calculate_saturation_pressure shown. An example file called 'example data.xlsx' is included with this manuscript.

517
You can load in your own data by first ensuring that your file is in the same folder  sion .xls or .xlsx. The file must be laid out in the same manner as the example file 532 'example data.xlsx'. The basic structure is also shown in Table 3. 533 Any extraneous columns that are not labeled as oxides or input parameters 534 will be ignored during calculations. The first column titled 'Label' contains sample 535 names. Note that the default assumption on the part of VESIcal is that this column 536 will be titled 'Label'. If no 'Label' column is found, the first non-oxide column name 537 will be set as the index column, meaning this is how samples can be accessed by 538 name (see Section 3.1.3). An index column can be specified by the user using the 539 argument label (see documentation below). The following columns must contain 540 compositional information as oxides. The only allowable oxides are: SiO 2 , TiO 2 , and CO 2 . Currently, VESIcal can only read these oxide names exactly as written 543 (e.g., with no leading or trailing spaces and with correct capitalization), but func-544 tionality to interpret variations in how these oxides are entered is planned (e.g., such 545 that "sio2. " would be understood as "SiO2"). All of these oxides need not be in-546 cluded; if for example your samples contain no NiO concentration information, you 547 can omit the NiO column. Omitted oxide data will be set to 0 wt% concentration.

548
If other oxide columns not listed here are included in your file, they will be ignored 549 during calculations. Notably, the order of the columns does not matter, as they are 550 indexed by name rather than by position. Compositions can be entered either in 551 wt% (the default), mol%, or mole fraction. If mol% or mole fraction data are loaded, 552 this must be specified when importing the tile.

553
Because VESIcal may misread column headings, we highly recommend that 554 users examine their data after loading into VESIcal and before performing calcula-555 tions. The user data, as it will be used by VESIcal, can be viewed at any time with 556 manuscript submitted to Earth and Space Science For the rest of this manuscript, data will be pulled from the example data.xlsx  The bulk composition stored in a dictionary, with values in wt%. In some cases, it may be desired to simply retrieve a sample composition and 702 use it elsewhere. In case normalization is desired, the 'norm' argument can be used.

703
To specify a normalization style, for example the 'fixedvolatiles' normalization rou-704 tine, the above code could instead be written as: The user may wish to print extracted sample composition to a terminal or 746 notebook cell to verify that the correct data was extracted. Using the current ex-                   Once a single composition is defined, conditions over which to calculate isobars 1158 and isopleths must be specified. The generated plot is isothermal, so only one tem- Next, the H 2 O and CO 2 dissolved in the melt at saturation is calculated at the 1226 specified temperature and over the range of specified pressures. Note that, because 1227 this function calculates two things (isobars and isopleths), two variable names must 1228 be given (below, "isobars, isopleths"). This calculation can be quite slow, and so it 1229 is recommended to set print status to True. Text S2). If so desired, this calculation can be performed for any initial pressure, 1268 but the default is the saturation pressure. If a pressure is specified that is above the 1269 saturation pressure, the calculation will simply proceed from the saturation pressure, each plotted line will be given the generic legend name of "Isobars n", with 1371 n referring to the nth isobars passed. Isobar pressure is given in parentheses.

1372
The user can pass their own labels as a list of strings. If more than one set 1373 of isobars is passed, the labels should refer to each set of isobars, not each 1374 pressure.