Dear all,
I want to simulate a long chain polymer. Since it is repeatable for its structure, I want to cut it into a smaller piece and use periodic boundaries to wrap both ends. But I am not sure if I joint the two ends, will LAMMPS be able to tell that they are just the end of one box and the beginning of another box? Or the program will still bend the chain to joint the both ends. Because on principle, by the definition of periodic boundaries, any distances larger than half box size should be wrapped back and make it smaller than that.
Thanks in advance for any suggestions or ideas.
Best wishes,
Andy
Hello Andy, the minimum-image convention will apply to each individual distance. The wrapped coordinates will appear like there is at least one very long bond per each chain, but for the sake of computing the bond potential and the other bonded terms, LAMMPS will compute distances between atoms modulo the box dimensions.
Giacomo
Hello Andy, the minimum-image convention will apply to each individual
distance. The wrapped coordinates will appear like there is at least one
very long bond per each chain,
That's true. I should have mentioned this. When I use VMD (for
example), I typically delete the bond near the boundary from the PSF
file, that I use for visualization. I do this to prevent it from
stretching that bond across the screen. Usually this looks okay for
short simulations. But if the polymer moves substantially during the
course of the simulation., and a different bond is straddling the
boundary, it can look weird. OVITO might work better than VMD for
visualzing these kinds of simulations, because it does not draw bonds
between atoms unless they appear to be within a (user-definable)
cutoff distance.
but for the sake of computing the bond
potential and the other bonded terms, LAMMPS will compute distances between
atoms modulo the box dimensions.
It is interesting. A few years ago most of the bonded potentials used
"domain->minimum_image()" to always select the periodic image that was
closest to the bond partner. That changed a little while. I think
LAMMPS is now using the neighborlist code to take care of this detail.
Since that change was made, sometimes LAMMPS will get this wrong:
(In that example, there were two short bonds, initially of almost the
same length, but one was longer than the "comm_modify" distance. Both
of these pairs of atoms (one pair per bond) were initially located on
the boundary. Forces from the bond whose length exceeded the
"comm_modify" distance was calculated using the wrong periodic image,
and LAMMPS thought that the bond was stretched across the simulation
box. After letting the system relax, the two atoms in that bond were
pulled away from the boundary to the center. Normally this would
cause a crash.)
So I'm mentioning this again if anyone runs into this problem and
googles this discussion thread. (Meanwhile, I should probably post a
feature request on github, as Axel suggested.)
Cheers
Andrew
Hello Andrew and Giacomo,
Now I run into another question. Will the NPT be sufficient to allow flexibility of the polymer? For instance, my polymer is lying on the z direction and I couple the x y pressures. Will my polymer behave like a rubber band on z direction? I question this because there should be solvents which are averaged to the total z direction pressure and I am afraid that the tension from polymer chain will not have much impact on the whole system.
Best wishes,
Andy
Hello Andrew and Giacomo,
Now I run into another question. Will the NPT be sufficient to allow
flexibility of the polymer? For instance, my polymer is lying on the z
direction and I couple the x y pressures. Will my polymer behave like a
rubber band on z direction?
You can definitely run constant pressure (or constant stress) simulations
with an infinite periodic polymer.
Here is an NPT simulation I ran which was periodic in the Z direction:
The polymer was initially straight. Only the Z dimension was allowed to
change.
The box compressed in the Z direction simply due to the tension in the
polymer pulling the two ends together. (If you look carefully, near the
end of the movie, you can see that it is a periodic system because there
are two copies of the "red" region of the polymer.)
I question this because there should be solvents which are averaged to the
total z direction pressure and I am afraid that the tension from polymer
chain will not have much impact on the whole system.
Is this a polymer under tension? If so, I assume you have applied tension
to the polymer by using a lower pressure in the Z direction compared to the
X and Y directions. If these pressures are significantly different I
supposes you should not couple X,Y,Z together.
Otherwise the polymer might have a negligible effect on the pressure in the
z direction, since as you say, the pressure would likely be dominated by
the solvent. Still, since your system is anisotropic, I would not couple
X,Y,Z together. (Fortunately, it does not sound like you are doing that.
You are only coupling XY together.)
I don't know if this was helpful. (But it got me thinking about what I
should do to make my own simulations more realistic.)
Cool.
Andrew
Thanks for your information.
What I do not know is how will the polymer changes its structure. So I made a smaller box which contains more solvents than it suppose to have and let the x y expand aggressively. (Maybe I should freeze the polymer and do this in the energy minimization part to avoid a extensively blow up at startup.) Originally, x y are asymmetry so I come up with the idea to couple them. For the z direction, I will let the dynamic control the system to see how the polymer behaves, but my plan was to use the same pressure as x and y.Although I have not started it, those are the way I plan to construct my system. I shall try it and see if any problems come up.
Best wishes,
Andy
