Newtonian Gravity on a Riemannian $3$-Manifold To solve the Poisson equation for the Newton Potential, say $\phi$, one can use the divergence theorem, such that
$$\int_U \nabla^2 \phi \sqrt{g}~ \mathrm dV= \int_{\partial U} \langle \nabla \phi,n\rangle  \sqrt{g_\Sigma}~ \mathrm d\Sigma=1,$$
where $\sqrt{g_\Sigma}$ represents the induced metric on the border.
Considering that $\phi$ is classic it will only depend on the relative position between the source and the point where we want to compute the potential, so, in order to use the divergence theorem, we need to find the set of points that are at the same distance of a given point (where we consider the source of the field to be) and that will be the border, $\partial U$, of the considered volume that contains the source.
Well, if space is the three dimensional Euclidean space then the volume that contains the source is straightforwardly the 3-ball whose center is at the source position and the border $\partial U$ the corresponding 2-sphere.
But if space is a different 3 Riemannian Manifold things are a "little" bit harder and finding the set of points that are at the same (geodesic) distance from the source is in general quite difficult or even impossible... But when searching for the geodesic curves on a manifold one can, in some circumstances, instead of searching for the curve that minimizes the distance
$$d~=~\int_a^b \left[g\left(\dot{c}(s),\dot{c}(s)\right)\right]^{1/2}~\mathrm ds,$$
search for the curve that minimizes the energy functional
$$E~=~\int_a^b \left[g\left(\dot{c}(s),\dot{c}(s)\right)\right] ~\mathrm ds,$$
since in those circunstances the curves that minimize the energy are the same as the curves that minimize the distance.
Well, finding the set of points that are at the same energy is much easier...I would then ask: what is the significance of the energy defined above and would it make sense, physically, to say that equipotential surfaces are not the set of points that at the same distance from the source but the set of points that have the same energy.
PS: Note that in the Euclidean Space, the set of points that are at the same distance is the same as the set of points that are at the same energy.
 A: You can use the energy functional to compute geodesics in any circumstances (even in GR!). The geodesic computed using the energy functional "differs" from a geodesic computed from the distance functional only in that it's automatically parametrised so that $g(\dot{c}(s),\dot{c}(s))$ is constant along the geodesic (this is automatic by conservation of energy if you consider the "energy" functional as the action for a particle). This last fact is very useful if you actually attempt to perform such a calculation by hand.
I do not think that trying to find points that are equidistant from some given point is going to help you solve the Poisson equation though; maybe if you are looking at a symmetric space you can compute the Newton kernel using the same trick physicists usually use to solve the Poisson equation for a point mass source, but it does not seem as that would help in full generality. I'd like to hear about any progress you make though.
EDIT: I forgot to specify how to compute the distance function... once you've solved the geodesic equation for curves joining some distinguished point $p^*$ to an arbitrary point $x$, the distance functional effectively becomes a distance function (of $x$) on the manifold so, in principle, you've solved the problem.
