A method for finding integrable Hamiltonian systems with one and a half degrees of freedom is introduced and applied to the problem of finding vacuum magnetic fields with rotational transform and a dense set of toroidal flux surfaces. In this method, the fixed points of the return map must be located, and the values of the corresponding residues calculated. The residues, being indicators of stochasticity, are systematically reduced by variation of the parameters of the system. Application of this method to the problem of finding integrable magnetic fields shows that it is practical, and it yields substantial improvement of the configuration. In particular, this method can result in an increase of the simple equilibrium‐&bgr; estimate, &igr;2&egr;, by a factor of 2 or more.