Time scalar & multi-dimensional, using new ERFA

With the new erfa routines, scalars are OK, and hence much of the code required in `Time` to keep track of scalar values could be removed. At the same time, arbitrary dimensions become possible:
```
In [1]: from astropy.time import Time; import numpy as np, astropy.units as u

In [2]: t = Time(50000.0, format='mjd')

In [3]: t2 = t + np.arange(6.).reshape(3,2) * u.day

In [4]: t2
Out[4]: 
<Time object: scale='utc' format='mjd' value=[[ 50000.  50001.]
 [ 50002.  50003.]
 [ 50004.  50005.]]>
```

This is work in progress mostly since tests of the multidimensional behaviour are still missing (though it passes all normal tests!). @taldcroft - please have a look!

p.s.  This relies on #3135; check the second commit to see just the changes in `Time`. EDIT -- this now stands on its own, as the required changes to ERFA were part of #3141