ENH: Add bounded linear least squares routine. Closes #3130.

This works, but I'd like to get some input from f2py experts (cc // @pv) on a few points. You can run the code with the below. The Fortran code is hard-coded for float32 precision. Presumably, that shouldn't be a big deal to change. However, I'm running into some other issues that are I think related to using Fortran 90. I've never wrapped Fortran 90 code before.

I can't change the intent of, say, a to be (in, out, copy). If I do, I receive an error like this from f2py. That's the main sticking point for me, I think. However, the wrapping code makes a copy, so maybe we can live with it?

```
gfortran:f77: /tmp/tmp_TCHXo/src.linux-x86_64-2.7/_bvls-f2pywrappers.f
/tmp/tmp_TCHXo/src.linux-x86_64-2.7/_bvls-f2pywrappers.f:27.31:

          real(kind=kind(one)), dimension(:,:),intent(in,out,copy)      
                               1
Error: Invalid character in name at (1)
/tmp/tmp_TCHXo/src.linux-x86_64-2.7/_bvls-f2pywrappers.f:39.16:

      call bvls(a, b, bnd, x, rnorm, nsetp, w, index_bn, ierr)          
                1
Error: Rank mismatch in argument 'a' at (1) (scalar and rank-2)
/tmp/tmp_TCHXo/src.linux-x86_64-2.7/_bvls-f2pywrappers.f:27.31:

          real(kind=kind(one)), dimension(:,:),intent(in,out,copy)      
                               1
Error: Invalid character in name at (1)
/tmp/tmp_TCHXo/src.linux-x86_64-2.7/_bvls-f2pywrappers.f:39.16:

      call bvls(a, b, bnd, x, rnorm, nsetp, w, index_bn, ierr)          
                1
Error: Rank mismatch in argument 'a' at (1) (scalar and rank-2)
error: Command "/usr/bin/gfortran -Wall -ffixed-form -fno-second-underscore -fPIC -O3 -funroll-loops -I/tmp/tmp_TCHXo/src.linux-x86_64-2.7 -I/usr/local/lib/python2.7/dist-packages/numpy/core/include -I/usr/include/python2.7 -c -c /tmp/tmp_TCHXo/src.linux-x86_64-2.7/_bvls-f2pywrappers.f -o /tmp/tmp_TCHXo/tmp/tmp_TCHXo/src.linux-x86_64-2.7/_bvls-f2pywrappers.o" failed with exit status 1
```


```
import numpy as np
from scipy import optimize

x = np.array([-1.7237128E+00,1.8712276E+00,-9.6608055E-01,
                  -2.8394297E-01,1.3416969E+00,1.3757038E+00,
                  -1.3703436E+00,4.2581975E-02,-1.4970151E-01,
                  8.2065094E-01])[:,None]
x = np.column_stack((np.ones_like(x), x))
np.random.seed(12)
y = np.dot(x, [1.8, 5.4]) + np.random.random(len(x))
bounds = None

print optimize.bounded_lstsq(x, y, bounds=bounds)

bounds = [(None, None), (None, 4)]
print optimize.bounded_lstsq(x, y, bounds=bounds)
```