Commit Graph

73 Commits

Author SHA1 Message Date
Gael Guennebaud
c44bbabdcc fix LU and QR solve when rank==0, fix LLT when the matrix is purely 0 2010-07-12 16:42:38 +02:00
Stuart Glaser
58fc972ed3 LU on limited-size matrices no longer allocates for temporaries. 2010-06-21 21:29:45 -07:00
Benoit Jacob
a9eabed421 Patch by 'Wolf' from the issue tracker:
Fix bug #90, missing type cast in LU, allow to use LU with MPFR.
2010-02-02 07:06:15 -05:00
Benoit Jacob
3e963ee69d EIGEN_ENUM_MIN ---> EIGEN_SIZE_MIN 2010-01-26 20:37:57 -05:00
Benoit Jacob
6cc9dc17f2 In LU / Inverse, decouple the output type from the input type.
This has long been done in the default branch
2010-01-26 18:45:23 -05:00
Benoit Jacob
6fd9248c09 add Intel copyright info 2009-12-15 08:43:31 -05:00
Benoit Jacob
4262117f84 backport 4x4 inverse changes:
- use cofactors
 - use Intel's SSE code in the float case
2009-12-15 08:16:48 -05:00
Benoit Jacob
8d31f58ea1 fix bug #70
Was trying to apply stupid invertibility check to top-left 2x2 corner.
2009-11-26 15:33:07 -05:00
Benoit Jacob
fc7f39980c backport improvement in 4x4 inverse precisio, and rigorous precision test 2009-11-23 10:27:10 -05:00
Benoit Jacob
de693cf34a remove extra ; 2009-10-28 10:04:13 -04:00
Benoit Jacob
21c4e0802d fix potential warning 2009-10-28 09:45:09 -04:00
Benoit Jacob
241b9d34a7 Hey, I was insomniac too ;)
This restores much of the performance benefit of Euler's method, without compromising accuracy (tested on 1e+7 matrices). Namely, my benchmark now runs in 1.5 s instead of 2.2 s. The same in the default branch runs in 1.08 s instead of 1.9 s, so the default branch benefits even more!
2009-10-28 03:50:29 -04:00
Benoit Jacob
9e15a6da2e Fix 4x4 matrix inversion. Applying Euler's trick is more tricky than what "high performance matrix inversion" websites would have you believe. Our 4x4 matrix inversion wasn't numerically stable, because in applying the Euler trick we didn't take the 2x2 block of biggest determinant. As a result, with float, we got relative errors above 1% every 1000 random matrices, and we got completely wrong results every 10000 matrices.
Note that this decreases the performance, but we're still significantly faster than the brutal cofactors approach.
2009-10-27 07:35:25 -04:00
Gael Guennebaud
c9f7a19053 make LU::solve() not to crash when rank=0
(transplanted from fe813911f2
)
2009-08-09 00:06:53 +02:00
Benoit Jacob
f5ae3a4b5a result of our experiments with LU tuning: implement very simple formula, that
turns out to be similar to Higham's formula already in use in LDLt
2009-05-07 20:35:26 +00:00
Benoit Jacob
3c3653b9de merge 953719: fix 4x4 inverse 2009-04-14 13:43:21 +00:00
Benoit Jacob
062d6bd47a documentation update/improvement 2009-01-26 16:22:40 +00:00
Benoit Jacob
b05c5dd1be compute the rank on the fly rather than at the end, and stop early
in the case of non-full rank (so big optimization in case the rank is low)
2009-01-26 16:14:20 +00:00
Benoit Jacob
a79deafc6d * mark Geometry as experimental
* install QtAlignedMalloc
* finish the renaming Regression->LeastSquares
* install LeastSquares directory (!!!)
* misc dox fixes
2009-01-26 13:53:43 +00:00
Benoit Jacob
00d7f8e567 * solveTriangularInPlace(): take a const ref and const_cast it, to allow passing temporary xprs.
* improvements, simplifications in LU::solve()
* remove remnant of old norm2()
2009-01-25 23:46:51 +00:00
Benoit Jacob
414ee1db4b Optimization in LU::solve: when rows<=cols, no need to compute the L matrix
Remove matrixL() and matrixU() methods: they were tricky, returning a Part,
and matrixL() was useless for non-square LU. Also they were untested. This is
the occasion to simplify the docs (class_LU.cpp) removing the most confusing part.
I think that it's better to let the user do his own cooking with Part's.
2009-01-25 16:33:06 +00:00
Benoit Jacob
81cb887baf fix bug in the computation of rank
very difficult to catch in unit-tests because this is very noisy
2009-01-20 16:21:56 +00:00
Benoit Jacob
15ca6659ac * the 4th template param of Matrix is now Options. One bit for storage
order, one bit for enabling/disabling auto-alignment. If you want to
disable, do:
Matrix<float,4,1,Matrix_DontAlign>
The Matrix_ prefix is the only way I can see to avoid
ambiguity/pollution. The old RowMajor, ColMajor constants are
deprecated, remain for now.
* this prompted several improvements in matrix_storage. ei_aligned_array
renamed to ei_matrix_array and moved there. The %16==0 tests are now
much more centralized in 1 place there.
* unalignedassert test: updated
* update FindEigen2.cmake from KDElibs
* determinant test: use VERIFY_IS_APPROX to fix false positives; add
testing of 1 big matrix
2009-01-04 15:26:32 +00:00
Benoit Jacob
d9e5fd393a * In LU solvers: no need anymore to use row-major matrices
* Matrix: always inherit WithAlignedOperatorNew, regardless of
vectorization or not
* rename ei_alloc_stack to ei_aligned_stack_alloc
* mixingtypes test: disable vectorization as SSE intrinsics don't allow
mixing types and we just get compile errors there.
2009-01-03 22:33:08 +00:00
Benoit Jacob
9e00d94543 * the Upper->UpperTriangular change
* finally get ei_add_test right
2008-12-20 13:36:12 +00:00
Benoit Jacob
f34a4fa335 * forgot to svn add 2 files
* idea of Keir Mierle: make the static assert error msgs UPPERCASE
2008-12-18 21:04:06 +00:00
Benoit Jacob
fabaa6915b * fix in IO.h, a useless copy was made because of assignment from
Derived to MatrixBase.
* the optimization of eval() for Matrix now consists in a partial
  specialization of ei_eval, which returns a reference type for Matrix.
  No overriding of eval() in Matrix anymore. Consequence: careful,
  ei_eval is no longer guaranteed to give a plain matrix type!
  For that, use ei_plain_matrix_type, or the PlainMatrixType typedef.
* so lots of changes to adapt to that everywhere. Hope this doesn't
  break (too much) MSVC compilation.
* add code examples for the new image() stuff.
* lower a bit the precision for floats in the unit tests as
  we were already doing some workarounds in inverse.cpp and we got some
  failed tests.
2008-12-18 20:36:25 +00:00
Benoit Jacob
c22d10f94c LU class:
* add image() and computeImage() methods, with unit test
* fix a mistake in the definition of KernelResultType
* fix and improve comments
2008-12-17 16:47:55 +00:00
Benoit Jacob
89f468671d * replace postfix ++ by prefix ++ wherever that makes sense in Eigen/
* fix some "unused variable" warnings in the tests; there remains a libstdc++ "deprecated"
warning which I haven't looked much into
2008-12-17 14:30:01 +00:00
Benoit Jacob
c1e2156d8a * Much better, consistent error msgs when mixing different scalar types:
- in matrix-matrix product, static assert on the two scalar types to be the same.
- Similarly in CwiseBinaryOp. POTENTIALLY CONTROVERSIAL: we don't allow anymore binary
  ops to take two different scalar types. The functors that we defined take two args
  of the same type anyway; also we still allow the return type to be different.
  Again the reason is that different scalar types are incompatible with vectorization.
  Better have the user realize explicitly what mixing different numeric types costs him
  in terms of performance.
  See comment in CwiseBinaryOp constructor.
- This allowed to fix a little mistake in test/regression.cpp, mixing float and double
- Remove redundant semicolon (;) after static asserts
2008-12-03 21:01:55 +00:00
Benoit Jacob
00f89a8f37 Update e-mail address 2008-11-24 13:40:43 +00:00
Gael Guennebaud
5066fe8bbe * sparse LU: add extraction of L,U,P, and Q, as well as determinant
for both backends.
* extended a bit the sparse unit tests
2008-10-20 17:03:09 +00:00
Gael Guennebaud
765219aa51 Big API change in Cholesky module:
* rename Cholesky to LLT
 * rename CholeskyWithoutSquareRoot to LDLT
 * rename MatrixBase::cholesky() to llt()
 * rename MatrixBase::choleskyNoSqrt() to ldlt()
 * make {LLT,LDLT}::solve() API consistent with other modules

Note that we are going to keep a source compatibility untill the next beta release.
E.g., the "old" Cholesky* classes, etc are still available for some time.
To be clear, Eigen beta2 should be (hopefully) source compatible with beta1,
and so beta2 will contain all the deprecated API of beta1. Those features marked
as deprecated will be removed in beta3 (or in the final 2.0 if there is no beta 3 !).

Also includes various updated in sparse Cholesky.
2008-10-13 15:53:27 +00:00
Gael Guennebaud
75649551c2 compilation fixes with MSVC 2008-09-03 11:26:19 +00:00
Gael Guennebaud
7e8aa63bb7 * Add Hyperplane::transform(Matrix/Transform)
* Fix compilations with gcc 3.4, ICC and doxygen
* Fix krazy directives (hopefully)
2008-08-31 13:32:29 +00:00
Gael Guennebaud
f2f48b6560 * remove LargeBit and related stuff
* replaced the Flags template parameter of Matrix by StorageOrder
  and move it back to the 4th position such that we don't have to
  worry about the two Max* template parameters
* extended EIGEN_USING_MATRIX_TYPEDEFS with the ei_* math functions
2008-08-23 17:11:44 +00:00
Gael Guennebaud
a6d387a359 Various compilation fixes for MSVC 9. All tests compile but some
still fail at runtime in ei_aligned_free() (even without vectorization).
2008-08-19 11:06:40 +00:00
Gael Guennebaud
95dd09bea6 * revert the previous interface change in solveTriangular (pointer vs reference)
* remove the cast operators in the Geometry module: they are replaced by constructors
  and new operator= in Matrix
* extended the operations supported by Rotation2D
* rewrite in solveTriangular:
  - merge the Upper and Lower specializations
  - big optimization of the path for row-major triangular matrices
2008-08-18 22:17:42 +00:00
Gael Guennebaud
8a3e6b1ee2 change solveTriangularInPlace() to take a pointer as input (as discussed on IRC).
extended the documentation of the triangular solver.
2008-08-12 07:49:59 +00:00
Benoit Jacob
f04c1cb774 Complete LU documentation 2008-08-11 21:26:37 +00:00
Benoit Jacob
17ec407ccd * add Regression module, from eigen1, improved, with doc and unit-test.
* fix .normalized() so that Random().normalized() works; since the return
type became complicated to write down i just let it return an actual
vector, perhaps not optimal.
* add Sparse/CMakeLists.txt. I suppose that it was intentional that it
didn't have CMakeLists, but in <=2.0 releases I'll just manually remove
Sparse.
2008-08-11 02:25:40 +00:00
Gael Guennebaud
b13148c358 renamed inverseProduct => solveTriangular 2008-08-09 20:06:25 +00:00
Benoit Jacob
d6e88f8155 * add LU unit-test. Seems like we have very good numerical stability!
* some cleanup, and grant me a copyright line on the determinant test.
2008-08-09 19:26:14 +00:00
Benoit Jacob
becbeda50a * reimplement the general case of inverse() on top of LU. Advantages:
- removes much code
  - 2.5x faster (even though LU uses complete pivoting contrary to what
inverse used to do!)
  - there _were_ numeric stability problems with the partial pivoting
approach of inverse(), with 200x200 matrices inversion failed almost
half of the times (overflow). Now these problems are solved thanks to
complete pivoting.

* remove some useless stuff in LU
2008-08-09 15:50:07 +00:00
Benoit Jacob
681e944670 *implement LU solver (solves any rectangular system)
*in test/CMakeLists : modify EI_ADD_TEST so that 2nd argument is
additional compiler flags. used to add -O2 to test_product_large so it
doesn't take forever.
2008-08-09 15:23:54 +00:00
Benoit Jacob
a41f2b4216 * fix bug in SwapWrapper : store the wrapped expression by reference
* optimize setIdentity: when the matrix is large enough it is better to
  setZero() and overwrite the diagonal
* start of LU solver, disabled for now
2008-08-09 04:37:09 +00:00
Benoit Jacob
9bbe396939 forgot to remove that old code 2008-08-07 21:50:29 +00:00
Benoit Jacob
5f350448e1 - add kernel computation using the triangular solver
- take advantage of the fact that our LU dec sorts the eigenvalues of U
in decreasing order
- add meta selector in determinant
2008-08-07 21:48:21 +00:00
Benoit Jacob
58ba9ca72f LU: remove partial-pivoting path (moderately useful since it's does
not allow to easily get the rank), fix a bug (which could have been
triggered by matrices having coefficients of very different
magnitudes).
Part: add an assert to prevent hard to find bugs
Swap: update comments
2008-08-07 04:31:05 +00:00
Benoit Jacob
09ef7db9d9 Add partial pivoting runtime option to LU.
Note: in fact, inverse() always uses partial pivoting because the algo
currently used doesn't make sense with complete pivoting. No num
stability issue so far even with size 200x200. If there is any problem
we can of course reimplement inverse on top of LU.
2008-08-05 15:43:11 +00:00