New user self-registration is disabled due to spam. Please email eigen-core-team @ lists.tuxfamily.org if you need an account.
Before reporting a bug, please make sure that your Eigen version is up-to-date!
Bug 1535 - Allow extraction of eigenvalues from a Schur matrixT()
Summary: Allow extraction of eigenvalues from a Schur matrixT()
Status: DECISIONNEEDED
Alias: None
Product: Eigen
Classification: Unclassified
Component: Eigenvalues (show other bugs)
Version: unspecified
Hardware: All All
: Normal API Change
Assignee: Nobody
URL:
Whiteboard:
Keywords:
Depends on:
Blocks:
 
Reported: 2018-04-03 09:00 UTC by Jan van Dijk
Modified: 2018-04-04 09:29 UTC (History)
5 users (show)



Attachments
Implement RealSchur::eigenvalues() (7.31 KB, patch)
2018-04-04 09:10 UTC, Jan van Dijk
no flags Details | Diff
Use RealSchur::eigenvalues() (2.37 KB, patch)
2018-04-04 09:13 UTC, Jan van Dijk
no flags Details | Diff

Description Jan van Dijk 2018-04-03 09:00:26 UTC
In my application I need the real Schur form of a matrix, as well as its eigenvalues (as provided by Lapack's DGEES).

Eigen's RealSchur will give me the Schur matrices U and T, but does not offer an interface to subsequently extract the eigenvalues from the 1x1 or 2x2 diagonal blocks of T. Such interface would be highly appreciated.

It seems that this would merely require a small factorization of the member compute() of the EigenSolver template: the part below the comment:

  // Compute eigenvalues from matT

in EigenSolver.h is exactly what I need --- but rather not duplicate.

IMHO, the cleanest solution is to move this code to a new public member RealSchur<>::eigenvalues(), that is also used by EigenSolver<>::compute().
Comment 1 Gael Guennebaud 2018-04-03 12:59:37 UTC
I see 3 options:

1 - call it RealSchur::eigenvalues(), and compute them at the end of RealSchur::compute (), similar to EigenSolver. Perhaps, with an option to compute() to request them?

2 - call it RealSchur::eigenvalues(), and compute them in a cached vector at the first call to eigenvalues(). The cache would be invalidated at every call to compute().

3 - no cache, and thus we would call it extractEigenvalues() to make it clear that each call to this function involves some computations.


To be more consistent with the other classes, I would go for either 1 or 2.
Comment 2 Jan van Dijk 2018-04-04 09:10:11 UTC
Created attachment 843 [details]
Implement RealSchur::eigenvalues()

This patch implements adds member eigenvalues() to the RealSchur template:

 * new private member m_eigenvalues hosts the eigenvalues
   of the matrix; these are calculated as part of the factorization.
   This vector can be retrieved through accessor eigenvalues();
 * new private member extractEigenvalues() resizes m_eigenvalues
   and extracts the eigenvalues from a previously calculated m_T.
   This member returns NumericalIssue if a !isfinite eigenvalue
   is encountered, Success otherwise;
 * computeFromHessenberg() calls extractEigenvalues() after doing the
   factorization. If the factorization generated an error, that
   error is kept for later inspection (Convergence failure, ...),
   otherwise m_info is overwritten by the result of extractEigenvalues();
 * compute() calls computeFromHessenberg() on a scaled matrix (as before).
   After that function returns, the (inverse) scaling factor is applied to
   m_eigenvalues (as it is to m_T itself);
 * in the special case of a near-zero matrix m_eigenvalues if manually
   set to zero (as is m_T).

NOTE (1):
Even if eigenvalues() is not called, there is one user-visible change:
m_info is set to NumericalIssue if non-finite eigenvalues are encountered;
before these went undetected.

NOTE (2):
This patch does not contain additional tests. This new functionality will be
indirectly tested by the EigenSolver tests once EigenSolver starts to use
RealSchur::eigenvalues(). (That is part of a patch 2/2).
Comment 3 Jan van Dijk 2018-04-04 09:13:54 UTC
Created attachment 844 [details]
Use RealSchur::eigenvalues()

This patch changes EigenSolver to use RealSchur::eigenvalues() to obtain the eigenvalues.
Comment 4 Jan van Dijk 2018-04-04 09:29:45 UTC
The two patches implement Gaels' proposition #1. I checked that all eigensolver_generic* and real_schur* tests pass after applying these patches.

Preparing the patches cost me a surprising amount of time, mostly because the intended error checking semantics of EigenSolver are not always clear (to me). As an example, let us assume that the computation of the Schur decomposition fails:
 * the EigenSolver object is marked to be initialized, although not all matrices are set (or even resized in accordance with the size of the incoming matrix).
 * m_eigenvectorsOk is set to computeEigenVectors, while m_matU is not calculated (or even resized) in case of a Schur failure.

I am happy to provide a documentation or code fix, but need some guidance here: what is the intended state of the RealSchur and EigenSolver objects after compute is called but does not return successfully? Should matrixT() and eigenvalues() (and maybe matrixU()) be resized at least? Or are no guarantees provided whatsoever in that case, and is it possible that the sizes are even inconsistent (with each other or with the matrix that was passed to compute)?

Note You need to log in before you can comment on or make changes to this bug.