8.3 Validate NVE timestep and energy conservation
These are unexecuted teaching inputs and starting settings to test. Original figures are schematics, not computed results. Use licensed VASP and PAW data, replace every placeholder, record the executable version and validate convergence.
8.3.1 Model, units and provenance
Use eV for energy, Å for length and eV/Å for force; 1 kbar = 0.1 GPa. State normalization per atom, molecule, primitive cell or simulation cell. Record PAW identifiers, release, ZVAL, ENMAX and permitted hashes; never redistribute POTCAR. SCF convergence addresses the chosen electronic problem; convergence of the target property requires separate tests.
Original schematic. Curves explain concepts; blank data areas await verified learner results. No calculation is claimed.
8.3.2 Worked case: procedure, interpretation and checks
Intuition
In a fixed-cell, isolated classical system, potential and kinetic energy continually exchange. The sum should not steadily gain or lose energy simply because the simulation clock advances. A thermostat can remove numerical heating and make an unreliable integrator look healthy. An NVE test therefore acts like removing the training wheels: it exposes timestep and force-quality problems. It does not prove that the functional describes the material accurately.
Prerequisites. Start from an equilibrated configuration with its velocities, validated masses and fixed cell. Keep the electronic Hamiltonian and smearing unchanged across all comparison branches. For metals at finite smearing, use the energy consistent with the force formulation and read VASP’s energy definitions rather than substituting E0 blindly. Do not test an energy-conservation claim while deliberately heating the system or changing the model.
A simple thermostat-free route is Andersen dynamics with zero collision probability. Use it as an isolated branch, removing contradictory thermostat tags. MDALGO ensemble choices
A controlled numerical experiment
Create a matrix of branches rather than editing one directory repeatedly. All begin from the same verified physical state. For example, compare timesteps 1.0, 0.5 and 0.25 fs, if appropriate to your system; these are test choices, not approved timesteps. Use the same physical duration, so smaller Δt requires more steps. For hydrogen-containing stiff modes, begin conservatively; for high-temperature collisions, retest at the highest intended temperature.
Then repeat a chosen Δt with a tighter SCF tolerance. EDIFF controls electronic energy convergence and is not itself a direct bound on force error; a smaller number helps only if SCF actually reaches it. Compare selected snapshot forces against a tighter reference. EDIFF
For each branch, calculate
\(\delta e(t)=\frac{E_{\mathrm{conserved}}(t)-E_{\mathrm{conserved}}(0)}{N_{\mathrm{atoms}}}\)
and fit \(\delta e(t)=a+bt\) over a documented interval. If energy is in eV, time in ps, and divided by atom count, b has units eV atom⁻¹ ps⁻¹. Plot the full trace as well as the slope. A bounded oscillation can be compatible with a stable finite-step integrator; a monotonic trend is a different symptom. A single slope can hide a jump or alternating transients. Report drift per atom, fluctuation range, maximum step jump and the SCF-failure count.
Use the scale of the scientific observable to set an acceptance criterion. A small absolute drift may still matter for a delicate transition or long trajectory. Do not advertise one drift threshold as correct for every compound and timescale. Compare whether halving the timestep or tightening electronic settings changes the quantity eventually reported, not merely whether the movie continues.
Diagnostic decision tree
- Drift improves markedly when Δt is halved: integration resolution was important; assess the smaller step at target temperature.
- Drift improves mainly when SCF is tightened: electronic force noise or convergence was important; inspect difficult snapshots.
- A single energy jump coincides with a restart: suspect lost velocities, changed state, changed settings, or incompatible restart information.
- The first step is abnormal: examine initialization, geometry and how the state was transferred before judging the entire integrator.
- T drifts while potential energy changes during structural relaxation: the system may be converting latent/configurational energy, so inspect the conserved sum and phase evolution together.
- Neither change helps: investigate the energy definition, force consistency, precision, projector approximations, constraints and file parsing.
Novice traps. Fitting only potential energy; including thermostat energy from a different ensemble without understanding it; comparing equal step counts rather than equal time; leaving a friction or velocity-rescaling mechanism active; accepting a run because its final energy happens to equal its initial energy; forcing identical trajectories as a reproducibility criterion for chaotic dynamics.
Exercise. Build a 3×2 test plan for timestep and SCF tolerance. Before running, write how you will identify every energy field and set the physical duration. After running, rank settings by cost and drift. Explain why a faster run with poor conservation is not a valid optimization.
8.3.3 Unexecuted inputs and analysis scaffolds
These are unexecuted teaching inputs and starting settings to test. Original figures are schematics, not computed results. Use licensed VASP and PAW data, replace every placeholder, record the executable version and validate convergence.
8.3.3.1 Input block 1
# UNEXECUTED NVE delta; continue positions and velocities appropriately
IBRION = 0
MDALGO = 1
ANDERSEN_PROB = 0.0
ISIF = 2
ISYM = 0
POTIM = 0.5
NSW = 2000
NBLOCK = 1
# No active Langevin, Nose or temperature-ramp settings in this branch.
8.3.4 Related learning paths
- 1.1 Four input files, one physical question
- 1.2 A convergence laboratory with an error budget
- 8.2 Prepare an AIMD calculation and decide when equilibration ends
- 8.4 Choose an NVT thermostat without confusing sampling and dynamics