GROMACS NVT and NPT Equilibration: A Step-by-Step Guide for Beginners

Equilibration prepares a solvated system for production MD in two phases. NVT runs at constant volume with a thermostat until the temperature settles at your target (often 300 K). NPT then runs with a barostat until pressure and density stabilise. Keep the solute fixed with position restraints during both, so the solvent relaxes around a structure that has not yet moved.

You minimised the energy, the potential dropped, and you are ready to hit run. Do not. Starting production molecular dynamics straight after energy minimisation is the fastest way to watch a simulation blow up on its first picosecond. The solvent is still packed at whatever density the setup tool guessed, the atoms have no sensible velocities, and the pressure is nowhere near your target. Equilibration fixes all three, in a fixed order, before a single frame of real dynamics is recorded. This guide from the StemSkills Lab team, who have spent 10+ years running production MD across drug-discovery and materials projects, walks through both equilibration phases with the exact commands and settings, then shows you how to prove each one worked. New to the engine itself? Start with our GROMACS molecular dynamics guide and return here when you reach the equilibration step.

What is equilibration in a GROMACS simulation, and why are there two phases?

Equilibration is the bridge between a static starting structure and a stable production run. It relaxes the system to the temperature and pressure you actually want to study, so that production MD samples a realistic equilibrium ensemble instead of the artefacts of your setup. GROMACS splits this into two ensembles because temperature and pressure are best settled one at a time:

• NVT (constant Number of particles, Volume, Temperature): also called the isothermal-isochoric or canonical ensemble. The box size is frozen while a thermostat brings the system to your target temperature and gives the atoms a proper thermal velocity distribution.
• NPT (constant Number of particles, Pressure, Temperature): the isothermal-isobaric ensemble. Now a barostat is switched on. It lets the box volume adjust until the pressure fluctuates around your reference value and the density reaches its correct equilibrium.

Running NVT first is deliberate. You want the temperature stable before you allow the volume to change, otherwise the barostat is chasing a moving target and the box can swing wildly. GROMACS is the open-source MD engine described in Abraham et al., “GROMACS: High performance molecular simulations through multi-level parallelism from laptops to supercomputers,” SoftwareX 1-2 (2015): 19-25 (doi:10.1016/j.softx.2015.06.001), and this two-phase equilibration is the standard workflow taught in its reference tutorials.

Why do you keep position restraints on during equilibration?

During equilibration you want the water (and ions) to reorganise around your solute, but you do not yet want the solute itself to move, because its starting coordinates came from a crystal structure or a model you trust. Position restraints solve this by attaching a stiff harmonic spring to each heavy atom of the solute, holding it near its starting position while everything else relaxes.

In GROMACS you activate them with a preprocessor flag in the .mdp file:

define = -DPOSRES

This tells gmx grompp to include the posre.itp file that gmx pdb2gmx generated when you built the topology. You also pass a reference structure to grompp with the -r flag, which is the coordinate set the restraints pull each atom back toward. Skip either the -DPOSRES define or the -r file and the restraints silently do nothing, which is one of the quieter ways an equilibration goes wrong.

How do you run the NVT equilibration phase in GROMACS?

NVT takes your minimised structure and, using an nvt.mdp parameter file, brings it to temperature. Assemble the run input file with grompp, then run it:

gmx grompp -f nvt.mdp -c em.gro -r em.gro -p topol.top -o nvt.tpr
gmx mdrun -deffnm nvt

Here -c em.gro is your energy-minimised structure, -r em.gro supplies the position-restraint reference, and -deffnm nvt tells mdrun to name every output file nvt.*. The settings that matter inside nvt.mdp are:

• gen_vel = yes, gen_temp = 300, gen_seed = -1: generate starting velocities. The GROMACS .mdp options reference states these are drawn “according to a Maxwell distribution” at gen_temp, which is why NVT is where velocities are born.
• tcoupl = V-rescale: the velocity-rescaling thermostat. Set tc-grps = Protein Non-Protein, tau_t = 0.1 0.1, and ref_t = 300 300 to couple solute and solvent to a 300 K bath with a 0.1 ps coupling time.
• pcoupl = no: no barostat yet. The volume stays fixed, which is what makes this the NVT ensemble.
• continuation = no: this is the first dynamics after minimisation, so constraints are generated fresh rather than continued.

A typical NVT phase is 100 ps (50000 steps at a 2 fs timestep). The V-rescale thermostat comes from Bussi, Donadio, and Parrinello, “Canonical sampling through velocity rescaling,” Journal of Chemical Physics 126, 014101 (2007) (doi:10.1063/1.2408420); it is the recommended thermostat because, unlike the older Berendsen scheme, it produces a correct canonical ensemble.

How do you run the NPT equilibration phase?

NPT starts from the end of NVT and adds a barostat so the density can settle. The key difference is that you now carry the velocities and state forward from NVT with a checkpoint:

gmx grompp -f npt.mdp -c nvt.gro -r nvt.gro -t nvt.cpt -p topol.top -o npt.tpr
gmx mdrun -deffnm npt

The -t nvt.cpt flag hands NPT the checkpoint from NVT, so the run continues with the thermal velocities you just equilibrated instead of generating new ones. Inside npt.mdp, the changes from the NVT file are:

• gen_vel = no and continuation = yes: you are continuing, not restarting, so keep the existing velocities and constraints.
• pcoupl = C-rescale: switch on the barostat. Use pcoupltype = isotropic, tau_p = 2.0, ref_p = 1.0 (1 bar), and compressibility = 4.5e-5 for water.

The compressibility of 4.5e-5 bar^-1 is the standard value for liquid water at ambient conditions and tells the barostat how strongly volume should respond to a pressure difference. C-rescale is the modern stochastic barostat; the classic Lemkul lysozyme tutorial uses Parrinello-Rahman here instead, which also works but oscillates badly if you start it far from equilibrium, so many practitioners prefer C-rescale (or Berendsen) for the equilibration step and save Parrinello-Rahman for production. NPT is usually run for another 100 ps or longer.

How do you confirm each equilibration phase actually converged?

Never assume equilibration worked because mdrun exited without an error. Extract the relevant quantity from the energy file with gmx energy and look at the plot. For NVT, check the temperature:

gmx energy -f nvt.edr -o temperature.xvg

At the prompt, type the number next to Temperature, then 0 to finish. A converged NVT run shows temperature climbing from its initial value and then fluctuating tightly around 300 K, with no upward or downward drift. For NPT, check pressure and density from npt.edr the same way:

gmx energy -f npt.edr -o pressure.xvg
gmx energy -f npt.edr -o density.xvg

Pressure in an MD box is genuinely noisy, so do not expect a flat line; you want the running average to sit near 1 bar. Density is the cleaner signal: it should rise or fall from its starting value and then plateau. For a system solvated in a common water model the plateau lands near 1000 kg/m³ (SPC/E water equilibrates close to 1008 kg/m³). A stable density plateau is the single clearest sign your system is ready for production. Building this habit of verifying every stage is exactly the kind of workflow discipline we sequence in the computational biology skills roadmap.

NVT versus NPT: what actually differs?

How do you fix the classic equilibration failures?

Most equilibration crashes come down to a handful of causes. Here is what each looks like and the fix.

• LINCS warnings and “system is exploding”. Atoms move impossibly far in one step, usually because energy minimisation was skipped or did not converge, the starting structure has clashes, or the timestep is too large before the system is relaxed. Fix: confirm minimisation brought the maximum force below your emtol, keep position restraints on, and if it persists, run the first equilibration with a smaller dt (for example 1 fs) before moving to 2 fs.
• “Water molecule starting at … can not be settled.” A SETTLE constraint failure, the same class of blow-up as above, pointing at a badly placed or overlapping water. Fix: re-run energy minimisation, and check that solvation and ion addition did not leave an atom on top of another.
• The box keeps shrinking or density is far off. Almost always a barostat problem, often Parrinello-Rahman applied to a system that was not yet close to equilibrium. Fix: equilibrate the pressure with C-rescale or Berendsen first, and confirm your compressibility and ref_p are physically sensible.
• Temperature never reaches the target. Usually the thermostat coupling is wrong. Fix: couple solute and solvent as separate groups (tc-grps = Protein Non-Protein) rather than the whole system as one, and check ref_t matches your intended temperature on every group.

Get equilibration right and production MD starts from a properly relaxed system instead of crashing on the first nanosecond. Get it wrong and no amount of production time will rescue the trajectory. That is why we treat it as a skill in its own right, not a box to tick.

Frequently asked questions

Do I always need both NVT and NPT before production?

For a solvated biomolecular system, yes. NVT gives the atoms a correct velocity distribution at your target temperature, and NPT lets the density relax to the right value at your target pressure. Skipping NPT leaves you at the setup tool’s guessed density, which is rarely correct.

How long should each equilibration phase run?

The common starting point is 100 ps for NVT and 100 ps for NPT, but the real answer is “until the property converges.” Watch temperature for NVT and density for NPT; if either is still drifting at the end, extend that phase rather than moving on.

What is the difference between the V-rescale thermostat and the C-rescale barostat?

V-rescale is a thermostat: it controls temperature by rescaling velocities and is used in both NVT and NPT. C-rescale is a barostat: it controls pressure by rescaling the box volume and is only switched on in NPT. They do different jobs and are often used together.

Why did my restraints not seem to hold the protein in place?

Two settings must both be present. The .mdp file needs define = -DPOSRES so grompp includes posre.itp, and the grompp command needs -r pointing at a reference structure. Miss either and the restraints are inactive even though the run completes.

Can I combine NVT and NPT into a single equilibration run?

You can, but separating them is clearer and easier to debug. Running them apart lets you confirm the temperature is stable before you let the volume change, and it gives you two clean plots to verify instead of one tangled one.