kaylai / vesical Goto Github PK

Explanations and examples of how to combine more than one model (e.g. one model for H2O, one for CO2), probably also with some text on why the code was written to allow this and some comments on extensibility.

Example:

m.calculate_dissolved_volatiles(sample=bulk_comp, temperature=900.0, pressure=2000.0, X_fluid=0.5, model='Dixon').result

Check all pressure units in bars for inputs and outputs, in MPa for MELTS calcs

Penny had the genius idea of using one of the fast models as an initial guess and then searching from that point.

Penny's comment:
"How much of the computation time is spend doing this step:

"The function works by calculating the equilibrium state of the given
melt at very high pressure (2,000 MPa) and then decreasing the pressure
in steps of 100 MPa until the mass of vapor is >0 grams."

If a lot, could it be sped up at least in some systems (e.g. basalts) by
using another one of the faster models to start at a more intelligent
pressure - E.g., if Shishkina returns 100 Mpa, Eguchi returns 200 Mpa,
could you start at one of those and work up and down to find the first
point where the mass of vapour is >0 grams?"

Will put this in "Future Features" for now, since it doesn't have to be implemented by v0.9, but it'd be nice if it were.

I am thinking ahead here and wondering if, when we really get down to comparing models, some of the differences will be due to the way preprocessing is done differently per model. Some may force normalization of the sample, others not. We should at least discuss and think about this such that the "differences" in model results we see are not partly an artifact of these choices.

Easy, right?

Doc string for calculate_equilibrium_fluid_comp is copy-pasted from calculate_saturation_pressure.

MooreWater is sort of a hybrid between a pure fluid model and mixedfluids model. In reality, when it calculates the fluid composition, for example, it's calculating what proportion of the fluid is H2O. The rest could be anything, I suppose (Gordon?).

Right now I am treating MooreWater as it's own style of class... not MixedFluids but not truly pure fluid either. It's fine, but a more elegant solution would be to standardize it to better integrate into the rest of the code. Not high priority since it all works fine for now.

Large values calculated for XH2Ofl... thought I fixed this yesterday...

Some of the MixedFluid (and the pure fluid methods it calls) will fail if the sample does not have both CO2 and H2O. Batch processing should be robust to this. The easiest way to control for this, I think, is to check for the presence of CO2 or H2O, and if only one is present, use a pure fluid model in place of the MixedFluid model. Additional consideration must be given where methods incorporate the concentration of the other volatile species, even in the pure fluid model.

When calculating volatile solubilities of a mixed fluid, the calculate_dissolved_volatiles methods of each model are called separately. In some models, the volatile solubility depends on other volatiles dissolved, e.g. Iacono-Marziano. The MixedFluids.calculate_dissolved_volatiles method should be able to deal with this.

Should check pressure, temperature, AND composition! We could write this to simply warn the user to please not use the models outside of their range (but allow them to). At some extreme point we can make the code simply refuse to do the calc! Thoughts?

Write code to export batch processed data.

Same structure as plot_isobars_and_isopleths()

I think an easy way of doing this is to have a bool statement before thermoengine is imported, which someone can configure in the python module file directly, if they are using it without thermoengine. I don't think this will cause any errors, unless thermoengine is specifically called.

However, I suspect there's a much more clever way of doing this.

The reason I think this is useful is for if people want to incorporate solubility code into other models (I've done this in the past, and is what my code was originally developed for...).

What are your thoughts? No urgency to doing this.

Low pressures are a source of problems. Their origins need to be constrained (perhaps scipy.optimize?) and something robust implemented to deal with it, particularly when calling functions like calculate_degassing_paths. At the moment I've hacked an inelegant solution here, but I think there will be a more robust generic way of dealing with this.

This is relatively low-priority. In the vast majority of realistic cases this will cause no problems.

If some samples in an ExcelFile object are outside of the calibration range for a particular model passed, how do we deal with this? One option is to return NaN and print something saying "Hey, by the way, we didn't calculate this sample because it's outside of the calibration range. If you really want to calculate it, call it as a single-sample".

Penny, I'm "assigning" this to you simply for visibility. We have been assigning each other to issues when we feel the other's input would be useful! --Kayla

The example_data.xlsx file should have a more diverse array of sample types! Right now it's a bunch of rhyolites...

Write code to perform calibration check on all core calculations.

This would maintain consistency with the MagmaSat method.

Simon, based on your email it sounds like the best way forward is for you to reformat your degassing paths methods to be more like MagmaSat. I'll assign this issue to you for now, but please let me know if I can tweak things, too.

Since we only allow a single degassing path to be calculated at once, maybe the methods should be named "calculate_degassing_path" instead of the current plural?

Once the mass balance equation has been assigned a number in the manuscript, the docstring should be updated to reference it correctly.

Dixon can return negative volatile concentrations, others probably can, too.

Standardize the option for arg mode='wtper' or mode='molfrac' for ExcelFile() class methods that return fluid compositions. Make molfrac the default.

Add mode arg to single-sample calculations.

Which models can perform which of the five core calculations? Which don't? Standardize.

Make a figure showing how reported data in the literature compares to that calculated using our code, to prove that it's working properly.

Small negative mass fractions are reported for the fluid at the start of the calculation, presumably from normalization discrepancies.

Currently non-MagmaSat models require two values passed as, e.g., X_fluid=(0.5,0.5). I assume this is (H2Ofluid, CO2fluid)? If so, isn't this redundant and we could simply pass one float? The MagmaSat model currently takes X_fluid as one float (X_fluid=0.5).

This becomes really important when comparing model results and plotting.

I have assigned this to all of us, as I think it will be fun for us to iterate over this until we find some really interesting results to show. I will take on the task of writing the first draft. We can tag up once all sections have a first draft and decide how to press on from there.

Implement fuzzywuzzy library to handle fuzzy string matching on oxide names for user input data.

See: https://docs.python.org/3.1/library/warnings.html

The MixedFluid.calculate_equilibrium_fluid_comp method assumes the volatile concentrations passed in sample are the system total volatile concentrations, not the concentrations dissolved in the fluid at the conditions of interest (and why it also requires a pressure to be chosen). To make this consistent with MagmaSat it should calculate assuming the concentrations are dissolved.

There could be an option for which way the function should be run, to keep this feature. Before changing, check that other methods don't rely on this particular way of implementing the method.

1. Imported data points (H2O, CO2) plus isobars, isopleths
2. Isobars, isopleths plus degassing path(s)
3. Others?

At the moment the calibration checks are 'hard coded' for each model, there's no reporting for the user anywhere of what the calibration limits are. This also means that the warnings will only tell you that you are outside a calibration range, but won't specify what it is...

e.g. the warning code produces this warning:
"Pressure is outside the calibration range of the CO2 Fugacity Model. Temperature is outside the calibration range of the CO2 Fugacity Model. Pressure is outside the calibration range of the H2O Fugacity Model. Temperature is outside the calibration range of the H2O Fugacity Model. "

This can be fixed by creating some kind of structure for storing the calibration limits, and then making the code check these. This will (potentially) limit the complexity of the possible checks, but will be much more useful for the user!

Make a table with all variable names, meaning, and usage.

We should discuss which license we want to use. Currently it is GPL, but this is self propagating, and I’m not sure we want that. MIT might be a better “standard” license to use. Thoughts?

How about a function within each model called get_calibration_range() that returns some statement about what kinds of samples should be used within that model?

The MagmaSat single-sample calculations return dictionaries. This was chosen because of the usage of the "verbose" argument. When verbose==True, a dictionary with lots of useful data is returned. Due to the large number of values returned, a dictionary makes the most sense for this case. When verbose==False, a simple dictionary is returned (consistent with verbose=True).

Non MagmaSat methods currently return series (value, value). Can we make these return dicts? Is this the best solution?