[lammps-users] Is it possible to use bond/react with rigid, linear molecules?

Hello!

I am trying to use fix bond/react to model a system which has some reacting rigid molecules. I am having a hard time with fix bond/react because some of my molecules are linear.

To make the molecules rigid, I could use something like SHAKE or RATTLE or RIGID. But for linear (180 deg angle) molecules, SHAKE and RATTLE do not work, as described in the fix-shake doc page. So RIGID seems to be the only option for making molecules rigid.

When using fix bond/react, I first set up the reaction:

fix rxns all bond/react stabilization yes statted_grp .03 &

react reaction_stp1 all 1 0.0 3.0 mol.pre mol.post my.map

If I use RIGID/SMALL for time integration, with something like:

fix nvt_react statted_grp_REACT rigid/nvt/small molecule temp 300 300 100

I get an error because the group ‘statted_grp_REACT’ created by bond/react is a dynamic group (containing all the non-reaction-site particles), which it seems rigid/nvt/small cannot use (error: “Fix nvt_react does not allow use with a dynamic group”)

To recap my logic, I believe the only choice for rigid linear molecules is to use RIGID, but time integration with RIGID is not possible because the group created by bond/react is dynamic.

Is this train of logic correct? Or is there some way to simulate rigid, linear molecules using bond/react?

Thanks for your ideas!

Evan

I cannot talk much about fix bond/react, but I have some insight on the reasons for the issues you are seeing.
First off, using fix shake or fix rattle wouldn’t make a difference. The reason is that for both, shake/rattle constraints and rigid bodies, some information about the constrained/rigid objects is collected at the beginning of a run and then used throughout for efficiency reasons. That collides with applying them to dynamic groups. Instead, there are ways that fixes can tell fix shake or rigid/small that a rigid object has been added or removed (and those have to be defined as a molecule object). So that creates significant limitations.

So the only way that I can think of to implement a simulation like you are describing would be to use fix halt with some custom compute style that detects if a condition for a reaction is true. That would stop the simulation, you could manipulate your system as needed, and then continue while re-initializing the rigid bodies at the beginning of the run continuation.

But perhaps Jake has a better idea.

Axel.

Hi Evan,

This issue has come up before, and you summarized it well. I briefly tried to add ‘dynamic group’ support to rigid bodies but ran into the complications Axel mentioned, and I don’t see a straightforward way to solve it from the bond/react side. Perhaps I’ll take another look eventually, but for now reactive rigid bodies and bond/react don’t mix well.

Jake

Hello again, thanks to both of you for your ideas! I had a feeling this would be the case. I’ll think through some work-arounds.

Best,

Evan