## Bug #375

### constraints = hangles is broken

**Description**

The "constraints=hangles" option is broken. The problem is that it generates too many constraints (more than the number of degrees of freedom available to be removed), leading to a singular constraint matrix. For example, consider the atoms

...-C-CH3

It generates a total of 10 distance constraints for those five atoms, one for each pair. Five atoms have 15 degrees of freedom, minus 10 constraints leaves only 5 residual degrees of freedom, which is too few. Nine constraints would be exactly right to make them move as a rigid body with six degrees of freedom. By adding one constraint too many, it makes the constraint matrix singular with the result that both LINCS and SHAKE fail to converge.

This can be fixed by reducing the number of constraints added by this option. There is no benefit to constraining **all** angles involving a hydrogen in any case. You really only care about ones of the form H-X-H or X-O-H. Those are the only ones needed to allow a 4-5 fs time step, and it avoids the problems described above while also not reducing flexibility as much.

The following code change implements this behavior. Starting at line 136 in topshake.c (the 4.0.5 version), it currently reads

`if ((nshake == eshALLANGLES) ||`

(count_hydrogens(info,3,ang->a) > 0)) {

Change that to

`int numhydrogens = count_hydrogens(info,3,ang->a);`

if ((nshake eshALLANGLES) ||

(numhydrogens > 1) ||

(numhydrogens 1 && toupper(**(info[ang->a[1]]))=='O')) {

### History

#### #1 Updated by Erik Lindahl over 10 years ago

Fixed+committed in master & release-4-0-patches.