Add superacids category: HF/SbF5 mixture densities and HF structure factor (Ref #379)

[mattiaperr]

Pre-review checklist for PR author

• I've confirmed the contribution guidelines.
• I have reviewed and understand all AI-generated code in this PR.
• I have added human-written tests for the new logic.
• I have properly cited any upstream algorithms or libraries the AI utilized.
• I have disclosed significant AI tool usage in the PR description.

Summary

test 1) HF/SbF5 mixture densities: NPT molecular dynamics at 288.65 K and 1 atm for three
compositions (pure HF, 10 mol % SbF5, pure SbF5), ~200 atoms each, 100 ps per system.
The density is taken from the average volume of the second half of the run, and the
metric is the MAPE against experimental densities (Shair and Schurig, 1951).

test 2) HF structure factor: a single NPT run of 100 HF molecules at 296 K and 1.2 bar for
50 ps. S(q) iis computed with the python library "MDANSE". Two metrics: the R-factor against experiment
(McLain et al., 2004) up to 4 A^-1, and the error on the position of the first peak.

Linked issue

Resolves #379

Progress

• Calculations
• Analysis
• Application
• Documentation

todo:

• HF_structure.zip has to be uploaded to S3

Testing

mace-mp-0a, mace-mp-0b3 and mace-mpa-0

New decorators/callbacks

Two changes to plot_scatter in ml_peg/analysis/utils/decorators.py, in their own commit:

1. highlight_range was being redrawn once per model, so the shaded band stacked up and the annotation was printed on top of itself. It is now drawn once per figure.
2. A new hlines option to draw dashed horizontal reference lines, used to mark the experimental target density on the density-vs-time plots. It defaults to None, so no existing plot is affected.

Other notes

This benchmark adds mdanse as a dependency, used for the structure factor.

Unit tests were added in tests/test_superacids.py covering the restart logic, the density calculation and the R-factor.

AI disclosure: I used an AI assistant (Claude) while working on this PR, mainly for refactoring and docstrings. All the physics of the simulations and the design of the workflow in the tests are mine, and I have reviewed every line.

[ElliottKasoar]

Thanks for this, @mattiaperr, this looks great so far! I've left a few minor comments.

If you can share the input data, it would be great to start testing this!

[joehart2001]

As this test involves running mutliple MDs in series, it could be a good idea to change parameterise this like in the liquid densities test

when you add the docs could you give an estiamte of how long the MD takes on GPU? if its not too expensive then it may not be woth doing this.

Otherwise, the calc, analysis and app seem to work well. For the analysis good and bad thresholds, is a 0% good value realistic due to error associated with the experimental measurements? For the app, we could consider some structure visualisation. as this is a trajectory we could visualise this, but this depends on what you think would be best for your test. Thanks!

[mattiaperr]

thanks!

I estimate 6 GPU hours for the pure HF system and 2 GPU hours for the pure SbF5 system, for a model that includes dispersion.
What do you think?

I wanted to follow up on this to discuss the structure of the test
I have developed, and currently use, a routine to estimate the HF structure factor and compare it with experimental data, in the style of the test: #445 .
Can I integrate such test into these same scripts, or should I create a separate calc/analysis/application within my same benchmark?
This would consist in a single NVT simulation at experimental density, which is much faster (25 ps to obtain the structure factor versus 100 ps for the density).

Given the concern about long simulation times, the whole benchmark could be reformulated as follows:

• if the HF density is too computationally expensive, use the structure factor as a metric (whose first peak position correlates with the density, and this fact can be taken into account as done in Add WiSE electrolyte benchmark (density, X-ray S(q), Li-O RDF) #445); and then perform only a single NPT simulation of the HF/Sbf5 density at a single point on the diagram.

What do you think?

[ElliottKasoar]

Hi @mattiaperr, thanks again for all the work on this!

In addition to the above comments, please can you take a look at our new filtering guidelines: https://ddmms.github.io/ml-peg/developer_guide/filter.html.

The main things to check are that:

1. The calculation is robust to errors (or this is suitably handled in analysis)
2. The analysis saves an info.json file along with tables etc
3. The analysis handles missing/incomplete data from models (and usually gives them a score of NaN/None)
4. The app reads the info.json file

I also wondered if you think there's anything else that would be interesting to visualise, beyond the scatter plot?

[ElliottKasoar]

thanks!

I estimate 6 GPU hours for the pure HF system and 2 GPU hours for the pure SbF5 system, for a model that includes dispersion. What do you think?

I wanted to follow up on this to discuss the structure of the test I have developed, and currently use, a routine to estimate the HF structure factor and compare it with experimental data, in the style of the test: #445 . Can I integrate such test into these same scripts, or should I create a separate calc/analysis/application within my same benchmark? This would consist in a single NVT simulation at experimental density, which is much faster (25 ps to obtain the structure factor versus 100 ps for the density).

Given the concern about long simulation times, the whole benchmark could be reformulated as follows:

• if the HF density is too computationally expensive, use the structure factor as a metric (whose first peak position correlates with the density, and this fact can be taken into account as done in Add WiSE electrolyte benchmark (density, X-ray S(q), Li-O RDF) #445); and then perform only a single NPT simulation of the HF/Sbf5 density at a single point on the diagram.

What do you think?

Sorry for such a delayed response. I've been working to introduce the element filtering I mention above, and wanted that in place before revisiting things to avoid things getting out of sync, and then I was away for quite a lot of June.

Personally I don't think it hurts to estimate the runtime, however quick it is, but I'd say that's definitely long enough that it's worth being clear.

Regarding the structure factor, I think it's somewhat up to you. In most cases, I'd say if it requires a different calculation (e.g. an entirely new MD simulation), rather than further analysis on the same data, it would usually be a different benchmark (so separate calc/analysis/app). If it's looking at exactly the same structure, we do have some examples of doing multiple types of calculation (e.g. iron_properties), so as long as it's set up in a way where we can run them individually if necessary (which might just be because one is labelled very slow and one isn't), we could do them all together.

I don't think we're overly concerned by the computational cost, we have tests which take days on GPU, I think @joehart2001 was making sure it's set up to run smoothly given the cost and that we document it, so as long as the cost tells us something interesting, I'd be inclined to keep it (but we can structure factor too).

[mattiaperr]

Thank you for your review!

The superacids category now has the density benchmark plus a separate one for the structure factor.

Main things implemented since the last review:

• restart handling for the MD runs
• errors are caught rather than raised
• element filtering: the analyses write info.json and both apps read it

As we discussed a while back, I used the Python MDANSE library to compute S(q). It is declared in pyproject.toml.
I am attaching HF_structure.zip that contains the input for runnining the structure factor benchmark.
HF_structure.zip

One thing to flag on shared code (analysis/utils/decorators.py): I added an hlines option to plot_scatter to draw the experimental target on the density-vs-time plots. It defaults to None, so no existing plot changes. The edits to this shared code are on a separate commit

Tested on mace-mp-0a, mace-mp-0b3 and mace-mpa-0. If this all looks fine, let me know whether you want to run the full set of models yourselves or whether I should do it.

I am also attaching the app JSONs, so you can look at the results: extract it into ml_peg/app/data/ and it creates the superacids/ folder.
superacids.zip