I’ve been doing a bulk liquid simulation near a silica surface. The only important information is that the liquid is hydrogen bonded to the surface and creates a silonal group. The hydrogen is bonded to the silica and I’m using Fix shake to contained the bond length. I have noticed that due to the frozen degree of freedom the temp calculation is off. This problem can be fixed by subtracting the frozen degree of freedom. My question is to know how similar is the time integration in using Fix shake or I just simply increase the spring constant so that essentially the bond length is almost frozen. I do know that Fix shake allows the use of larger time step and that is not my concern at the moment. My concern is about the temp of these hydrogen atoms and their velocity which is affected by the frozen degree of freedom.
I’m suing June 2022 version of Lammps on Ubuntu.
The only way to know is to set up tests and quantify this yourself.
If you don’t care about the simulation time, why don’t you just apply a realistic (i.e. properly parameterized) force constant? The point of using a stiff force constant is that you would not have your geometry change too much from the constrainted lenghts, e.g. during a minimization where SHAKE cannot be applied… and if the simulation time should eventually matter, you can always consider using r-RESPA to reduce the cost of force computation for those components of the system where the forces change more slowly.
You may want to check if your temperature calculation is really off because of fix shake or because of some other issue in your simulation. Most temperature computes in LAMMPS “know” how to account for frozen DoFs from SHAKE. (Of course, some temperature computes don’t, and if you are using one of those and know it you can safely ignore this post.)
Stree, you’re bringing up a good point. Thanks.
I’ve been doing some try and error with using multiple thermostat versus using a single thermostat, and the results of temperature at equilibrium is different. Does this mean that there is a problem with my simulation or in fact, using multiple thermostat is a common practice when it’s necessary to enforce temperature convergence of different type of atoms?
Achieving equipartitioning (i.e. each subsystem has the same temperature) is part of reaching equilibrium. If you don’t reach this, that may have some reason you need to understand. This can be due to bad coupling between the subsystems.This would be highly dependent on the specific setup of your simulation. Enforcing equipartitioning with multiple thermostats can be a way to reach equlibrium, but it is not guaranteed. You also have to check, that this property is preserved when you turn off thermostatting or switch to a Nose-Hoover thermostat.
For a slab configuration with a bulk? (more likely a thin film?) liquid in contact, the situation can be quite complex depending on how you set up and equilibrate the slab and the liquid. There has been numerous discussions of this subject in the past, so you can try digging through the archives or study the materials and methods sections of relevant publications.
If there is a problem with the slab, (most likely a geometrical problem since fix shake diverges) then why this issue is hidden when I apply a thermostat for the bulk and surface? Whereas when I separate the surface thermostat and leave the rest (bulk) under another thermostat, then fix shake diverges.
I guess the important question is why Lammps does not detect the problem when a thermostat is applied versus separating the thermostats?
Lots of issues can be hidden with changing how you thermalize your system. The question is whether your choices make sense from the perspective of generating a physically meaningful system and trajectory. The fact that you have a simulation that does not crash is in no way proof of correct physics.
Then probably some other settings are not suitable.
It is not the job of LAMMPS to tell you whether your simulation is meaningful, it is yours to reason this out.
For starters, it should be possible to run without a thermostat at all (after guiding the system to equilibrium, that is). Good energy conservation without a thermostat is often a helpful quality measure about the reliability of a simulation setup.
Combinations of thermostats are not the same as a single thermostat.
Work through the math of this. This is not a LAMMPS issue, it is a science issue.
One simple (but significant) example: if a (deterministic) thermostat is applied to all particles in a simulation, and they all obey Newton’s Third Law, and the total momentum of the system is initially zero, then it is conserved at zero. If separate thermostats are applied to separate parts of a system, total system momentum is generally not conserved.