Initial velocities for MD non-restart simulations

[trdurrant]

Question

For dynamics there is the following code (which will end up being called by md.run():

    def _run_dynamics(self) -> None:
        """Run dynamics and/or temperature ramp."""
        [...]
        # Set velocities to match current temperature
        if not self.restart:
            self._set_velocity_distribution()

i.e. if the calculation is not a restart, set a Maxwell-Boltzmann distribution of velocities at self.temp

But this seems unexpected to me, as this means that user supplied velocities are discarded. This has a few effects:

• If I start a new MD run from the positions and velocities of an old run, the velocities will be re-initialised by this line
• For NVE calculations, it is not possible to for a user to set their own initial distribution of velocities, (if no temperature is directly set, then the default temperature of 300 is used)

In contrast, in the initialiser for the NVT_CSVR class, there is the following code:

        # Velocity distribution must be non-zero before dynamics is set
        if np.isclose(self.struct.get_kinetic_energy(), 0.0, rtol=0, atol=1e-12):
            if self.logger:
                self.logger.warning("Velocities modified during MD initialisation")
            self._set_velocity_distribution()

This is more in line with what I would have expected, but still in this case the velocity distribution set here will be replaced by the first code block, regardless of if this code triggers or not.

I would like to run NVE with a custom velocity distribution, but I don't see how to do this at the moment. Does it make sense to remove the initial block of code, and instead try and include something like the second block in each of the MD classes? Or perhaps I haven't understood this code correctly.

Platform

No response

Version

No response

[alinelena]

hi @trdurrant,

you are right if the user provides velocities they shall be used... @oerc0122 @ElliottKasoar @harveydevereux thoughts?

[ElliottKasoar]

Good question. Perhaps it should just be an option (probably defaulting to True to preserve the current behaviour) that (re)sets the initial velocity distribution? So if set to False, we use the velocities read in

[trdurrant]

I don't have a problem with doing it that way, although if the velocity distribution is non-zero and is being reset, I do think it should raise a warning (as in the second snippet)

[ElliottKasoar]

That warning is primarily because it's unusual that we modify the atoms on initialisation, rather than run, which is particularly surprising.

If we have an explicit option, hopefully this behaviour isn't too surprising, but I'd definitely happily include a logging.info note that it's happened.