Add Eigenvalue Range Functionality for Symmetric Eigenvalue Problems #6510
Conversation
- Utilizing prefixes - Changing intent on vl and vu from hide to in
After exposing the parameters vl and vu in the signature for the fortran subroutines, the routines are wrapped to preserve the original ordering of arguments and append vl and vu at the end of the list.
Tests make sure the previously calling format is maintained and vl and vu are pushed to the end of the arg list.
|
I really don't get why these changes aren't passing. I just did the following:
So I really don't get what's going on. Hopefully it can be worked out. This would be a useful feature to have. REMOVED |
|
A google search for "*** glibc detected *** python: free(): invalid next size (normal)" turns of a lot of results from ~4-6 years ago. Sounds like it has to do with some old libraries or Linux builds. Not sure why it's failing here and not elsewhere. |
|
It's typically a sign of memory corruption occurring. As said earlier, these bugs don't always manifest if you're lucky enough nothing important is in the place where stuff gets overwritten. The backtrace it gives is useless, because the crashing free() is a symptom, and the cause is elsewhere. |
scipy/linalg/flapack.pyf.src
Outdated
| integer intent(hide),dimension(2*m) :: isuppz | ||
| integer intent(in),depend(n) :: lwork=18*n | ||
| complex intent(hide),dimension(lwork) :: work | ||
| <ftype2c> intent(hide),dimension(lwork) :: work | ||
| integer intent(hide),depend(n) :: lrwork=24*n | ||
| real intent(hide),dimension(lrwork) :: rwork |
|
iow, it's a sign of something wrong, and in this case it seems the cause could be found by carefully reading through the change made, as I wrote in the other pull request. |
|
Shoot. I'm sorry about this. I was poring over this code so much already I was blind to that. Thanks for the help. |
A new parameter `eigrng` is added to `eigh` which allows the specification of computing eigenvalues only within the half-open interval set by `eigrng`. Checks are done to ensure `eigrng` and `eigvals` are both not set and `eigrng` is a valid range.
Adding tests for the new parameter `eigrng` to accompany previous tests of `eigh`. Also checking properly raised errors.
apapanico
changed the title
Removing extra blank lines
| (vl, vu) = eigrng | ||
| w_full, v_full, info = evr(a1, uplo=uplo, jobz=_job, range="V", vl=vl, | ||
| vu=vu, overwrite_a=overwrite_a) | ||
| ix = nonzero(w_full)[0].tolist() |
There was a problem hiding this comment.
What if there are zero eigenvalues?
Is the number of eigenvalues found returned by the LAPACK routine in the M parameter?
If not, how is it possible to know how many were found?
What is the purpose of ISUPPZ parameter in the dsyevr?
The nonzero and tolist are probably not needed.
| @@ -7,6 +7,7 @@ | |||
| # additions by Bart Vandereycken, June 2006 | |||
| # additions by Andrew D Straw, May 2007 | |||
| # additions by Tiziano Zito, November 2008 | |||
| # additions by Alex Papanicolaou, August 2016 | |||
There was a problem hiding this comment.
You don't necessarily need to add a comment here --- who modified what is kept track of by the version contorl.
In particular, the file has been modified by many people since 2008.
| @@ -240,7 +241,7 @@ def eig(a, b=None, left=False, right=True, overwrite_a=False, | |||
|
|
|||
| def eigh(a, b=None, lower=True, eigvals_only=False, overwrite_a=False, | |||
| overwrite_b=False, turbo=True, eigvals=None, type=1, | |||
| check_finite=True): | |||
| check_finite=True, eigrng=None): | |||
There was a problem hiding this comment.
I would prefer to write eig_range or just range instead of abbreviations.
| @@ -404,6 +404,73 @@ | |||
| _flapack.sgegv = sgegv | |||
| _flapack.zgegv = zgegv | |||
|
|
|||
| # Patching EVR Methods to handle value ranges | |||
There was a problem hiding this comment.
Should say that this is used to retain backward compatibility.
| wrapper.__module__ = 'scipy.linalg._flapack' | ||
| wrapper.__doc__ = _evr_doc.format( | ||
| name=name, w_type=type.lower(), type=type) | ||
| wrapper.__name__ = name |
There was a problem hiding this comment.
Should also copy the _cpointer attribute for full backward compat.
| for typ in v['dtype']: | ||
| for overwrite in v['overwrite']: | ||
| for eigenrange in v['eigrng']: | ||
| for lower in v['lower']: |
There was a problem hiding this comment.
Note that itertools.product can be used to avoid deeply nested loops:
dims = (DIM,)
dtypes = ('f', 'd', 'F', 'D')
...
parts = [dims, dtypes, ...]
for dim, typ, overwrite, eigenrange, lower in itertools.product(*parts):
yield (...)
| ssyevr = evr_decorator(ssyevr, 'ssyevr', 'f') | ||
| dsyevr = evr_decorator(dsyevr, 'dsyevr', 'd') | ||
| cheevr = evr_decorator(cheevr, 'cheevr', 'F') | ||
| zheevr = evr_decorator(zheevr, 'zheevr', 'D') |
There was a problem hiding this comment.
You need also replace the function in _flapack:
_flapack.ssyevr = ssyevr
| @@ -27,6 +27,10 @@ | |||
| from scipy.linalg.lapack import get_lapack_funcs | |||
| from scipy.linalg.blas import get_blas_funcs | |||
|
|
|||
| # Special import for monkey patched functions | |||
| from scipy.linalg.lapack import dsyevr, ssyevr, cheevr, zheevr | |||
There was a problem hiding this comment.
The change in lapack.py needs to be made so that this is not necessary.
In particular, get_lapack_funcs should return the same functions.
|
Hi @apapanico while investigating #8205 it turned out that lwork, liwork parameters were the culprit. I've done something similar and only then I've seen this PR and the related discussion. From the look of it, it covered quite some distance. But out parameters such as Are you still planning to finish this by any chance? Otherwise, I can take over if you wish. |
Updated, renamed, and cleaned up. @apapanico 2016-08-25
Add Eigenvalue Range Functionality for Symmetric Eigenvalue Problems
This PR does the following:
-evrLAPACK methods to allow for computing eigenvalues only in a half-open interval specified by the parametersvlandvu.-evrwrappers to reorder the arguments such thatvlandvucome at the end of the optional argument list. Inscipy/linalg/lapack.py, the Fortran routines are wrapped after import and provided a new doc string to match the appropriate signature. Moreover, a bug is addressed whereilandiuare NOT properly ignored in-evrmethods when not settingrange='I'despite the claims of the LAPACK documentation. For example, ifilandiuare set to 1 andrange='A', then all eigenvalues are returned but only 1 eigenvector. This is fixed by settingil=1andiu = n(the size of the matrix) whenrange='A'orrange='V'so that there are no missing eigenvectors. This approach is an improvement on Updating LAPACK -syevr methods #6503 which could not ensure a truly preserved ordering of input arguments.eigrng, is added toscipy.linalg.eighthat allows specifying this mode of use.eigrngmust be a tuple of length 2 and the first value must be less than the second value. It cannot be specified in conjunction witheigvals. Finally, the return values are all eigenvalues, and eigenvectors if desired, within the range.eighfunctionality is as it should be, and that this new mode is incompatible with the mode specified by the parametereigvalswhich computes eigenvalues corresponding to an index range.