Not yet reviewed by Pith; the record is open.

In this work, we develop energy-preserving iterative schemes for the (non-)linear systems arising in the Gauss integration of Poisson systems with quadratic Hamiltonian. Exploiting the relation between Gauss collocation integrators and diagonal Pad\'e approximations, we establish a Krylov-subspace iteration scheme based on a $Q$-Arnoldi process for linear systems that provides energy conservation not only at convergence --as standard iteration schemes do--, but also at the level of the individual iterates. It is competitive with GMRES in terms of accuracy and cost for a single iteration step and hence offers significant efficiency gains, when it comes to time integration of high-dimensional Poisson systems within given error tolerances. On top of the linear results, we consider non-linear Poisson systems and design non-linear solvers for the implicit midpoint rule (Gauss integrator of second order), using the fact that the associated Pad\'e approximation is a Cayley transformation. For the non-linear systems arising at each time step, we propose fixed-point and Newton-type iteration schemes that inherit the convergence order with comparable cost from their classical versions, but have energy-preserving iterates.