The main result of this paper is the discretization of Hamiltonian systems of the form $\ddot x = -K \nabla W(x)$, where $K$ is a constant symmetric matrix and $W\colon\mathbb{R}^n\to \mathbb{R}$ is a polynomial of degree $d\le 4$ in any number of variables $n$. The discretization uses the method of polarization and preserves both the energy and the invariant measure of the differential equation, as well as the dimension of the phase space. This generalises earlier work for discretizations of first order systems with $d=3$, and of second order systems with $d=4$ and $n=1$.