We discuss a variety of factors determining the behavior of a combined Fourier analysis plus Alternating Direction Implicit (ADI) implementation to solve Poisson’s equation in cylindrical coordinates. The boundary value of the potential is found through integration of Green’s functions over the density distribution. We have varied (1) the distance to the grid boundary, (2) the number of terms used in the Green’s functions and (3) the number of ADI iterations in simulations of a detached binary initially on a circular orbit. Numerical errors resulting from the Poisson solver are reflected in epicyclic and apsidal motion of the binary and cumulative numerical errors result in motion of the system center of mass over long term evolutions. In light of constraints on computing resources, we explore optimal parameters for the Poisson solver for long term, direct numerical simulations.