I am working on implementing MD simulations with multiple thermostats to control the temperature of individual slabs. To achieve this, I divide the system into multiple slabs along the y-direction. This approach helps mitigate large local temperature deviations caused by latent heat release during free solidification of a pure metal. Below is the relevant portion of my script:
velocity all create ${temp} ${rand}
label loop
variable i loop 28
variable y1 equal -28+(${i}-1)*2
variable y2 equal -28+(${i}+0)*2
variable ss equal 1000+${i}
region reg${i} block -3 3 ${y1} ${y2} -3 3 units box
group grp${i} region reg${i}
fix t${i} grp${i} langevin ${temp} ${temp} 0.1 ${ss}
next i
jump SELF loop
############################################################
fix 2 all npt temp ${temp} ${temp} 0.1 y 0 0 1.0
run ${Nstep}
This script runs successfully. Is this script correct? Also, how to choose the best damping time for the Langevin thermostats?
Can you please be more specific? The script runs successfully, so what details of the script do you think are incorrect? My main concern is that the damping time may not be appropriate, which may affect the atom dynamics.
Or could I use multiple NPT thermostats for these slabs? But LAMMPS does not seem to support this.
The fact that an input runs to completion is by no means an indication that it is correct. It can be performing a completely unphysical simulation and still run to completion. It only means that you have no syntax errors.
You are thermostatting atoms multiple times. That is completely unphysical and bogus.
There is lots of literature about thermostats, including langevin thermostats. You would need to study those to learn what is appropriate for your circumstances. There is no global “this is the right choice”.
That said, you entire setup of applyting multiple thermostats to multiple layers makes no sense with a thermostat like fix langevin, since that is applied to individual atoms anyway. So just one thermostat would suffice.
That makes even less sense.
For a good reason, you cannot have multiple fix commands changing the same box dimension at the same time. That is completely unphysical.
Thank you for the detailed explanation and prompt response.
My goal is to simulate a free solidification process at a temperature slightly below the melting point. I already have a script using a global NPT thermostat, where the system size is free to change in the solidification direction:
fix 2 all npt temp ${temp} ${temp} 0.1 y 0 0 1.0
This script works well and produces meaningful results. However, this method is not good enough because solidification releases latent heat, which may raise the local temperature near the solid-liquid interface well above the target value. There are several methods to address this issue, and the most effective one I found in the literature (see below) involves applying multiple thermostats to control the temperature of small regions within the system.
J. Monk, et al. Determination of the crystal-melt interface kinetic coefficient from molecular dynamics simulations." Modelling and Simulation in Materials Science and Engineering 18, no. 1 (2009): 015004.
Is there a convenient way in LAMMPS to implement this multiple-thermostat strategy?
Applying the additional langevin thermostats is in my personal opinion turning your simulation of a real physical behavior into a computer animation. The release of the heat is a physical effect, but using a thermostat to suppress that has no physical equivalent.
Feel free to ignore my opinion, but don’t expect any further comments from me.
Consider @akohlmey first answer and have a look at your log file. It will contain warnings that hopefully will make clear(er) to you why your setup is bogus.
Consider reading fix langevin and NPT doc pages as well as the How-to thermostat and barostat pages. They’ll show you how to use langevin correctly.
The first concerns the physics. I agree that suppressing heat release is unphysical in principle. However, when the goal is to measure the kinetic coefficient—specifically, the slope of the solidification velocity versus temperature—it becomes important to accurately capture the local temperature near the solid–liquid interface. While using a global thermostat is more physically justified, the global temperature does not necessarily reflect the local conditions at the interface. Based on the methodology in the published paper I referenced, it seems reasonable to assume that we could suppress heat release without significantly affecting the atomic dynamics relevant to the kinetic coefficient.
The second layer is the technical challenge of implementing multiple thermostats in LAMMPS. I am not fixed on using the Langevin thermostat—that was just one approach I tested. I have tried several alternatives, but none worked. I initially thought this would be a straightforward task in LAMMPS, but it appears to be more complicated.
@Germain, regarding the log file: I couldn not find any warnings related to the Langevin thermostat in log.lammps; it appears to be working perfectly. My main concern was about selecting an appropriate damping time, which may affect the atom dynamics. I have also attached my LAMMPS script run.in for reference. I will certainly continue reviewing the documentation for more insights.