Skip to content

Commit

Permalink
FIX: Rewrapped ?gels and ?gels_lwork routines (scipy#8333)
Browse files Browse the repository at this point in the history
* FIX: Rewrapped ?gels and ?gels_lwork

and added fat/tall matrix size args test
  • Loading branch information
ilayn authored Feb 13, 2018
1 parent 9581f52 commit 27dd818
Show file tree
Hide file tree
Showing 2 changed files with 43 additions and 78 deletions.
104 changes: 26 additions & 78 deletions scipy/linalg/flapack.pyf.src
Original file line number Diff line number Diff line change
Expand Up @@ -1464,102 +1464,50 @@ interface

end subroutine <prefix2c>gesvd_lwork

subroutine <prefix2>gels(trans,m,n,nrhs,a,lda,b,ldb,work,lwork,info)
subroutine <prefix>gels(trans,m,n,nrhs,a,lda,b,ldb,work,lwork,info)

! lqr,x,info = gels(a,b,lwork=..,overwrite_a=False,overwrite_b=False)
! Solve Minimize 2-norm(A * X - B).

fortranname <prefix2>gels
callstatement (*f2py_func)(trans,&m,&n,&nrhs,a,&lda,b,&ldb,work,&lwork,&info)
callprotoargument char*,int*,int*,int*,<ctype2>*,int*,<ctype2>*,int*,<ctype2>*,int*,int*
callprotoargument char*,int*,int*,int*,<ctype>*,int*,<ctype>*,int*,<ctype>*,int*,int*

character optional,intent(in),check(*trans=='N'||*trans=='T'):: trans = 'N'
integer depend(b),intent(hide):: nrhs = shape(b,1)
integer depend(a),intent(hide):: lda = shape(a,0)
character optional,intent(in),check(*trans=='N'||*trans==<'T',\0,'C',\2>):: trans = 'N'
integer intent(hide),depend(a):: m = shape(a,0)
integer intent(hide),depend(a):: n = shape(a,1)
integer depend(b),intent(hide):: ldb = shape(b,0)
integer optional,intent(in),depend(nrhs,m,n),&
check(lwork>=1||lwork==-1) &
:: lwork=MAX( MIN(m,n) + MAX(MIN(m,n), nrhs ), 1 )
<ftype2> dimension(m,nrhs),intent(in,out,copy,out=x) :: b
<ftype2> dimension(m,n),intent(in,out,copy,out=lqr) :: a
<ftype2> depend(lwork),dimension(MAX(lwork,1)),intent(hide,cache):: work
integer intent(hide),depend(b):: nrhs = shape(b,1)
<ftype> dimension(m,n),intent(in,out,copy,out=lqr):: a
integer intent(hide),depend(a):: lda = MAX(1,shape(a,0))
<ftype> intent(in,out,copy,out=x),depend(trans,m,n),dimension(MAX(m,n),nrhs),check(shape(b,0)==MAX(m,n)) :: b
integer depend(b),intent(hide):: ldb = MAX(1,MAX(m,n))
integer optional,intent(in),depend(nrhs,m,n),check(lwork>=1||lwork==-1)::lwork=MAX(MIN(m,n)+MAX(MIN(m,n),nrhs),1)
<ftype> depend(lwork),dimension(MAX(1,lwork)),intent(hide,cache):: work
integer intent(out)::info

end subroutine <prefix2>gels

subroutine <prefix2>gels_lwork(trans,m,n,nrhs,a,lda,b,ldb,work,lwork,info)
end subroutine <prefix>gels

! work,info = gels_lwork(m,n,nrhs,trans='N')
! Query for optimal lwork size
subroutine <prefix>gels_lwork(trans,m,n,nrhs,a,lda,b,ldb,work,lwork,info)
! ?GELS LWORK Query for optimal block size

fortranname <prefix2>gels
fortranname <prefix>gels
callstatement (*f2py_func)(trans,&m,&n,&nrhs,&a,&lda,&b,&ldb,&work,&lwork,&info)
callprotoargument char*,int*,int*,int*,<ctype2>*,int*,<ctype2>*,int*,<ctype2>*,int*,int*

character optional,intent(in),check(*trans=='N'||*trans=='T'):: trans = 'N'
integer intent(in) :: m
integer intent(in) :: n
integer intent(in) :: nrhs
integer depend(m),intent(hide):: lda = m
integer depend(m),intent(hide):: ldb = m
integer intent(hide) :: lwork = -1
<ftype2> intent(hide) :: a
<ftype2> intent(hide) :: b
<ftype2> intent(out) :: work
integer intent(out)::info

end subroutine <prefix2>gels_lwork

subroutine <prefix2c>gels(trans,m,n,nrhs,a,lda,b,ldb,work,lwork,info)

! lqr,x,info = gels(a,b,lwork=..,overwrite_a=False,overwrite_b=False)
! Solve Minimize 2-norm(A * X - B).

fortranname <prefix2c>gels
callstatement (*f2py_func)(trans,&m,&n,&nrhs,a,&lda,b,&ldb,work,&lwork,&info)
callprotoargument char*,int*,int*,int*,<ctype2c>*,int*,<ctype2c>*,int*,<ctype2c>*,int*,int*

character optional,intent(in),check(*trans=='N'||*trans=='T'):: trans = 'N'
integer depend(b),intent(hide):: nrhs = shape(b,1)
integer depend(a),intent(hide):: lda = shape(a,0)
integer intent(hide),depend(a):: m = shape(a,0)
integer intent(hide),depend(a):: n = shape(a,1)
integer depend(b),intent(hide):: ldb = shape(b,0)
integer optional,intent(in),depend(nrhs,m,n),&
check(lwork>=1||lwork==-1) &
:: lwork=MAX( MIN(m,n) + MAX(MIN(m,n), nrhs ), 1 )

<ftype2c> dimension(m,n),intent(in,out,copy,out=lqr) :: a
<ftype2c> dimension(m,nrhs),intent(in,out,copy,out=x) :: b
<ftype2c> depend(lwork),dimension(MAX(lwork,1)),intent(hide,cache):: work
integer intent(out)::info

end subroutine <prefix2c>gels

subroutine <prefix2c>gels_lwork(trans,m,n,nrhs,a,lda,b,ldb,work,lwork,info)
callprotoargument char*,int*,int*,int*,<ctype>*,int*,<ctype>*,int*,<ctype>*,int*,int*

! work,info = gels_lwork(m,n,nrhs,trans='N')
! Query for optimal lwork size
character optional,intent(in),check(*trans=='N'||*trans==<'T',\0,'C',\2>):: trans = 'N'
integer intent(in),check(m>=0) :: m
integer intent(in),check(n>=0) :: n
integer intent(in),check(nrhs>=0) :: nrhs

fortranname <prefix2c>gels
callstatement (*f2py_func)(trans,&m,&n,&nrhs,&a,&lda,&b,&ldb,&work,&lwork,&info)
callprotoargument char*,int*,int*,int*,<ctype2c>*,int*,<ctype2c>*,int*,<ctype2c>*,int*,int*
<ftype> intent(hide):: a
integer intent(hide):: lda = MAX(1,m)
<ftype> intent(hide):: b
integer intent(hide):: ldb = MAX(1,MAX(m,n))
integer intent(hide):: lwork=-1

character optional,intent(in),check(*trans=='N'||*trans=='T'):: trans = 'N'
integer intent(in) :: m
integer intent(in) :: n
integer intent(in) :: nrhs
integer depend(m),intent(hide):: lda = m
integer depend(m),intent(hide):: ldb = m
integer intent(hide) :: lwork = -1
<ftype2c> intent(hide) :: a
<ftype2c> intent(hide) :: b
<ftype2c> intent(out) :: work
<ftype> intent(out):: work
integer intent(out)::info

end subroutine <prefix2>gels_lwork
end subroutine <prefix>gels_lwork

subroutine <prefix2>gelss(m,n,minmn,maxmn,nrhs,a,b,s,cond,r,work,lwork,info)

Expand Down
17 changes: 17 additions & 0 deletions scipy/linalg/tests/test_lapack.py
Original file line number Diff line number Diff line change
Expand Up @@ -146,6 +146,23 @@ def test_clapack(self):
class TestLeastSquaresSolvers(object):

def test_gels(self):
seed(1234)
# Test fat/tall matrix argument handling - gh-issue #8329
for ind, dtype in enumerate(DTYPES):
m = 10
n = 20
nrhs = 1
a1 = rand(m, n).astype(dtype)
b1 = rand(n).astype(dtype)
gls, glslw = get_lapack_funcs(('gels', 'gels_lwork'), dtype=dtype)

# Request of sizes
lwork = _compute_lwork(glslw, m, n, nrhs)
_, _, info = gls(a1, b1, lwork=lwork)
assert_(info >= 0)
_, _, info = gls(a1, b1, trans='TTCC'[ind], lwork=lwork)
assert_(info >= 0)

for dtype in REAL_DTYPES:
a1 = np.array([[1.0, 2.0],
[4.0, 5.0],
Expand Down

0 comments on commit 27dd818

Please sign in to comment.