describing atom style variable using python

LAMMPS LAMMPS Mailing List Mirror
1

Hello all,

I need to perform a dynamical simulation where the z-direction forces on particles are sampled using a python function. The python function that I am utilizing reads a tabular data of pre computed z-position and fz-values on a 1D grid, and returns an interpolated value of fz (float) for a given value of z (float). Here is a pseudocode

def interpolated1dfz(input_zpos):
read_input_data;
f = perform_interpolation;
fzval = f(input_zpos);
return fzval

Now, I want to use this function to create an atom style variable like

compute peratom ION property/atom z
variable samplefz atom interpolated1dfz(v_peratom)

This variable would be used somewhere in addforce fix later in the script.

Ignore the incorrect syntax on the variable command, but is it possible define atom style variable using it this way? Or do I have to update the per atom fz inside of python, by looping over all atoms in the group?

Please advise.

2

Hello all,

I need to perform a dynamical simulation where the z-direction forces on
particles are sampled using a python function. The python function that I
am utilizing reads a tabular data of pre computed z-position and fz-values
on a 1D grid, and returns an interpolated value of fz (float) for a given
value of z (float). Here is a pseudocode

def interpolated1dfz(input_zpos):
    read_input_data;
    f = perform_interpolation;
    fzval = f(input_zpos);
    return fzval

Now, I want to use this function to create an atom style variable like

compute peratom ION property/atom z
variable samplefz atom interpolated1dfz(v_peratom)

This variable would be used somewhere in addforce fix later in the script.

Ignore the incorrect syntax on the variable command, but is it possible
define atom style variable using it this way?

​no. python style variables act like equal style variables, i.e. global
​not per-atom expressions. please see the documentation for the variable
command.

​axel.​

3

Thanks Axel, is there any way you can point me to achieve the desired effect, i.e., adding force to the particles based on their position using an interpolated data.

Based on my understanding, I can define a compute in lammps input file and extract it in python

eng = lmp.extract_compute(id,style,type)

Now I can calculate the force values from interpolated data. The issue is, how do I communicate these sampled forces back to an atom style variable, and perform this every time step.

Thanks,

4

Thanks Axel, is there any way you can point me to achieve the desired
effect, i.e., adding force to the particles based on their position using
an interpolated data.

​interpolated how? and from what?​

Based on my understanding, I can define a compute in lammps input file and
extract it in python

eng = lmp.extract_compute(id,style,type)

Now I can calculate the force values from interpolated data. The issue is,
how do I communicate these sampled forces back to an atom style variable,
and perform this every time step.

​i don't understand the logic of this. you can access computes directly in
an atom-style variable​. why this obsession with python?

axel.

5

I have a z-position and mean-force at that z-position data stored in a file from an all atom MD simulation. I intend to run a fictitious dynamical simulation where the z-force on the particles of interest is sampled from this dataset. To keep things simple, if the z location (z_loc) of the particle falls between z_k and z_k+1 (k and k+1 being the indices from the stored data), the force is linearly interpolated between fz_k and fz_k+1 as:

fz_loc = (fz_k+1 - fz_k) / (z_k+1 - z_k) * (z_loc - z_k)

Though I can access computes directly in atom-style variable, the available arguments for atom-style variable does not support interpolation (as far as I am aware). Therefore, I was hoping to use python. If there is a cleaner approach, I would like to learn about it as well.

6

I have a z-position and mean-force at that z-position data stored in a
file from an all atom MD simulation. I intend to run a fictitious dynamical
simulation where the z-force on the particles of interest is sampled from
this dataset. To keep things simple, if the z location (z_loc) of the
particle falls between z_k and z_k+1 (k and k+1 being the indices from the
stored data), the force is linearly interpolated between fz_k and fz_k+1 as:

fz_loc = (fz_k+1 - fz_k) / (z_k+1 - z_k) * (z_loc - z_k)

Though I can access computes directly in atom-style variable, the
available arguments for atom-style variable does not support interpolation
(as far as I am aware). Therefore, I was hoping to use python. If there is
a cleaner approach, I would like to learn about it as well.

​the cleanest approach would be to implement this feature as a new fix
style or make it an extension of fix addforce​ in a similar fashion to
tabulation for various force styles (but ideally generalized into a 3d-grid
with the option to do spline interpolation as an alternative to linear
interpolation).

another, somewhat simpler way, would be to implement this as a callback
function to fix external. while a new fix is the most generic solution,
this one would require the least amount of programming, and you don't need
to know much about how to write code for LAMMPS.

if you absolutely want to do this via python, you would need to have the
python package installed *and* the python wrapper to the LAMMPS library, as
you'd have to primarily use the python style variable as a callback
mechanism that is triggered by a fix (e.g. fix addforce), but won't return
a force, but rather - through the python library wrapper - applied the
force directly to the force array. this is - in my opinion - the messiest
and least efficient approach.

axel.

7

Here are some other ideas.

If you wanted to add a force between run commands, you could

simply have Python invoke the runs, then use the lib interface

in Python to get a ptr to the force array for all the atoms,

and adjust the forces from Python. For anything but forces

that could work, but forces are zeroed and reset at each

timestep during a run, so that won’t work.

So you need someway to get adjusted forces from an external program

at the correct point in the timestep. This is exactly what

fix external does. Your driver program would have to define

the callback function and trigger LAMMPS to perform a (long) run

that makes the callbacks. That’s certainly possible from a C/F++/Fortran

code - see examples/couple/lmpqst.cpp that does this for LAMMPS

plus a quantum code to compute forces. It’s probably also possible

from a Python script, using the Python wrapper on LAMMPS as a lib.

You’d have to figure out how to specify the callback to a Python

function in C-syntax. I think I did that in another code; I’m almost

positive it’s possible.

Another option is to use the python command in LAMMPS.

As Axel said, the python variables LAMMPS provides

are like equal-style variables, i.e. a single value.

However the Python function they call can do anything.

So if you used fix addforce with fz specified as a python-style

variable, it would invoke the Python function you provide.

That function could extract a pointer to the current coords

and force for each atom. It could do the interpolation

for each atom and adjust the force on each atom. The

function could then return 0.0 which will be what fix addforce

gets as a result. It would then increment all the fz forces by 0.0,

i.e. a no-op.

I think that would work, but obviously haven’t tried it.

Steve

8

Thanks Axel and Steve, for the inputs.

I went ahead and modified the examples/couple/lmpqst.cpp file to suit my needs and I am able to drive the lammps with the resultant executable with z-force interpolated from the position and force dataset. I would like to get some clarifications on the original lmpqst program and my modifications at this point.

The first issue is to be sure if my modified forces with fix external are current. I am adding the following lines after defining fquest

if (info->me == 0) lmp2qst->gather(&f[0][0],3,&fquest[0][0]);
else lmp2qst->gather(&f[0][0],3,NULL);

which I believe are copying the forces from LAMMPS to quest array (I do this because I need x and y forces unchanged). Then a loop is run over atom indices i to modify only z component as

fquest[i][2] = fzpart(zpos, fzpos, xquest[i][2]);

where fzpart takes the [zpos, fzpos] map and returns interpolated value at xquest[i][2]. xquest is copied from LAMMPS in a similar fashion as fquest above. Later, fquest is passed to LAMMPS.

I have a few questions at this point:

Are the forces which are assigned to fquest from LAMMPS at the “current” time? Specifically, is the operation of pairwise computation and kspace calculations (plus other relevant bonded interactions etc) already been done before fix external is in play? Or the forces are set only using fix external, and no routine calculations are performed previously? The output trajectory shows non-zero values for fx and fy, but I need to be very sure.

Are the forces recomputed after the positions are updated and then dumped into the dump file? If yes, is fix external also in action at this point? Or just the routine computations (pair, bonded, kspace etc).

Thanks in advance,

9

Fix external adds its forces after the normal pair/bond/kspace forces have been

calculated. If you dump forces (which occurs at the end of the timestep),

you can easily debug whether the forces you want to add from your

external program have been added correctly to each atom.

Steve

10

Thanks a lot, the computation works as intended.