reimplement pchip interpolator

Convert `pchip` interpolator into a wrapper around `BPoly`. The interpolation algorithm is the same.

Related PRs: gh-3260 (root finding), gh-3112 (Akima splines). I guess a relation between this algorithm and Akimas should be clarified at some stage.

A couple of points I'd like to flag:
* As it stands, there are two slightly different API flavors for derivatives. `Krogh`-based Interpolators in `polyint` define `interpolator.derivative(x)` to evaluate a *value* of a derivative at point `x`. Meanwhile, `PPoly.derivative()` returns a callable object representing the derivative. Which I personally find more attractive. So far in this PR I've implemented first version in the name of backwards compatibility, but I'd like to change it to PPoly way --- any strong opinions on this?
* So far I did not implement `derivatives`, which  are supposed to return a list of derivatives at a point: easy to do, but do we actually need this?
* How much do we want to expose an underlying `BPoly` object. An alternative implementation would be to just have a function to return `BPoly` object and get rid of the `PchipInterpolator` class altogether. In this case we'll need to add an `axis` parameter to PPoly/BPoly.