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 1493 - dense Q extraction and solve is sometimes erroneous for complex matrices
Summary: dense Q extraction and solve is sometimes erroneous for complex matrices
Status: RESOLVED FIXED
Alias: None
Product: Eigen
Classification: Unclassified
Component: QR (show other bugs)
Version: 3.4 (development)
Hardware: All All
: Normal Wrong Result
Assignee: Nobody
URL:
Whiteboard:
Keywords:
Depends on:
Blocks: 3.4
  Show dependency treegraph
 
Reported: 2017-12-20 00:29 UTC by Jeff Trull
Modified: 2018-04-17 09:41 UTC (History)
2 users (show)



Attachments
Demonstrating that Q is applied incorrectly in some cases by using an identity RHS (1.04 KB, text/x-c++src)
2017-12-20 00:29 UTC, Jeff Trull
no flags Details

Description Jeff Trull 2017-12-20 00:29:37 UTC
Created attachment 809 [details]
Demonstrating that Q is applied incorrectly in some cases by using an identity RHS

HouseholderQR and ColPivHouseholderQR apply Q incorrectly on the right side, and Q' incorrectly on the left side, for the case of complex numbers. The latter case impacts solving as well, because Eigen needs to calculate Q'*RHS.

A simple test is to calculate the decomposition of a small matrix, and then apply it to the left and right of an identity matrix of the appropriate size.  Both answers should be the same, but they are not. I've attached an example.
Comment 1 Jeff Trull 2017-12-20 16:52:11 UTC
My proposed fix is https://bitbucket.org/eigen/eigen/pull-requests/355
Comment 2 Christoph Hertzberg 2018-01-19 16:40:00 UTC
The current behavior indeed can't be right.
Most importantly the patch should include unit tests to prevent regressions (extent the tests in test/householder.cpp)
Comment 3 Christoph Hertzberg 2018-04-15 08:55:21 UTC
I finally pushed a fix (to the default branch only):
https://bitbucket.org/eigen/eigen/commits/25b342390facc827

It does slightly change the behavior for several .apply*OnThe{Left,Right}() methods (only for complex scalars), so I'm not 100% sure if this should be backported to 3.3
Comment 4 Gael Guennebaud 2018-04-16 13:01:00 UTC
Thank you for the fix.

The documentation of setTrans()/trans() are now incorrect, we should really rename them setReverseFlag/reverseFlag(). Ideally we would keep setTrans()/trans() for compatibility with only a deprecation warning for real scalar types, and a compilation error for complex types, because in the later case the user code is likely broken.
Comment 5 Christoph Hertzberg 2018-04-16 13:32:54 UTC
I agree, we should rename those (and the m_trans member) -- alternative encode the transposing into the type.
I'm not sure if we really need to keep the setTrans()/trans() methods. Those were protected members, so at least they were not easy to access.
The public interface have been the .transpose()/.adjoint() methods for a while.
Comment 6 Christoph Hertzberg 2018-04-17 09:41:30 UTC
Renaming trans to reverseFlag pushed here:

https://bitbucket.org/eigen/eigen/commits/6acd5a89a93cdc

Btw, are these constructs necessary?

typename internal::conditional<NumTraits<Scalar>::IsComplex,
  typename internal::remove_all<typename VectorsType::ConjugateReturnType>::type,
  VectorsType>::type,

It seems that the same logic already happens inside `VectorsType::ConjugateReturnType` -- everything works fine with, e.g.:

typedef HouseholderSequence<
  typename internal::remove_all<typename VectorsType::ConjugateReturnType>::type,
  typename internal::remove_all<typename CoeffsType::ConjugateReturnType>::type,
  Side
> ConjugateReturnType;

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