Issue/Question about bond/react

Hi,

I would like to use bond/react. I have looked through the examples, which worked fine, but I’m running into an issue in my case.

My system consists of linear polymer chains. I have bonds, angles, and dihedrals in the system. I want to change and ultimately break a bond in the linear chain.

For this case, I thought that I need to specify 8 atoms in the molecule files and map_file. This way the edge atoms are 3 bonds away from the bond ID atoms.
However, when I specify this the simulation crashes with an error of “bond atoms missing”. The reason for this is that the potential energy becomes really high. However, there is no reason the potential energy should be so high. In fact, even when the template-ID(pre-reacted) and template-ID(post-reacted) are the same, there is a huge increase in the pair-wise potential energy and the simulation crashes.

I tried only specifying 6 atoms in the map_file and the template-ID molecule files. If there is no change in atom types in the template-ID(post-reacted), this simulation runs. The fix bond/react shows the output of cumulative reactions occurring. However, even though the template-ID(post-reacted) had specified that a bond had changed type, no bonds changed type as seen by a data file written after the fix bond/react run a while. Thus, even though the thermo output of fix bond/react showed reactions happening, the bond change specified in the templates did not occur.

I tested that there was no issue in the example case when the template-ID(pre-reacted) and template-ID(post-reacted) were the same. I think there perhaps is something I’m missing with my linear polymer case.

I would like some guidance on this issue. I’m happy to provide any additional info/files that would be helpful in addressing this issue.

Thanks,

-Amulya

Hi Amulya,

You are correct, you must include atoms within three bonds of reacting atoms (when dihedrals are involved). Note that reacting atoms are all atoms that change type or connectivity, not just the ‘bond ID’ atoms responsible for triggering the reaction.

Bond/react should have thrown an error for the 6 atom case. It does not update topology too close to edge atoms. However, we did not think of the possibility that a bond would change type, but not the involved atoms’ types. Is this the case? (either way, I will add that additional check)

Please let me know if this helps. If that issue persists for the 8-atom-template system, you may need to send an example.

Thanks,

Jake

Thanks for your response. I’m actually still having trouble in the 8 atom case. I’ve attached relevant files here that contains this example. I specified an 8-atom topology in which the bond between atoms 4 and 5 changes to have a different bond constant. This for some reason results in bond atoms missing with the pair-wise potential energy increasing dramatically resulting in bad dynamics. I have no idea why this is happening and any suggestions would be very helpful.

Thanks,

rxn_unreacted.data_template (350 Bytes)

rxn_map (141 Bytes)

rxn_reacted.data_template (350 Bytes)

melt.data (93.4 KB)

in.melt_bond_react (1.69 KB)

Hi Amulya,

It appears we’ll need your supplementary files, for example containing the custom angle definitions, to test this example.

Jake