Recommendation for molecular motors in LAMMPS?

Hello.

I’ve been writing up simulations of a crosslinked, semiflexible polymer network modeled by bead-spring chains. Crosslinks are modeled as a distinct bond type between atoms rather than explicit molecules in themselves. I’d like to incorporate active forces as seen with molecular motors in cytoskeletal systems between pairs of atoms on separate polymers. Does anyone have any recommendations or solutions to modeling these effects in LAMMPS?

I see two options:

  1. Somehow select pairs of beads on separate filaments and use ‘addforce’ to create some force-extension function based on their separation
  2. Ideally, I’d like to be able to create a processive motor, effectively a harmonic bond that ‘walks’ along filaments. I would definitely need to code this up as a bond type that somehow has access to bond topology information that walks along a filament. However this may even require creating a new atom type that has some kind of topology that defines a polarity for the filaments, differentiating between bonding sites at the plus- and minus-ends. I’ve actually done this before with my own C++ code, but my MD integration was not implemented for HPC purposes and was not efficient at all.
  3. Model motors as discrete molecules/atom types that form bonds with ‘fix bond/create’ and ‘fix bond/break’ in simulation. Not sure how reliable this would be, but I’m planning on testing it anyway as it seems the easiest to write up.
    So I guess I’m asking for recommendations regarding what to implement, specifically which classes I should look to modify/customize in LAMMPS for options 1 and 2.

All help/pointers are appreciated.

Regarding your third option at least, fix bond/create and bond/break generally cannot be used at the same time. Instead, consider fix bond/react.

Regarding your third option at least, fix bond/create and bond/break generally cannot be used at the same time. Instead, consider fix bond/react.

I agree. For now, use fix bond/react. It is very powerful.

— avoid for now —

If you google, you may stumble upon examples of molecular motors simulated with LAMMPS on the internet, such as:
https://lammps.sandia.gov/workshops/Aug19/talk_jewett.pdf
https://www.youtube.com/watch?v=QO4LbHGAgxU
https://www.youtube.com/watch?v=EEbt07vZHew
https://www.youtube.com/watch?v=0cuEeCdyokU
(details here)
But I don’t currently recommend that you use the method described in these examples. The code that makes these examples run (“fix_bond_modify.cpp”, “fix_bond_new.cpp”) is not currently supported. Instead, fix bond/react can, in principle, do everything that is shown in these videos (and eventually more).

Andrew

P.S. I don’t think this is helpful for what you are currently trying to do, but for others who stumble upon this post, there are some examples of rotating motors using “fix twist”. (Note: “dihedral_style table” with the “linear” keyword can also be used to create rotating motors in LAMMPS as well, but it is more difficult to use)
https://www.youtube.com/watch?v=uAZ0R_oPv-4
https://www.youtube.com/watch?v=TzCJ4BTUQ20

I will support this code (“fix_twist.cpp”) and it will hopefully be included in future versions of LAMMPS.

For now, you can find the code for that here:
https://github.com/jewettaij/lammps/blob/fix_twist/src/USER-MISC/fix_twist.cpp
https://github.com/jewettaij/lammps/blob/fix_twist/src/USER-MISC/fix_twist.h

P.P.S. To make fix bond/react more convenient to use, in the future, I plan to add features to moltemplate.sh that will generate the “molecule template” files that fix bond/react needs for every step in the reaction. (Right now, moltemplate only creates LAMMPS data files, however the two file formats are pretty similar. Eventually, users will be able to describe each step in the reaction using syntax similar to the “MCA” syntax shown in the PDF file above.)