- berendsen_pscale(inputrec, pcoupl_mu, state->box, state->box_rel,
- start, homenr, state->x, md->cFREEZE, nrnb);
+ /* We should only scale after a step where the pressure (kinetic
+ * energy and virial) was calculated. This happens after the
+ * coordinate update, whereas the current routine is called before
+ * that, so we scale when step % nstpcouple = 1 instead of 0.
+ */
+ if (inputrec->nstpcouple == 1 || (step % inputrec->nstpcouple == 1))
+ {
+ berendsen_pscale(inputrec, pcoupl_mu, state->box, state->box_rel,
+ start, homenr, state->x, md->cFREEZE, nrnb);
+ }