Solution model with vacancies, Fe2+ and Fe3+ and site partitioning

Question originally from Stan Roozen (solution model with vacancies, Fe2+ and Fe3+ and site partitioning · Issue #525 · geodynamics/burnman · GitHub):

I am currently working on fitting a solid solution model of tourmaline with 14 endmembers to thermodynamic measurements obtained from natural solid solutions (where I know the site occupancies via XRD). I have used your polytope script to construct a set of independent endmembers that describe the solution. However, when I fit my natural compositions to this set of independent endmembers many of the indep endmembers’ fractions are negative. I want to do a multiple linear regression to get endmember thermodynamic properties (endmember and if significant W parameters). However, I am still unsure about the validity of using negative independent endmember fractions directly in the fitting of the thermodynamic properties.

So far I have not found a paper that uses the (negative) indep endmember fractions directly in the fitting. I can see how it can turn into problems, for example:

The mineral composition F-dravite would be a linear combination of independent endmembers F-uvite-uvite+dravite. F-dravite would weigh in the regression of uvite by -1 although it has no uvite at all in its formula. Also, the margules parameters would get different meanings when you take for example the interaction parameter wuvite-dravite in case of the f dravite composition would weight in the regression of w by w*-1*1.

I think therefore (correct me if I am wrong or if there is an alternative way) that I should do the data fitting of endmember and W parameters using a set of positive independent endmembers.

I am trying to use the example_polytopetools.py to generate a polytope with independent endmembers which all have positive endmember fractions based on a given input bulk composition of the solid solution (one of my natural samples)

However, my solid solution has vacancies, Fe2+ /Fe3+, and site partitioning of elements between different sites. I could not really find an example of a mineral solid solution class with vacancies, different valence, and site partitioning to use as an example.

so far I have ignored the partitioning of the elements among the different sites and used the following code to try to get the set of independent endmembers whose endmember fractions are all positive.

tm_model2_burnman_template.txt

However, It keeps complaining about my user-defined element names like Vc (vacancy) and Fe_f (ferric iron), especially during ordering in the IUPAC order. Could you please let me know what I have to change in the code so I can use user-defined elements?

Additionally, in my extended model, I have to work with site partitioning and site fractions with user-defined names like Fe_X , Fe_Y, Mg_X, Mg_Y, Fef_X, Fef_Y ect. How do I set up a solution model class which allows for the partitioning of elements? Is there an example somewhere?

I would like to do the composite fitting to positive indep endmembers also with a model with site fractions. Also, I would like to do composition fitting like:
composition_fitting.txt
but with a solution model with site fractions.

Hi Stan,

However, when I fit my natural compositions to this set of independent endmembers many of the indep endmembers’ fractions are negative.

This is completely normal.

I am still unsure about the validity of using negative independent endmember fractions directly in the fitting of the thermodynamic properties.

It is completely valid.

The mineral composition F-dravite would be a linear combination of independent endmembers F-uvite-uvite+dravite.

Yes, and the Gibbs energy of F-dravite might also incorporate some mixing terms (see appropriate section of Myhill and Connolly (2021), which summarises and builds on work by Holland, Powell, Diener etc.).

… the margules parameters would get different meanings when you take for example the interaction parameter wuvite-dravite in case of the f dravite composition would weight in the regression

Yes, also this.

Could you please let me know what I have to change in the code so I can use user-defined elements?

The endmember mineral formulae can only use real elements (those in the IUPAC list). The site formulae, however, can use user-defined species names (Fef, Vac). So in your example, the ffoit and fluorbuergerite endmember formulae would be “Fe2Al7B3Si6H4O31” and “NaFe3Al6B3Si6O30F”, and then in the Solution you would use the site formulae. No need to declare the elements in the Solution class.

Your site formulae currently look like this:

[[schorl(), ‘[Na][Fe3/9Al6/9]9[B]3[Si]6[H]3[H]O31’],
[dravite(), ‘[Na][Mg3/9Al6/9]9[B]3[Si]6[H]3[H]O31’],
[uvite(), ‘[Ca][Mg4/9Al5/9]9[B]3[Si]6[H]3[H]O31’],
[ffoit(), ‘[Vc][Fe2/9Al7/9]9[B]3[Si]6[H]3[H]O31’],
[olenite(), ‘[Na][Al]9[B]3[Si]6[O]3[H]O28’],
[fluoruvite(), ‘[Ca][Mg4/9Al5/9]9[B]3[Si]6[H]3[F]O3O’],
[fluorbuergerite(), ‘[Na][Fe_f3/9Al6/9]9[B]3[Si]6[O]3[F]O27’],
[tidravite(), ‘[Na][Mg1/9Ti2/9Al6/9]9[B]3[Si]6[O]3[O]O27’],
[bolenite(), ‘[Na][Al]9[B]3[Si3/6B3/6]6[H]3[H]O31’]]

If you only have one species on a site, you can ignore the site. It looks like boron is the only species on site 3, so you might as well remove it in that model. BurnMan only reads the bits of the site formulae in square brackets and the number or fraction immediately after each bracket.

Your last two sites currently have H listed as an independent species. I guess by “H” you mean “OH”? It might be better to call this species “Oh” for clarity.

How do I set up a solution model class which allows for the partitioning of elements? Is there an example somewhere?

Site distributions of species in solutions are defined by their independent endmember proportions. You can specify any site-site partitioning you like by varying those proportions under constant bulk-compositional constraints. How you do that is a personal decision. When I’m fitting models, I usually find the equilibrium composition under the bulk compositional constraint by using burnman.equilibrate on the mineral of interest. This may not yield the global minimum Gibbs energy if the starting guess is inappropriate (part of the fun of nonlinear thermodynamics) - so sometimes it is worth brute-forcing the problem via grid search. The problem is harder if you also don’t know something else about the speciation (the ferric-ferrous ratio, for example). This is why full characterisation of experimental samples is so valuable.