[lammps-users] running LAMMPS with forces from an external code

Is it possible to use lammps to perform dynamics with forces provided by an external (python) code?
If so, could you please suggest a viable strategy to implement it?

I have a machine learning code that provides the cartesian components of the force acting on each atom.

I am not interested in defining pair or bond terms.

Most of the computer time will be used by the machine learning code, so efficiency for the lammps part is not an issue. In other words, I’d rather have a simple solution than an efficient one.

Here is what I tried:

A) according to the manual, python codes can set string variables in lammps, so for N atoms I declare 3N ‘string’ variables of the form:

variable sf1x string 0.0 # string to contain fx for atom 1

and 3N ‘equal’ variables:

variable f1x equal v_sf1x # value of fx for atom 1

then I create N groups, one for each atom:

group a1 id 1

and activate a fix to set the force components

fix set1 a1 setforce v_f1x v_f1y v_f1z

next I declare a ‘python’ function to compute the forces:

python forces input 1 SELF format p file set_forces.py

here, the code in set_forces.py extracts the coordinates from lammps, computes the forces with a neural network, and sets the string variables with lines of the form:

lmp.set_variable(“sf1x”,"g" F[0,0])

the dynamics is then executed computing the forces at every step:

run 100 every 1 “python forces invoke”

This approach works only up to 31 atoms due to the MAX_GROUP limit of 32.
I tried recompiling after setting MAX_GROUP to a higher value, but groups past the 32nd are not properly constructed and contain too many atoms.
Is there a way to increase the number of groups, or to apply fixes directly to atoms?

B) If the forces are in vector f, then the python command

lmp.scatter_atoms(“f”,1,3,f)

will change the force values directly, i.e., without passing through the ‘string’ and ‘equal’ variables, the groups and the fixes.
Forces, however, are zeroed and recomputed inside lammps prior to each step, so this route does not work, at least without tampering with the ‘run’ command, which is probably beyond my capabilities.
Is there a simple way to interface python directly with a fix->post_force() function?

Thank you for any suggestion.

Cecco

Is it possible to use lammps to perform dynamics with forces provided by an external (python) code?

yes.

If so, could you please suggest a viable strategy to implement it?

please have a look at fix external. https://docs.lammps.org/fix_external.html
the callback function for fix external can also be a python function when it is registered from python using the LAMMPS python module.

this fix was meant to be used for the kind of usage scenario that you describe.

I have a machine learning code that provides the cartesian components of the force acting on each atom.

I am not interested in defining pair or bond terms.

Most of the computer time will be used by the machine learning code, so efficiency for the lammps part is not an issue. In other words, I’d rather have a simple solution than an efficient one.

Here is what I tried:

(deleted so that people don’t get the wrong idea)

This approach works only up to 31 atoms due to the MAX_GROUP limit of 32.
I tried recompiling after setting MAX_GROUP to a higher value, but groups past the 32nd are not properly constructed and contain too many atoms.
Is there a way to increase the number of groups, or to apply fixes directly to atoms?

no. this is a horrible approach on many levels, anyway. best you forget about it instantly.

B) If the forces are in vector f, then the python command

lmp.scatter_atoms(“f”,1,3,f)

will change the force values directly, i.e., without passing through the ‘string’ and ‘equal’ variables, the groups and the fixes.
Forces, however, are zeroed and recomputed inside lammps prior to each step, so this route does not work, at least without tampering with the ‘run’ command, which is probably beyond my capabilities.

modifying the flow of control in the MD loop is a bad idea as well.

Is there a simple way to interface python directly with a fix->post_force() function?

this is also possible using fix python/invoke: https://docs.lammps.org/fix_python_invoke.html

Thank you for any suggestion.

another option to couple two codes, without having to change LAMMPS or be “creative” in script programming is to use the MolSSI driver interface (MDI) and implement the MDI protocol into your python code. https://docs.lammps.org/Howto_mdi.html
this is a very new feature, so there may be some rough edges, but on the other hand, this is under active development and the authors are likely going to be interested in giving you a helping hand and have one more application to demonstrate the use of that feature.

axel.

That example is outdated. Please disregard it.

You should be either using the C library interface or the python module.
Both are documented online.

For the C library interface you only need to include library.h and for the python module you need no header. Both have a function to register a callback for fix external.

Dear Axel,

thank you for your reply. I tried to use ‘fix external’ but I’m having troubles passing the pointer for the callback function (pf/callback) or for the array (pf/array) to lammps.
I could not find instructions in the manual and the only examples with ‘fix external’ do this from C++.
For instance, the example examples/COUPLE/lammps_quest uses the C++ lines:

More specifically, please check out:

There is no nice instructive example at the moment, but you can look at the unit tests how these functions are used: