Change default qr subroutine to one that allows m<n. Update documentation for fat QR. Return permutation vector instead of matrix in pivoted QR.

This should address the issues in #4555, so now
```
julia> A=randn(2,4)
2x4 Array{Float64,2}:
 1.14841  0.148657   0.50509    0.870583
 0.5955   1.78331   -0.352758  -1.95542 

julia> qr(A)
(
2x2 Array{Float64,2}:
 -0.887745  -0.460335
 -0.460335   0.887745,

2x4 Array{Float64,2}:
 -1.29362  -0.952887  -0.286004   0.127291
  0.0       1.51469   -0.54567   -2.13667 )
```

@dmbates I browsed through the GLM code and couldn't find any call directly to `geqrt3` anymore. If I have missed it, this pull request might affect GLM because I have changed `qrfact` to call `geqrt`. They are very similar as `geqrt` is just a higher level interface to `geqrt3`.