I am currently using LAMMPS (version 8 Feb 2023) to perform shear deformation simulations on an amorphous system. In my setup, the simulation box is sheared up to 61 times its initial length. I am using the NPT ensemble along with shear deformation at a strain rate of 1e10/s, with the following commands:
fix npt all npt temp 77 77 100 iso 1 1 1000
fix deform all deform 1 yz final ${final_strain} remap x flip yes
Here, ${final_strain} is set to 61 times the initial length of the simulation cell.
While the density of the system remains steady, I observe a significant increase in the amplitude of stress fluctuations (highlighted in the red rectangle). This issue does not appear in the early stages of the simulation (green rectangle).
I would appreciate any insights into why these increasing stress fluctuations are occurring during the later part of the shear process.
To investigate this further, I tried the following:
I restarted the simulation from a point of high fluctuation using a restart file. Initially, the amplitude of the stress fluctuations returns to normal, but it builds up again after some time.
I took dump trajectories (containing only xyz values) from the fluctuating region and used the ārerun commandā on them using LAMMPS. The increase in stress fluctuation amplitude persisted, suggesting that the phenomenon is already imprinted in the trajectory data itself.
It is unlikely you get any relevant help if you do not provide more relevant information on the kind of system you attempt to simulate, what stress you attempt to compute and how you implemented your simulation. At the moment there is very little we can say that can be relevant to your issue and to determine if this is expected behavior or not.
Also is there a reason to use such an old version of the code?
Hi @Germain
Thank you so much for your reply. As noted earlier, Iām running a shear-deformation simulation on an amorphous system contained in a triclinic box with periodic boundaries in all directions. The simulation uses real units with a 5 fs timestep. Please let me know if you need any further details about the setup.
To impose a shear rate of 1 Ć 10Ā¹ā° sā»Ā¹, I first calculated ${final_strain} and then set the total number of timesteps accordingly. The simulation employs two fixes simultaneously:
fix npt all npt temp 77 77 100 iso 1 1 1000
fix deform all deform 1 yz final ${final_strain} remap x flip yes
The NPT fix keeps the pressure at 1 atm, after which the fix deform command shears the system along the yz plane.
In the stress plot I shared earlier, the y-axis shows P_yz (pressure component in yz), while the x-axis represents the timestep.
Iām running these calculations on our universityās HPC cluster, which is why Iām using an older LAMMPS version.
There is no need to repeat yourself if I mention that this is not enough information. What I meant is about the physical system you want to simulate and to which extent (the size of your box). A metallic glass does not behave the same as bulk polymer, yet both are amorphous. Amorphous is a state of matter. Many things can be amorphous.
Depending on your model, this can be quite a high integration time step. But without more information that is hard to tell with certainty.
What I mean by implementation is which commands you used in your script for the whole setup. Sometimes, errors lie in commands we donāt think about first. But anyway, is there any reason to use both NPT and deform? Why not equilibrate your simulation in NPT and stick to the equilibration volume when shearing? Barostats have to react to the instantaneous pressure of your system to let it relax and this takes some time. Note that with a very high shear rate as yours, the shear stress might not have time to relax at all. Also note that with your increase timestep, your Tdamp and Pdamp constants are quite low.
As with many HPC clusters, Iām pretty sure there are ways to download and compile your own version of LAMMPS on your user session, but you do you as the saying goes.
Iām simulating an amorphous glass with a coarse-grained water model and intend to shear the simulation box until it is 61 Ć its original length.
The reason is that the systemās density changes as it deforms, so a fixed-volume approach may not capture that behavior.
Would combining fix nvt/sllod with fix deform/pressure let me maintain the target pressure without artificially locking the volume? Iām not sure whether that setup can satisfy both requirements.
See, this is what I meant when asking for details. Water is a very particular kind of simulation on its own and coarse-grain models add a layer of complexity. Hence, this makes things difficult to discuss in a meaningful way. It is nearly impossible to comment on the correctness of the behavior youāre seeing if we have no idea about the physics we have to expect.
What you should do is refer to the literature about your model and the property you are trying to simulate and try to understand what is going on in your simulation. Have you looked at your simulationās trajectory? How stable are you structures? What happens when the glass dislocates? How is your box changing? Then it is up to you to interpret what you are seeing in physics terms and link it to the phenomenon you want to model. If you notice something incoherent, then you can come and explain why, and we will help if possible.
But for the moment it seems that you lack a detailed analyses on what is going on compared to what you were trying to achieve.
Note also that SLLOD equations are for liquid flows. Not for glass and solid deformations.
Not sure if that is going to sit well with LAMMPSā domain decomposition algorithm. It will be horribly inefficient for sure. Also, unless you are simulating a particularly large system, you will get some very thin slices of a cell, where atoms will interact with themselves through periodic boundaries and that will make your results very questionable.
As @Germain already stated, It is essentially impossible to make an accurate assessment of the situation without being involved in the research yourself. For sure, this is a problem not caused by LAMMPS directly, but more likely through the model and methods that are being implemented (and thus would likely happen the same way with other MD codes). Thus, I can only speculate. Two things would come to my mind, a) the unusual geometry and unexpected finite size effect due to the extreme shearing; and b) a bad choice of Nose-Hoover damping parameters for the barostat: you may have chosen a parameter that is close to resonant with the fictitious dynamics of the Nose-Hoover chain.
You should not be using a barostat with such a tricky system setup unless you know for sure that you expect your real-life system to respond to pressure changes on picosecond timescales. The more things you add to a non-equilibrium simulation, the more likely it will be that you mis-specify something or have some components interact poorly.
I notice you did not use SLLOD ā you should. As the box shears particles will pick up some amount of flow velocity which must be accounted for during thermostatting (otherwise, if this velocity is counted as thermal energy, the particles will be overcooled).
As far as I know, SLLOD equations of motion assume homogeneous planar shear flow. I donāt think this will be the case of sheared glass system even containing water. But this is just a guess given the few details we have on the simulated system and the simulation procedure.
Youāre right ā both that SLLOD most restrictively refers to laminar shear flow conditions and that there is not enough info for us to help constructively
