Generating comparison potential between OpenFF-2.0.0 and Bespoke fit fragment

[tlfobe]

Hello! I attended a bespoke fit workshop last October and I'm now using it for some of my simulations. I had a question about one of the figures presented in that workshop.

I wanted to generate a figure similar to this comparing the torsions between OpenFF and BespokeFit (possibly also with QM data), but it looks like this figure includes more than just the torsion potentials (since they look more asymmetric and might include clashes from atoms as the torsion is rotated). What energetic contributions are being used in this figure?

Some comparison plots I've generated have only looked at the energy of torsion, but are rarely asymmetric.

[j-wags]

Not a complete answer, but I think the fact that there's a chiral carbon involved in the torsion is one thing that can cause the asymmetry.

@jthorton Would you happen to have the original inputs for this?

[tlfobe]

Hi @j-wags! Thanks for the response! I guess I was wondering if there is a way to evaluate how well the bespoke fit is to the selected fragments. For now, I've mostly used this package as a black-box workflow but if possible, we would like to compare the QM torsion distribution to the bespoke fit FF and the original FF we used (openff-2.0) similar to this figure. I know the QM calculations are cached, so if we could access those to look at the output energies being fit?

[jthorton]

Hi @tlfobe thanks for using BespokeFit! The plots from the paper include the full FF potential energy so they can be compared to the QM data. A notebook to reproduce this style of plot is included in the supporting data here called plot_bespoke_torsiondrive.ipynb. However, this is aimed at pulling records from QCArchive for reference data rather than a local BespokeFit instance. As you point out the data is cached and so it could be made to work with BespokeFit if needed.

Another option would be to tell BespokeFit to keep the ForceBalance fitting output before running a parameter optimisation this can be done by exporting the variable BEFLOW_OPTIMIZER_KEEP_FILES=True before the run. Then plots comparing the initial and current FF potential energy surface with the QM reference for each iteration of the optimisation can be found in a folder called optimize.tmp.

The last other option if you have already run the fits would be to regenerate the ForceBalance fitting input files using the following code, this will generate all of the input needed to run a ForceBalance single-point evaluation of the optimised FF. Note that as discussed in the paper these are all different ways of comparing the FF. The first example runs a full torsiondrive from the seed conformers, the second two ForceBalance options run constrained geometry optimisations from the QM optimised geometries with optional restraints.

# the bespokefit server must be running in the background
from openff.bespokefit.executor import BespokeExecutor
from openff.bespokefit.executor.services.coordinator.storage import get_task_ids
from openff.bespokefit.optimizers.forcebalance.factories import ForceBalanceInputFactory

# get all tasks from the optimizer
task_ids = get_task_ids()

result = BespokeExecutor.retrieve(optimization_id=task_ids[0])
opt_schema = result.results.input_schema.stages[0]
# turn off the restaint if you want to fully relax
opt_schema.targets[0].extras = {'restrain_k': 0.0}
opt_schema.parameters = []
 # generate the FB inputs with the bespoke FF in a folder called "test_bespoke"
ForceBalanceInputFactory.generate(root_directory='test_bespoke', schema=result.results.input_schema.stages[0], initial_force_field=result.bespoke_force_field)

you will then need to edit the optimize.in file and change the jobtype to single before running.

As you can see this is not the easiest to do currently, so let us know how you get on and any ideas you have to help make this easier!
I hope this helps!

[tlfobe]

Thank you so much for the help! When I ran bespoke-fit a month ago, it looks like some of the optimizations kept the optimize.tmp
directories for the fragments we were interested in and these had a set of comparison plots for each fragment that was used in the optimization.