def explicit_jacobian_solve(matvec, b):
    return lax.stop_gradient(jnp.linalg.solve(api.jacobian(matvec)(b), b))
def explicit_jacobian_solve(matvec, b):
    return lax.stop_gradient(jnp.linalg.solve(api.jacobian(matvec)(b), b))