Performance and stability improvements to wcs.all_world2pix

Significant (typically 1000-2000 times faster) improvements of the
performance and stability of `wcs.all_world2pix` function. The used
algorithm is a completely vectorized implementation of the method of
consecutive approximations. While this method is a low-order method, it
is more than adequate for image distortions that are typical in
astrophysical images (smooth, slowly varying polynomials). The method
is somewhat similar to the existing `HSTWCS.all_sky2pix`, though the
choice of the function, initial approximation, and implementation,
makes the method used here to converge in half of the iterations of
`HSTWCS.all_sky2pix`.

This implementation no longer requires users to have SciPy installed in order to run.

While there is a check to see if WCS contains higher-order distortions
(and return the exact solution if it does not), the method itself does
not *require* this check and will converge to the exact solution after
one iteration. This check is mostly for performance purposes to avoid
unnecessary iterative solution for images that do not contain
distortions (and pix->sky requires only WCS-based, i.e., CD-only
transformations).

The proposed replacement for the current implementation of
`all_world2pix` avoids unnecessary use of `wcs_world2pix` and
`all_pix2world` such as `wcs_world2pix(all_pix2world(pix))` by
replacing this with `pix2foc`. This avoids unnecessary CD-matrix based
transformation CD(CD^{-1})=I which in itself results in about 50%
performance improvement.

The function has additional code to check for *divergence* condition
(as opposite to *slow convergence*) when the consecutive approximations
*diverge*. This can happen when the pixel coordinates are *far* away
from the reference pixel (distortions are "severely" extrapolated), or
for poorly-defined geometric distortions.

A "new and improved" :) `NoConvergence` exception has attributes that allows users to retrieve best solution reached by the function, and indices of *slowly converging* solutions (points) that *did not converge* in the specified maximum number of iterations as well as induced of *diverging* solutions (likely an indication of poorly described distortions or very large distortions).

Finally, the new test function `test_all_world2pix_direct` (proposed as
a replacement to the old test) provides a much more improved
(versatile) and *thorough* test of the `all_world2pix` function.

These improvements address the issues raised in #1977 and #2294 (in
addition to performance improvements).