create/destroy bonds on the fly

Hi all,
i'm new to LAMMPS and i'm working on its connection to our own python
program - we use LAMMPS to calculate forces for (parts of) our system.
Therefore, i want to use LAMMPS python connection (as a shared
library). Work-flow of our program is like this:

1) specify type of job (can be QM, QM/MM, adaptive QM/MM,
metadynamics....) & do initialization

2) if needed separate system into subsystems (to e.g. quantum part &
classical part) and pass coordinates & settings to external programs
that calculate forces

3) collect forces (and other interesting properties) from the external
programs

4) calculate forces combining forces from 3), do magic...

5) propagate system
  -> repeat from 2)

Therefore, for some type of jobs we need to call external programs fast
and that is why i'm working on connection to LAMMPS. While execution of
commands will be fast via the "lammps_command()" function from
library.cpp (also i added some code so i can pass input as a string so
I don't need to do it command at a time) I'm having problems with
system description, because especially bonds can be broken/created
during the simulation - this is what we do in the 4) "do magic" step -
we can break/create bonds there and we do new connection list in our
program and need to pass it to external MD programs. So my questions
are:
1) can i change bonds of already created system by specifying new
connection list in a way you do by read_data command in "bonds"
section, or in a way you change coordinates of atoms (
lammps.scatter_atoms('x',1,3, NewCoordinates) )

if 1) is not possible:
2) i could create new lammps() object each time step & destroy the old
one. However, is it possible to pass system description in the input
file (e.g. not by a command "read_data data.file" inside input
script file), or I'll need to dig with read_data.cpp?

i'm sorry if these questions (especially the second one) are part of
basic LAMMPS knowledge, but searching google & mailing list gave me no
answers.
thank you.
best,
stanislav.

I think you are saying that different invocatinos of LAMMPS,

e.g. to run a few steps of dynamics or energy minimization

are with essentially different systems. Maybe same # of atoms

but a new (changed) set of bonds.

In which case I think you need to re-create the system each

time with a data file that stores the topology. You don’t need

to create/destroy an instance of LAMMPS to do that. You

can just use the “clear” command, then read in a new system.

Note that the lib interface can take a file of LAMMPS commands

as input (i.e. an input script), so if your driver code can

write that file (and the data file it reads in), each time you

want to create a new system, it should work.

That is what I would do first. If there is some speed bottleneck

in writing/reading these files (doubtful if you are doing a reasonably

substantial computation with LAMMPS), and you are doing

something incremental (e.g. breaking one or two bonds, creating

one or two bonds), then it wouldn’t be hard to add some

hooks in the lib interface to allow you to tweak the topology

of those atoms. But I wouldn’t start that way …

Steve

Dear Steve,
thank you for answer.

I think you are saying that different invocatinos of LAMMPS,
e.g. to run a few steps of dynamics or energy minimization
are with essentially different systems. Maybe same # of atoms
but a new (changed) set of bonds.

Yes, there can be situations that subsequent invocations of LAMMPS
would need to work with system with changed bonds.

In which case I think you need to re-create the system each
time with a data file that stores the topology. You don't need
to create/destroy an instance of LAMMPS to do that. You
can just use the "clear" command, then read in a new system.

Note that the lib interface can take a file of LAMMPS commands
as input (i.e. an input script), so if your driver code can
write that file (and the data file it reads in), each time you
want to create a new system, it should work.
That is what I would do first. If there is some speed bottleneck
in writing/reading these files (doubtful if you are doing a
reasonably
substantial computation with LAMMPS),

Yes, our code is able to create input script & data file. I don't know
exactly how big will the overhead be for the LAMMPS, but for NAMD we're
having very big increase in wall time and writing inputs and reading
outputs is significant part of this increase (i'm comparing classical
MD simulation that is ran by namd vs the same simulation ran by our
code using namd to calculate forces). Therefore i would like to
directly pass connectivity table or even whole "data file" as a
string/matrix/whatever to some method of LAMMPS instead of writing it
to hard drive and then calling lammps_command("read_data datafile") -
we're doing this each time step...
As I mentioned - it's possible to execute commands directly -i.e.
without writing input script to HDD, by passing command to
"lammps_command()".
However, I don't know if I can (and how) pass "bonds", "angles"...
sections of data file in a similar way.

thank you.
best,
stanislav.

Therefore i would like to
directly pass connectivity table or even whole “data file” as a
string/matrix/whatever to some method of LAMMPS instead of writing it
to hard drive and then calling lammps_command(“read_data datafile”) -
we’re doing this each time step…
However, I don’t know if I can (and how) pass “bonds”, “angles”…
sections of data file in a similar way.

As I said, you can’t currently. You’d have to extend the

lib interface to use that kind of info. At a minimum

the system would have to re-neighbor if the bond

connectivity changed, since it affects the neighbor lists.

And the special exclusion lists for each atom would

have to be re-computed.

Until you have actual timing evidence that the simpler

solution (via files) is slow, I don’t suggest going

down this path.

Steve

OK, thank you.
best,
stanislav.