Add wrapper to linear solve.

[gileshd]

In cases where the matrix A should be is positive semidefinite, replace jnp.linalg.solve with wrapper in utils/utils.py so that we have control over the solver used.

For the time being the solver itself is not changed.

Note that the wrapper is only used is scenarios where I was reasonably confident that the matrix A is psd.

[slinderman]

This test failure happens from time to time and is tracked in an open issue. Merging anyway.