"fix bond/break" or how to delete a specified bond

LAMMPS LAMMPS Mailing List Mirror
1

Dear everyone,

I believe this post is closely related to the problem I am encountering (a bit naive perhaps, sorry I am new to lammps).

Any answer/comment/advice is well appreciated. Thanks a lot.

Brief description of the problem (input and read file attached):

I want to sample the energy surface of 2 coarse grained interacting systems.

System 1 is composed of 3 atoms bonded between them (bond 1, bond 2 and angles 1).

System 2 is a single macroatom, interacting with a pair potential with the other 3 atoms in system 1.

To move system 2 around, I use a for loop over distance (no MD at the moment, in the long run I want to study transport properties).

So far everything is good.

I’d like to allow the bonded interaction in system 1 to be relaxed if the 4th atom comes too close to one of the atoms in system 1 (say < 4.5Ang from atom 2).

As suggested above, I try to define a new type (ID?) of bond, “3” in this case. The read_file.inp is modified as follow:

  • 2 bonds -> 2 bonds # no change

  • 2 bond types -> 3 bond types

  • Bond Coeffs # +1 line
    1 200 2.0 -> 1 200 2.0

2 200 3.0 -> 2 200 3.0

#nothing -> 3 100 3.0

In the lammps.inp I would then insert:

fix bond_2_close 2 bond/break 1 2 4 4.5

fix bond_2_far 3 bond/create 1 2 4 4.5

but I cannot fix the group ID properly (already in the first line).

"ERROR: Could not find fix group ID (…)

Last command: fix bond_2_close 2 bond/break 0 2 4 4.5"

So the question is what is the group ID of the bond I am going to break and the one I would create afterwards?

Given that I am not interested in dynamics at this point, I also thought to nest a if condition on the value of the bond in the for loop over distance, but I would ultimately end up at the same point.

Thank you.

Kind regards,

Marco

lammps.inp (1.37 KB)

read_file.in (457 Bytes)

2

Hi Marco,

LAMMPS is probably not the right tool for sampling the energy space of a 4-atom system, without running dynamics.

That said: The group ID must be the ID of a previously defined group, using the group command: https://lammps.sandia.gov/doc/group.html

You should review some example input scripts provided with the LAMMPS distribution. In your case, sounds like the group ID should be ‘all.’

Finally, you should not use fix bond/break and bond/create at the same time. For this, consider fix bond/react.

Thanks,

Jake

3

Hello Jake,

thank you for the mail.

LAMMPS is probably not the right tool for sampling the energy space of a 4-atom system, without running dynamics.

At this stage I want to optimize some coarse grained parameters to fit DFT results, before running transport calculations with MD. I am open to discuss on this point.

Finally, you should not use fix bond/break and bond/create at the same time. For this, consider fix bond/react.

To simplify the problem, let’s have only one bond (bond 1) to be modified when system 2 gets too close.

I now have 2 molecule templates (sys1.txt and sys2.txt attached) plus a map.txt file.

The only change to be made is bond1.

read_data read_file.txt
molecule sys1_hard sys1.txt
molecule sys1_weak sys1.txt
molecule sys2 sys2.txt

I could not insert “Bonds” in the molecular template (ERROR: Invalid atom ID in Bonds section of molecule file) so they are still in the read_file.txt. So, the main question is how to tell lamps that what I want to change is the bonds and not exchange atoms?

After some digging and hands on I ended up with more questions. I definitely need to re-iterate with you after this:

  1. when I fix bond/react as:

fix my_23 all bond/react stabilization no react my_rxn_23 all 1 4.0 7.0 sys1_weak sys1_hard map.txt

it only works for Rmin greater than Rmax, otherwise I get a "ERROR: Fix bond/react cutoff is longer than pairwise cutoff "

  1. To understand what is going on with the radius, I commented out the fix bond/react line and tried to print bonds with a compute but I get a “ERROR: Variable bond_dist: Mismatched compute in variable formula”:

this compute works

compute 1 all com

variable center_of_mass equal c_1[1]

fix print_com all print 1 "${center_of_mass} " file com.txt screen no title “com”

this does not work

compute 2 all bond/local force #dist engpot
variable bond_dist equal c_2[1]
fix print_bond all print 1 "${bond_dist} " file bond.txt screen no title “b1” ## this line crashes

So why is fix print_com working while print_bond is not?
(I also tried some dump command but I could not get any of this to work.)

  1. last, I would like to know if I can visualize sys.txt and read_file.txt with vmd or any other visualization tool.

All files are attached.

Many thanks

Marco

sys1.txt (139 Bytes)

sys2.txt (70 Bytes)

map.txt (72 Bytes)

read_file.txt (492 Bytes)

lammps.inp (1.4 KB)

4

first, you absolutely can specify bonds in a molecule template, so troubleshoot this first

you did not attach the molecule templates you used for bond/react. note, they must contain the same number of atoms

you should review the type of vector output from ‘compute bond/local.’ error messages explained here: https://lammps.sandia.gov/doc/Errors_messages.html

for visualizing molecule templates, I recently submitted an extension to Axel’s topotools: https://github.com/akohlmey/topotools
otherwise, you’ll have to convert to another format first