Clarify how to compute ideal gas contribution to pressure

[g-bauer]

This issue was discussed in #58 (comment) and targets computation of the ideal gas pressure:

impl Compute for PressureAtTemperature {
    type Output = f64;
    fn compute(&self, system: &System) -> f64 {
        :
        let natoms = system.size() as f64;
        return natoms * K_BOLTZMANN * self.temperature / volume + virial / (3.0 * volume);
             //  ^--- natoms, degrees of freedom, number of molecules?
    }
}

[Luthaf]

So I am pretty sure this should be the number of degrees of freedom. It is what LAMMPS declare using in it's documentation, and it is what RASPA (another MD and MC capable software) uses for both MC and MD computation.

It also goes to the number of atoms in fully flexible systems, and to the number of molecules in rigid molecules system.

Finally, in the Allen and Tildesley book, they derive the instantaneous pressure using the instantaneous temperature (eq 2.55), and the instantaneous temperature uses the number of degrees of freedom (eq 2.50).

You can find the two corresponding pages from the book, containing the whole derivation.

[Luthaf]

Proposition to solve this: we can add a Molecule::degree_of_freedom() -> usize function, that is used to get the number of degree of freedom for each molecule, from 3 if the molecule is rigid to 3*natoms if the molecule is flexible, going through 3*natoms - n_constraints for constrained MD (RATTLE/SHAKE algorithms).

I think the main issue here is that we don't know whether a molecule is simulated a rigid or as flexible. We could start by assuming that MC means rigid and MD means flexible, but this would break quickly as we add more functionalities.

If/when we implement #137, we can have a nice way to resolve this. We can add a rigid key (true/false/[list of constraints]) in the molecule section, and use this in the Molecule::degree_of_freedom function.

What do you think?

[g-bauer]

I think the main issue here is that we don't know whether a molecule is simulated a rigid or as flexible.

I agree, we should think about how we can introduce constraints in a way that makes sense for both MC and MD as well as solves the degrees of freedom problem.

Often, MC codes introduce specific functions e.g. for constant bond lengths and angles. I think for a code that has both MC and MD, introducing method-exclusive functions makes little sense. I'd favor for example a constraints section in the input that takes inputs such as angles or/and bonds or rigid - that is what you also propose, from my understanding.

We should discuss this in a separate issue, since there are several aspects with regards to the implementation that are important for things like containing algorithms or configurational bias MC.

[Luthaf]

I think I have an idea for this, using roughly the same trick as we use for temperature computation: we could ask the propagator for the number of simulated degrees of freedom, and use this to compute the pressure. I'll try to send a PoC soon, so that we can discuss with the code =)