Setting an atom type rigid when certain conditions are met

Hello,

I’d like to be able to run an energy minimisation where atom type x is held rigid ONLY when it is bonded to atom type y. I have managed to do this in GROMACS, but haven’t found a way to do it in LAMMPS. Any suggestions?

Thanks,

Greig

The options available to LAMMPS for keeping bonds rigid are listed in the documentation of fix shake/rattle.

There is no direct way to do what you request to do, but there could be workarounds:

  • The simplest would be to write a script to process your data file and assign a new bond type to the bonds that you want to keep rigid based on your criteria and then adjust the number of bond types accordingly and copy over the BondCoeffs parameters.
  • A more complex, on the fly operation could be (this is untested) to use the python command to load a custom Python script and then use python invoke to execute that script. The script would then use the LAMMPS python module to create a group, gather the list of bonds, loop over them, extract atom types, and add all atoms where the desired conditions are met into that group and apply fix shake to that group only.

Thanks a lot.

To be clear, there already is a bond type associated with the atoms I’d like to hold rigid. So is there a way I could hold all atoms that are joined by a certain bond type rigid?

Your nomenclature is not clear. You cannot “hold atoms rigid”, only bonds or angles or entire groups of atoms (hence fix shake or fix rigid or one of their variants will be applied for that purpose).

People new to LAMMPS often confusingly refer to immobilizing atoms as “keeping them fixed” or “keeping them rigid”. The terms “rigid” and “fix(ed)” in the context of LAMMPS have specific meanings and thus should not be used to refer to immobilizing atoms. Please clarify.

Sorry, I did mean keeping the atoms immobile.

I have an atom type which I do want to be mobile, except when bonded to another specific atom type. The bond between these 2 atom types has a bond type associated with it. Can I use fix shake or fix rigid to keep all atoms involved in this bond type immobile? The sticking point is that it doesn’t look like there’s a way to create a group based on bond type?

My two suggestions for what you can do still apply with minor modifications.

  • You can still write a script that processes your data file and then rather writes out the atom IDs you want to keep immobile with a value of 1 to a file compatible with the atomfile variable style. If you then read the file into such a variable and use that variable to define a group, you have what you want.
  • You can do the same thing on-the-fly using the LAMMPS python command and module where you can grab the list of bonds via the gather_bonds() method or its NumPy variant and then compose a group command based on the extracted atom IDs.