Hello,
I am trying to tether the CoM (centers-of-mass) of groups of atoms to specific points (fed as inputs). “Fix spring tether” does this. But, I would like to be unhindered by the 32 group limit.
The command “Fix spring/chunk” appears to be an alternative. However, this only tethers the instantaneous CoM to the initial CoM (and not to specified points).
I feel the easiest solution is to modify “Fix spring/chunk” to allow inputting a list of the coordinate of these specific points and bypassing the initial CoM calculation.
Please may get your thoughts on whether there are any better solutions I may be missing?
Also, I am new to LAMMPS development. So, any recommendations or templates for programming would be greatly appreciated. Looking at “Fix setforce” has helped. Is there any ‘correct’ way to input a long list of coordinates into a fix command?
Thank you!
This is usually done in LAMMPS by reading a file. See fix rigid, for example, which can read a list of center of mass, and moments of inertia from a file to override the values for rigid bodies built from overlapping extended particles.
I suggest to associate the position data with molecule IDs, which you can then also use to group the atoms with compute chunk/atom.
I would try this: add a dummy “anchor” particle to each molecule in the desired position. Make them all of a new type that isn’t used by any existing atom type so it’s easy to group them up into an anchors group
Then this code will set up a per-atom compute telling each atom the displacement between its molecular COM and its molecular “anchor”:
group atoms subtract all anchors
compute cmol all chunk/atom molecule ids once
compute anc anchors com/chunk cmol
compute com atoms com/chunk cmol
compute anc_atom all chunk/spread/atom cmol anc[*]
# you might be able to avoid recalculating with
# variable anc0 atom c_anc_atom
# variable anc atom v_anc0
# and then use v_anc instead of c_anc_atom
compute com_atom all chunk/spread/atom cmol com[*]
You can then impose the desired forces per-atom, for example
variable k equal 0.1
variable fx atom v_k*(c_com_atom[1]-c_anc_atom[1])
variable fy atom v_k*(c_com_atom[2]-c_anc_atom[2])
variable fz atom v_k*(c_com_atom[3]-c_anc_atom[3])
fix molecular_anchor atoms addforce v_fx v_fy v_fz
Then finally remember to not move the anchors, ideally by omitting them from your eventual integrator fix.
Of course you must think about how you want to “spread” the force on a molecule back to its individual atoms. But I think this would work. You wouldn’t have to write any new C++ code.
As a bonus, you can record your anchor positions as part of your trajectory and visualize them naturally in your desired software, instead of coming up with a new input format and then having to process that during visualization as well.
Dr. Kohlmeyer,
Thank you very much for the pointers. I can create an optional argument then in fix spring/chunk to allow for a center-of-mass file input.
Dr. Tee,
Thank you so much for elucidating the power of LAMMPS input file commands in constructing new fixes. For spreading the forces back onto the atoms, I would likely use a mass weighted distribution as is done in fix spring. If I am not wrong, the above provides the same force onto all atoms within a chunk? I will try to better understand your example. It is very interesting. Thank you!
Yes, the script section I wrote imposes identical force on all the atoms in a molecule. You can write a per-atom variable “telling” each atom its molecular mass with
variable m atom mass
# make sure you exclude the anchors from the following compute!
compute molmass atoms reduce/chunk cmol sum v_m
compute molmass_atom all chunk/spread/atom cmol c_molmass
and then you will be able to access both the per-atom masses and the molecular masses in the formula you wish to build, possibly
Dr. Kohlmeyer,
I am trying to do similar to what Sankha has mentioned.
I am trying to tether chunks to their initial COM (beginning of first run) for series runs having fix spring/chunk defined separately. Is there a way to keep same COM for subsequent runs or a way to input COM parameters from a file?
Please see the documentation about fix spring/chunk as it explains how it behaves in restart scenarios.
The challenge is not so much restarting the fix but having the exact same chunks defined. That depends strongly on information that you have not shared. Please see this recent discussion for an example of doing that for cases where just redefining compute chunk/atom will create different chunks: Grouping atoms across restarts
Thanks for the reply.
I am not having any problem in defining chunks, since in my case chunks are separate molecules (more than 32). With jump command, I am using fix spring/chunk multiple times in a simulation having different “k” values. I am trying to use same “R0m” values for all the sequential steps having fix spring/chunk. What is the best way to use same “R0m”?