Ice Ih at 60K

[alberto santonocito]

Dear MBX Team,

I am currently performing PIMD simulations of Ice Ih using the MBX software with i-pi. I have attached my input and output files for a simulation run at 60 K (using 109 beads) for your reference.

I am analyzing the calculated thermodynamic quantities and would appreciate your feedback on the following points to ensure my interpretation of the MB-pol model is correct.

1. Kinetic Energy and Total Enthalpy:

I calculated an average Kinetic Energy of roughly 33.4 kJ/mol per molecule. Using a harmonic ZPE estimate (1/4 of the sum of intramolecular frequencies: 2 stretching + 1 bending), I approximated the intramolecular contribution to be roughly 28 kJ/mol per molecule. I attribute the remaining ~5.4 kJ/mol in my result to the kinetic contribution of intermolecular ZPE (librations).

Consequently, my Total Enthalpy (<K> + <U> + PV) is positive, as the kinetic energy (~33.4 kJ/mol per molecule) outweighs the potential energy (~ -25.7 kJ/mol per molecule). Could you confirm if this is the expected behavior for MB-pol, or is this anomalous?

2. Potential Energy Magnitude:

Regarding the potential energy, I obtained a value of roughly -25.7 kJ/mol per molecule. I assume this is higher than the experimental cohesive energy (~ -58 kJ/mol per molecule) because MB-pol is fully flexible, and the intramolecular strain adds a large positive term to the potential. Is this interpretation correct?

3. Heat Capacity Convergence:

 Finally, I aim to calculate the specific heat capacity (Cv). I am concerned that with 109 beads, the statistical noise in the virial estimator for total energy fluctuations might be too large to obtain a converged Cv. Would you recommend trusting the fluctuation method here, or is numerical differentiation of Energy vs Temperature the preferred approach for MB-pol?

Thank you very much for your time and assistance.

Best regards,

Alberto Santonocito