# Bond Constraints with SHAKE and V-rescale thermostat

**URL:** <https://gromacs.bioexcel.eu/t/bond-constraints-with-shake-and-v-rescale-thermostat/9167>\
**Category:** User discussions\
**Created:** [May 17, 2024, 2:54pm UTC](https://gromacs.bioexcel.eu/t/bond-constraints-with-shake-and-v-rescale-thermostat/9167 "2024-05-17T14:54:17Z")\
**Posts on this page:** 2\
**Page:** 1

<div class="post-metadata">

**Author:** ![delysermr](https://avatars.discourse-cdn.com/v4/letter/d/47e85d/32.png) [@delysermr](https://gromacs.bioexcel.eu/u/delysermr)\
**Post date:** [May 17, 2024, 2:54pm UTC](https://gromacs.bioexcel.eu/t/bond-constraints-with-shake-and-v-rescale-thermostat/9167/1 "2024-05-17T14:54:17Z")

</div>

GROMACS version: 2019.4  
GROMACS modification: No

I am running a large CG simulation with a martini-like force field. Some of the bonds are constrained to a set distance, others are allowed to fluctuate. I am running in the NPT ensemble at 1 bar (Parrinello-Rahman barostat) and 300K. Time step is 10 fs, thermostat coupling constant is 2 ps, barostat coupling constant is 5 ps.

When I use the v-rescale thermostat and SHAKE to constrain the bonds, the average temperature comes out to around 275K. Further, when I use gmx energy to output the temperature, the output is:

Energy Average Err.Est. RMSD Tot-Drift  
'-------------------------------------------------------------------------------  
Temperature 5519.72 1600 31223 11013.5 (K)

yet the temperature.xvg file created by gmx energy gives values around 275K for the temperature at each output step.

I tried using Nose-Hoover & SHAKE and I get the right temperature. I tried using v-rescale & LINCS and I get the right temperature. I tried changing the constrained bonds to harmonic bonds and I get the right temperature.

Is/was there a known issue when combining the SHAKE algorithm with the v-rescale thermostat?

Thanks,

Michael DeLyser

---

<div class="post-metadata">

**Author:** ![hess](https://dub1.discourse-cdn.com/flex017/user_avatar/gromacs.bioexcel.eu/hess/32/416_2.png) [@hess](https://gromacs.bioexcel.eu/u/hess)\
**Post date:** [May 17, 2024, 3:01pm UTC](https://gromacs.bioexcel.eu/t/bond-constraints-with-shake-and-v-rescale-thermostat/9167/2 "2024-05-17T15:01:37Z")

</div>

My guess is that your SHAKE parameters are not sufficiently accurate, or maybe cannot be sufficiently accurate, and the constraint errors pump energy into the system faster than the v-rescale thermostat can pump out with a coupling time of 2 fs. If my guess is correct, you should not use SHAKE, or at least not with the SHAKE parameters you are using.
