NEMD Water Simulation

Hi,
I am trying to run a Non-Equilbrium molecular dynamic simulation which creates a steady-state temperature profile across the domain with the water being maintained at a constant pressure of 1.0 atm. I have been referencing the example in.spce.ehex located in the HEAT examples folder. How would one need to alter this code to achieve the desired results mentioned above?
Regards,
Peter

The examples/KAPPA folder has scripts which initialize and maintain a temperature profile
across a system, for the purpose of computing a thermal conductivity. Something similar
can be done for molecular systems. You don’t say what kind of water model you want
to use. Fix ehex is a variant of SHAKE for a rigid-body model.

Steve

Hi Steve,
Thank you for the response. Sorry to exclude the water model type. I would like to use the SPC/E water molecule. If I am understanding correctly, a similar approach seen in the KAPPA example can be applied to the HEAT SPC/E water molecule example to maintain a specific temperature profile across the domain. Would replacing the fix NVE with a fix NPH be the optimal approach to maintaining a constant pressure within the domain as well? I preferably would like the volume not to change during the simulation, is this possible when trying to control the pressure? Lastly, what recommendations would you make to reduce the oscillations of the temperature gradient in the KAPPA example? The gradient is relatively stable, oscillating about the set point, but reduced oscillations would be preferable for my simulation. Thanks in advance for your help.

Regards,
Peter

Hi Steve,
Thank you for the response. Sorry to exclude the water model type. I would like to use the SPC/E water molecule. If I am understanding correctly, a similar approach seen in the KAPPA example can be applied to the HEAT SPC/E water molecule example to maintain a specific temperature profile across the domain. Would replacing the fix NVE with a fix NPH be the optimal approach to maintaining a constant pressure within the domain as well?

I preferably would like the volume not to change during the simulation, is this possible when trying to control the pressure?

that is not possible, since the constant pressure integrator will change the volume in response to the pressure. if you want to keep the volume fixed, then you should gradually change the volume with the change_box command, and then equilibrate for a bit, until the average pressure is close (enough) to what you desire.

Lastly, what recommendations would you make to reduce the oscillations of the temperature gradient in the KAPPA example? The gradient is relatively stable, oscillating about the set point, but reduced oscillations would be preferable for my simulation.

temperature (and pressure) fluctuations depend on the system size. simulate a larger system and the fluctuations will be smaller.

axel.

Hi Axel,
Thank you for the email. My last follow up question is as follows,
Is it possible to fix the temperature in a region on either end of the domain and a heat flux in the middle to create a NEMD? I am currently defining three regions. One on either end which are the heat sources and one in the middle which is a heat sink with periodic boundary conditions over the entire domain. I am fixing the temperatures with a langevin thermostat in the source regions and the heat flux in the central region with ehex. I am running a fix nve on all of the molecules in order to integrate the equations of motion. Does this approach make sense?
Regards,
Peter

Hi Axel,
Thank you for the email. My last follow up question is as follows,
Is it possible to fix the temperature in a region on either end of the domain and a heat flux in the middle to create a NEMD? I am currently defining three regions. One on either end which are the heat sources and one in the middle which is a heat sink with periodic boundary conditions over the entire domain. I am fixing the temperatures with a langevin thermostat in the source regions and the heat flux in the central region with ehex. I am running a fix nve on all of the molecules in order to integrate the equations of motion. Does this approach make sense?

typically, you would have 4 regions. e.g. one on the side that is a heat source/sink, one in the middle that is a heat sink/source and then then two in-between regions with a temperature gradient (one up, the other down). otherwise, you would have direct heat transfer from the heat source to the sink due to periodic boundaries.
please see the files examples/HEAT/in.spce.hex and examples/HEAT/in.spce.ehex and their

axel.