Compare commits

..

83 Commits

Author SHA1 Message Date
Benoit Jacob
ed6eb5a625 bump 2010-02-11 21:39:41 -05:00
Benoit Jacob
9488a12125 work around brain dead ICC 2010-02-11 19:32:56 -05:00
Piotr Trojanek
7b44957c4b std:: namespace fixup for more restricive compilers such as QNX's QCC 2010-02-10 22:27:35 +01:00
Hauke Heibel
743ad75595 BenchTimer backport (clock_gettime & QueryPerformanceCounter). 2010-02-03 21:55:01 +01: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
cd34a1d351 backport bug fix by Jitse. 2010-01-28 14:00:09 -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
Gael Guennebaud
7852a48a2f fix matrix product with EIGEN_DEFAULT_TO_ROW_MAJOR 2010-01-25 21:56:01 +01:00
Benoit Jacob
d209120180 * Introduce EIGEN_DEFAULT_TO_ROW_MAJOR tests option
---> Now only product_large fails with EIGEN_DEFAULT_TO_ROW_MAJOR.
* Fix EIGEN_NO_ASSERTION_CHECKING tests option
* Fix a crash in Tridiagonalization on row-major matrices + SSE
* Fix inverse test (numeric stability noise)
* Extend map test (see previous fixes in MapBase)
* Fix vectorization_logic test for row-major
* Disable sparse tests with EIGEN_DEFAULT_TO_ROW_MAJOR
2010-01-25 14:00:02 -05:00
Thomas Capricelli
55c0707b1d fix the script again (definitely?) + cleaning 2010-01-22 19:28:33 +01:00
Benoit Jacob
72044ca925 fix a super nasty bug: on row-major expressions that are NOT vectors but that
do have LinearAccess, the MapBase::coeff(int) and MapBase::coeffRef(int)
methods were broken.
2010-01-21 23:33:20 -05:00
Benoit Jacob
c2b8ca7493 if EIGEN_DONT_ALIGN then don't try to vectorize (was giving a #error later on). 2010-01-21 22:32:16 -05:00
Gael Guennebaud
018cb8975a fix plugin doc 2010-01-17 19:55:08 +01:00
Benoit Jacob
3ab280ce4e add missing semicolon in the example 2010-01-17 12:40:19 -05:00
Benoit Jacob
b40030753b Added tag 2.0.11 for changeset 5f73a8df20 2010-01-10 11:30:40 -05:00
Benoit Jacob
5f73a8df20 bump 2010-01-10 11:30:10 -05:00
Thomas Capricelli
8a6d5f10dc backport from tip : actually stop on compile failure 2010-01-06 17:17:40 +01:00
Benoit Jacob
ba6ed5fa5f Fix CoeffReadCost in Part: it must account for the cost of the
conditional jump. This makes Part considered an "expensive" xpr
 that must be evaluated in operations such as Product.
This fixes bug #80.
2010-01-02 13:04:04 -05:00
Benoit Jacob
e4c88c14ec clarify docs as requested on the forum 2010-01-02 12:54:55 -05:00
Benoit Jacob
74207a31fa backport the fix to bug #79, and the unit test 2010-01-02 12:45:49 -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
Gael Guennebaud
b581cb870c fix #74: sparse triangular solver for lower/row-major matrices 2009-12-14 10:20:35 +01:00
Gael Guennebaud
72fc81dd9d backport quaternion slerp precision fix 2009-12-05 18:28:17 +01:00
Gael Guennebaud
f36650b00a fix MSVC10 compilation issues 2009-12-02 19:34:37 +01: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
a161a70696 Added tag 2.0.10 for changeset 8f1ce52e76 2009-11-25 08:54:17 -05:00
Benoit Jacob
8f1ce52e76 bump 2009-11-25 08:46:42 -05:00
Benoit Jacob
268df314f1 make the inverse_4x4 test pass more consistently 2009-11-25 08:43:20 -05:00
Benoit Jacob
522022ebfc wow, restore Gael's changeset 5455d6fbe8
that I had accidentally undone in my changeset c64ca6870e
.
2009-11-25 08:31:25 -05:00
Thomas Capricelli
d048d7e712 fix bitbucket url after recent change 2009-11-24 23:08:16 +01:00
Benoit Jacob
cd3c8a9404 typo 2009-11-23 11:23:30 -05:00
Benoit Jacob
ec8f37ac75 improve precision test 2009-11-23 11:20:48 -05:00
Benoit Jacob
fc7f39980c backport improvement in 4x4 inverse precisio, and rigorous precision test 2009-11-23 10:27:10 -05:00
Gael Guennebaud
70af59c455 an attempt to fix compilation with recent MSVC 2009-11-23 10:29:40 +01:00
Benoit Jacob
f4dd399499 fix warnings 2009-11-16 14:15:47 -05:00
Benoit Jacob
153557527e backport: init-by-zero option: resize with same size must be a NOP 2009-11-16 13:47:02 -05:00
Benoit Jacob
6aad7f80ff avoid infinite loop, optimization not important, so a plain for loop is the safe way 2009-11-12 14:09:53 -05:00
Benoit Jacob
e3f6c3115a backport the initialize-by-0 option 2009-11-12 12:53:24 -05:00
Benoit Jacob
a2c838ff8f fix PowerPC platform detection 2009-11-11 10:52:00 -05:00
Thomas Capricelli
1e2f56c35a backport of b53c2fcc99
: fix install dir for *.pc

Ingmar Vanhassel <ingmar@exherbo.org>:
Packages that don't install architecture-specific files should install
their pkg-config file to datadir, not libdir.
2009-11-11 15:35:12 +01:00
Hauke Heibel
808c4e9581 Fixed the packport of 62 - Packet4f/d/i does not exist in 2.0. 2009-11-05 10:49:49 +01:00
Hauke Heibel
65331c3884 backporting3979f6d8aad001174160774b49b747430a7686b5
: fixed bug #62
2009-11-04 17:49:34 +01:00
Benoit Jacob
e158cdd61d fix Matrix::Map/MapAligned documentation, and rephrase the tutorial on Map 2009-10-31 14:45:50 -04:00
Benoit Jacob
c64ca6870e this explicit keyword can't hurt 2009-10-31 11:49:20 -04:00
Benoit Jacob
6a90f6c5f0 * default MatrixBase ctor: make it protected, make it a static assert, only do the check when debugging eigen to avoid slowing down compilation for everybody (this check is paranoiac, it's very seldom useful)
* add private MatrixBase ctors to catch cases when the user tries to construct MatrixBase objects directly
2009-10-31 11:48:33 -04:00
Gael Guennebaud
22dd13fdb9 backporting fix of #65 2009-10-29 14:26:38 +01:00
Gael Guennebaud
5455d6fbe8 backporting fix of #65 2009-10-29 14:26:00 +01: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
Benoit Jacob
3d365a75cd Added tag 2.0.9 for changeset 38bc82a6f7 2009-10-24 19:38:58 -04:00
Benoit Jacob
38bc82a6f7 bump 2009-10-24 16:35:46 -04:00
Benoit Jacob
6173eb67ff really fix pkgconfig support and make install.
* mistake: was using the install dir instead of binary dir
* was also using INCLUDE_INSTALL_DIR before it was set, so on a first cmake run, the pkgconfig file was bad
2009-10-24 16:16:48 -04:00
Benoit Jacob
9f89431cea install NewStdVector 2009-10-23 19:58:37 -04:00
Benoit Jacob
79e392472a Added tag 2.0.8 for changeset e1c96f3fe0 2009-10-23 18:47:33 -04:00
Benoit Jacob
e1c96f3fe0 bump 2009-10-23 18:37:05 -04:00
Benoit Jacob
46f0fe3b4b SVD:
* implement default ctor, users were relying on the compiler generater one and reported issues on the forum
* turns out that we crash on 1x1 matrices but work on Nx1. No interest in fixing that, so just guard by assert and unit test.
2009-10-23 18:05:55 -04:00
Benoit Jacob
e17e4f3654 just this is enough 2009-10-23 17:51:32 -04:00
Benoit Jacob
2b006ae430 fix "make install"
aaargh
my shiny birthday release, it's broken already!
2009-10-23 17:47:12 -04:00
Benoit Jacob
df6923fd2b Added tag 2.0.7 for changeset 21d081c6da 2009-10-22 16:15:15 -04:00
Benoit Jacob
21d081c6da bump 2009-10-22 16:14:03 -04:00
Benoit Jacob
81fd1e9060 support gcc 3.3 2009-10-22 15:53:23 -04:00
Benoit Jacob
be8ae0d45f * CMakeLists: only pass -Wextra if it's supported (it's not on gcc 3.3)
* MapBase: put static first (fix gcc 3.3 warning)
* StdVector: add missing newline at end
2009-10-22 15:20:02 -04:00
Hauke Heibel
fc8b54c142 Code cleanup. 2009-10-22 20:06:05 +02:00
Hauke Heibel
76d578fb99 This does fix #61. See the comments in #61 for details. 2009-10-22 09:29:59 +02:00
Benoit Jacob
9af431e78e tentative fix for bug #61 2009-10-21 18:55:54 -04:00
Hauke Heibel
a4a3e511d0 Added missing resize case for m_outerSize==0. 2009-10-15 23:30:15 +02:00
Hauke Heibel
ffee27bf72 Fixed uninitialized variables in unit tests.
Fixed /W4 warning C4512 (LU was left out on purpose).
2009-10-14 08:33:36 +02:00
Benoit Jacob
5e24fbbead add assert for M>=N 2009-10-13 14:17:25 -04:00
Benoit Jacob
09364c8d05 fix compilation with gcc 3.3 2009-10-12 12:29:07 -04:00
Gael Guennebaud
3b8938ee1a backport d97d307fcf 2009-10-06 10:36:59 +02:00
Gael Guennebaud
e43d630d80 fix #10: the reallocateSparse function was half coded
(transplanted from 55de162cf6
)
2009-06-08 14:05:23 +02:00
Thomas Capricelli
eb1df142a3 backport changes in tip related to eigen_gen_docs 2009-10-04 03:38:13 +02:00
Thomas Capricelli
746d8b7ce9 update URL for adol-c (backport from tip) 2009-09-27 02:00:45 +02:00
Benoit Jacob
656c8faeb8 merge 2009-09-25 09:09:49 -04:00
Benoit Jacob
8084dbc86a copy the Memory.h file from the devel branch and remove some added trailing spaces.
This is now very harmless to do as the big change (EIGEN_ALIGN preprocessor stuff and the body of ei_aligned_malloc) was already introduced in 2.0.6.

Should address Björn's issue, and also improve FreeBSD platform detection.
2009-09-25 09:09:14 -04:00
Thomas Capricelli
79ebba4f52 clean tags... 2.0.6 was tagged three times, with two different values.
The best way to "push" a tag is to edit the .hgtags instead of using 'hg
tag xx' several times with the same value 'xx'.
2009-09-25 03:19:28 +02:00
Benoit Jacob
b362b45cff update .hgignore to ignore files created by eigen_gen_credits 2009-09-24 07:06:31 -04:00
Rhys Ulerich
c67b8b7ce0 Added pkgconfig support 2009-09-23 22:05:46 -04:00
Benoit Jacob
bcd621fcd5 Added tag 2.0.6 for changeset de88fb67d6 2009-09-23 12:11:26 -04:00
68 changed files with 949 additions and 313 deletions

View File

@@ -20,4 +20,5 @@ CMakeCache.txt
tags
.*.swp
activity.png
gmon.out
*.out
*.php*

View File

@@ -1,7 +1,13 @@
project(Eigen)
cmake_minimum_required(VERSION 2.6.2)
set(EIGEN_VERSION_NUMBER "2.0.6")
set(INCLUDE_INSTALL_DIR
"${CMAKE_INSTALL_PREFIX}/include/eigen2"
CACHE PATH
"The directory where we install the header files"
FORCE)
set(EIGEN_VERSION_NUMBER "2.0.12")
set(EIGEN_VERSION "${EIGEN_VERSION_NUMBER}")
set(CMAKE_MODULE_PATH ${PROJECT_SOURCE_DIR}/cmake)
@@ -12,6 +18,9 @@ if(NOT WIN32)
option(EIGEN_BUILD_LIB "Build the binary shared library" OFF)
endif(NOT WIN32)
option(EIGEN_BUILD_BTL "Build benchmark suite" OFF)
if(NOT WIN32)
option(EIGEN_BUILD_PKGCONFIG "Build pkg-config .pc file for Eigen" ON)
endif(NOT WIN32)
if(EIGEN_BUILD_LIB)
option(EIGEN_TEST_LIB "Build the unit tests using the library (disable -pedantic)" OFF)
@@ -21,7 +30,12 @@ set(CMAKE_INCLUDE_CURRENT_DIR ON)
if(CMAKE_COMPILER_IS_GNUCXX)
if(CMAKE_SYSTEM_NAME MATCHES Linux)
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wnon-virtual-dtor -Wno-long-long -ansi -Wundef -Wcast-align -Wchar-subscripts -Wall -W -Wpointer-arith -Wwrite-strings -Wformat-security -Wextra -fno-exceptions -fno-check-new -fno-common -fstrict-aliasing")
include(CheckCXXCompilerFlag)
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wnon-virtual-dtor -Wno-long-long -ansi -Wundef -Wcast-align -Wchar-subscripts -Wall -W -Wpointer-arith -Wwrite-strings -Wformat-security -fno-exceptions -fno-check-new -fno-common -fstrict-aliasing")
check_cxx_compiler_flag("-Wextra" has_wextra)
if(has_wextra)
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wextra")
endif()
if(NOT EIGEN_TEST_LIB)
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -pedantic")
endif(NOT EIGEN_TEST_LIB)
@@ -71,6 +85,13 @@ endif(EIGEN_TEST_NO_EXPLICIT_VECTORIZATION)
include_directories(${CMAKE_CURRENT_SOURCE_DIR} ${CMAKE_CURRENT_BINARY_DIR})
if(EIGEN_BUILD_PKGCONFIG)
configure_file(eigen2.pc.in eigen2.pc) # uses INCLUDE_INSTALL_DIR
install(FILES ${CMAKE_CURRENT_BINARY_DIR}/eigen2.pc
DESTINATION share/pkgconfig
)
endif(EIGEN_BUILD_PKGCONFIG)
add_subdirectory(Eigen)
add_subdirectory(unsupported)

View File

@@ -1,4 +1,6 @@
set(Eigen_HEADERS Core LU Cholesky QR Geometry Sparse Array SVD LeastSquares QtAlignedMalloc StdVector)
set(Eigen_HEADERS Core LU Cholesky QR Geometry
Sparse Array SVD LeastSquares
QtAlignedMalloc StdVector NewStdVector)
if(EIGEN_BUILD_LIB)
set(Eigen_SRCS
@@ -20,12 +22,6 @@ if(CMAKE_COMPILER_IS_GNUCXX)
set(CMAKE_CXX_FLAGS_RELWITHDEBINFO "${CMAKE_CXX_FLAGS_RELWITHDEBINFO} -g1 -O2")
endif(CMAKE_COMPILER_IS_GNUCXX)
set(INCLUDE_INSTALL_DIR
"${CMAKE_INSTALL_PREFIX}/include/eigen2"
CACHE PATH
"The directory where we install the header files"
FORCE)
install(FILES
${Eigen_HEADERS}
DESTINATION ${INCLUDE_INSTALL_DIR}/Eigen

View File

@@ -26,7 +26,7 @@
#define EIGEN_SSE2_BUT_NOT_OLD_GCC
#endif
#ifndef EIGEN_DONT_VECTORIZE
#if !defined(EIGEN_DONT_VECTORIZE) && !defined(EIGEN_DONT_ALIGN)
#if defined (EIGEN_SSE2_BUT_NOT_OLD_GCC) || defined(EIGEN_SSE2_ON_MSVC_2008_OR_LATER)
#define EIGEN_VECTORIZE
#define EIGEN_VECTORIZE_SSE

View File

@@ -1,10 +1,23 @@
#ifndef EIGEN_QTMALLOC_MODULE_H
#define EIGEN_QTMALLOC_MODULE_H
#if (!EIGEN_MALLOC_ALREADY_ALIGNED)
#ifdef QVECTOR_H
#error You must include <Eigen/QtAlignedMalloc> before <QtCore/QVector>.
#endif
#ifdef Q_DECL_IMPORT
#define Q_DECL_IMPORT_ORIG Q_DECL_IMPORT
#undef Q_DECL_IMPORT
#define Q_DECL_IMPORT
#else
#define Q_DECL_IMPORT
#endif
#include "Core"
#if (!EIGEN_MALLOC_ALREADY_ALIGNED)
#include <QtCore/QVector>
inline void *qMalloc(size_t size)
{
@@ -26,4 +39,11 @@ inline void *qRealloc(void *ptr, size_t size)
#endif
#ifdef Q_DECL_IMPORT_ORIG
#define Q_DECL_IMPORT Q_DECL_IMPORT_ORIG
#undef Q_DECL_IMPORT_ORIG
#else
#undef Q_DECL_IMPORT
#endif
#endif // EIGEN_QTMALLOC_MODULE_H

View File

@@ -144,4 +144,4 @@ class vector<T,DummyAlloc,true>
#endif // EIGEN_STDVECTOR_MODULE_H
#endif // EIGEN_USE_NEW_STDVECTOR
#endif // EIGEN_USE_NEW_STDVECTOR

View File

@@ -43,6 +43,8 @@ struct ei_scalar_add_op {
inline const PacketScalar packetOp(const PacketScalar& a) const
{ return ei_padd(a, ei_pset1(m_other)); }
const Scalar m_other;
private:
ei_scalar_add_op& operator=(const ei_scalar_add_op&);
};
template<typename Scalar>
struct ei_functor_traits<ei_scalar_add_op<Scalar> >
@@ -138,6 +140,8 @@ struct ei_scalar_pow_op {
inline ei_scalar_pow_op(const Scalar& exponent) : m_exponent(exponent) {}
inline Scalar operator() (const Scalar& a) const { return ei_pow(a, m_exponent); }
const Scalar m_exponent;
private:
ei_scalar_pow_op& operator=(const ei_scalar_pow_op&);
};
template<typename Scalar>
struct ei_functor_traits<ei_scalar_pow_op<Scalar> >

View File

@@ -133,6 +133,8 @@ struct ei_member_redux {
inline result_type operator()(const MatrixBase<Derived>& mat) const
{ return mat.redux(m_functor); }
const BinaryOp m_functor;
private:
ei_member_redux& operator=(const ei_member_redux&);
};
/** \array_module \ingroup Array
@@ -290,6 +292,9 @@ template<typename ExpressionType, int Direction> class PartialRedux
protected:
ExpressionTypeNested m_matrix;
private:
PartialRedux& operator=(const PartialRedux&);
};
/** \array_module

View File

@@ -222,7 +222,7 @@ class Block<MatrixType,BlockRows,BlockCols,PacketAccess,HasDirectAccess>
class InnerIterator;
typedef typename ei_traits<Block>::AlignedDerivedType AlignedDerivedType;
friend class Block<MatrixType,BlockRows,BlockCols,PacketAccess==AsRequested?ForceAligned:AsRequested,HasDirectAccess>;
friend class Block<MatrixType,BlockRows,BlockCols,PacketAccess==int(AsRequested)?ForceAligned:AsRequested,HasDirectAccess>;
EIGEN_INHERIT_ASSIGNMENT_OPERATORS(Block)

View File

@@ -116,6 +116,9 @@ struct CommaInitializer
int m_row; // current row id
int m_col; // current col id
int m_currentBlockRows; // current block height
private:
CommaInitializer& operator=(const CommaInitializer&);
};
/** \anchor MatrixBaseCommaInitRef

View File

@@ -178,6 +178,9 @@ template<typename ExpressionType> class Cwise
protected:
ExpressionTypeNested m_matrix;
private:
Cwise& operator=(const Cwise&);
};
/** \returns a Cwise wrapper of *this providing additional coefficient-wise operations

View File

@@ -47,11 +47,11 @@ struct ei_traits<DiagonalCoeffs<MatrixType> >
typedef typename ei_unref<MatrixTypeNested>::type _MatrixTypeNested;
enum {
RowsAtCompileTime = int(MatrixType::SizeAtCompileTime) == Dynamic ? Dynamic
: EIGEN_ENUM_MIN(MatrixType::RowsAtCompileTime,
: EIGEN_SIZE_MIN(MatrixType::RowsAtCompileTime,
MatrixType::ColsAtCompileTime),
ColsAtCompileTime = 1,
MaxRowsAtCompileTime = int(MatrixType::MaxSizeAtCompileTime) == Dynamic ? Dynamic
: EIGEN_ENUM_MIN(MatrixType::MaxRowsAtCompileTime,
: EIGEN_SIZE_MIN(MatrixType::MaxRowsAtCompileTime,
MatrixType::MaxColsAtCompileTime),
MaxColsAtCompileTime = 1,
Flags = (unsigned int)_MatrixTypeNested::Flags & (HereditaryBits | LinearAccessBit),

View File

@@ -109,6 +109,9 @@ template<typename ExpressionType, unsigned int Added, unsigned int Removed> clas
protected:
ExpressionTypeNested m_matrix;
private:
Flagged& operator=(const Flagged&);
};
/** \returns an expression of *this with added flags

View File

@@ -279,6 +279,8 @@ struct ei_scalar_multiple_op {
EIGEN_STRONG_INLINE const PacketScalar packetOp(const PacketScalar& a) const
{ return ei_pmul(a, ei_pset1(m_other)); }
const Scalar m_other;
private:
ei_scalar_multiple_op& operator=(const ei_scalar_multiple_op&);
};
template<typename Scalar>
struct ei_functor_traits<ei_scalar_multiple_op<Scalar> >
@@ -294,6 +296,8 @@ struct ei_scalar_quotient1_impl {
EIGEN_STRONG_INLINE const PacketScalar packetOp(const PacketScalar& a) const
{ return ei_pmul(a, ei_pset1(m_other)); }
const Scalar m_other;
private:
ei_scalar_quotient1_impl& operator=(const ei_scalar_quotient1_impl&);
};
template<typename Scalar>
struct ei_functor_traits<ei_scalar_quotient1_impl<Scalar,true> >
@@ -306,6 +310,8 @@ struct ei_scalar_quotient1_impl<Scalar,false> {
EIGEN_STRONG_INLINE ei_scalar_quotient1_impl(const Scalar& other) : m_other(other) {}
EIGEN_STRONG_INLINE Scalar operator() (const Scalar& a) const { return a / m_other; }
const Scalar m_other;
private:
ei_scalar_quotient1_impl& operator=(const ei_scalar_quotient1_impl&);
};
template<typename Scalar>
struct ei_functor_traits<ei_scalar_quotient1_impl<Scalar,false> >
@@ -323,6 +329,8 @@ template<typename Scalar>
struct ei_scalar_quotient1_op : ei_scalar_quotient1_impl<Scalar, NumTraits<Scalar>::HasFloatingPoint > {
EIGEN_STRONG_INLINE ei_scalar_quotient1_op(const Scalar& other)
: ei_scalar_quotient1_impl<Scalar, NumTraits<Scalar>::HasFloatingPoint >(other) {}
private:
ei_scalar_quotient1_op& operator=(const ei_scalar_quotient1_op&);
};
// nullary functors
@@ -335,6 +343,8 @@ struct ei_scalar_constant_op {
EIGEN_STRONG_INLINE const Scalar operator() (int, int = 0) const { return m_other; }
EIGEN_STRONG_INLINE const PacketScalar packetOp() const { return ei_pset1(m_other); }
const Scalar m_other;
private:
ei_scalar_constant_op& operator=(const ei_scalar_constant_op&);
};
template<typename Scalar>
struct ei_functor_traits<ei_scalar_constant_op<Scalar> >

View File

@@ -66,11 +66,11 @@ template<typename Derived> class MapBase
inline const Scalar* data() const { return m_data; }
template<bool IsForceAligned,typename Dummy> struct force_aligned_impl {
AlignedDerivedType static run(MapBase& a) { return a.derived(); }
static AlignedDerivedType run(MapBase& a) { return a.derived(); }
};
template<typename Dummy> struct force_aligned_impl<false,Dummy> {
AlignedDerivedType static run(MapBase& a) { return a.derived()._convertToForceAligned(); }
static AlignedDerivedType run(MapBase& a) { return a.derived()._convertToForceAligned(); }
};
/** \returns an expression equivalent to \c *this but having the \c PacketAccess constant
@@ -99,7 +99,7 @@ template<typename Derived> class MapBase
inline const Scalar coeff(int index) const
{
ei_assert(Derived::IsVectorAtCompileTime || (ei_traits<Derived>::Flags & LinearAccessBit));
if ( ((RowsAtCompileTime == 1) == IsRowMajor) )
if ( ((RowsAtCompileTime == 1) == IsRowMajor) || !int(Derived::IsVectorAtCompileTime) )
return m_data[index];
else
return m_data[index*stride()];
@@ -108,7 +108,7 @@ template<typename Derived> class MapBase
inline Scalar& coeffRef(int index)
{
ei_assert(Derived::IsVectorAtCompileTime || (ei_traits<Derived>::Flags & LinearAccessBit));
if ( ((RowsAtCompileTime == 1) == IsRowMajor) )
if ( ((RowsAtCompileTime == 1) == IsRowMajor) || !int(Derived::IsVectorAtCompileTime) )
return const_cast<Scalar*>(m_data)[index];
else
return const_cast<Scalar*>(m_data)[index*stride()];

View File

@@ -54,7 +54,7 @@ template<> inline int machine_epsilon<int>() { return 0; }
inline int ei_real(int x) { return x; }
inline int ei_imag(int) { return 0; }
inline int ei_conj(int x) { return x; }
inline int ei_abs(int x) { return abs(x); }
inline int ei_abs(int x) { return std::abs(x); }
inline int ei_abs2(int x) { return x*x; }
inline int ei_sqrt(int) { ei_assert(false); return 0; }
inline int ei_exp(int) { ei_assert(false); return 0; }
@@ -67,7 +67,7 @@ inline int ei_pow(int x, int y) { return int(std::pow(double(x), y)); }
template<> inline int ei_random(int a, int b)
{
// We can't just do rand()%n as only the high-order bits are really random
return a + static_cast<int>((b-a+1) * (rand() / (RAND_MAX + 1.0)));
return a + static_cast<int>((b-a+1) * (std::rand() / (RAND_MAX + 1.0)));
}
template<> inline int ei_random()
{

View File

@@ -25,6 +25,11 @@
#ifndef EIGEN_MATRIX_H
#define EIGEN_MATRIX_H
#ifdef EIGEN_INITIALIZE_MATRICES_BY_ZERO
# define EIGEN_INITIALIZE_BY_ZERO_IF_THAT_OPTION_IS_ENABLED for(int i=0;i<base().size();++i) coeffRef(i)=Scalar(0);
#else
# define EIGEN_INITIALIZE_BY_ZERO_IF_THAT_OPTION_IS_ENABLED
#endif
/** \class Matrix
*
@@ -233,7 +238,14 @@ class Matrix
&& (RowsAtCompileTime == Dynamic || RowsAtCompileTime == rows)
&& (MaxColsAtCompileTime == Dynamic || MaxColsAtCompileTime >= cols)
&& (ColsAtCompileTime == Dynamic || ColsAtCompileTime == cols));
m_storage.resize(rows * cols, rows, cols);
#ifdef EIGEN_INITIALIZE_MATRICES_BY_ZERO
int size = rows*cols;
bool size_changed = size != this->size();
m_storage.resize(size, rows, cols);
if(size_changed) EIGEN_INITIALIZE_BY_ZERO_IF_THAT_OPTION_IS_ENABLED
#else
m_storage.resize(rows*cols, rows, cols);
#endif
}
/** Resizes \c *this to a vector of length \a size
@@ -243,10 +255,17 @@ class Matrix
inline void resize(int size)
{
EIGEN_STATIC_ASSERT_VECTOR_ONLY(Matrix)
ei_assert(SizeAtCompileTime == Dynamic || SizeAtCompileTime == size);
#ifdef EIGEN_INITIALIZE_MATRICES_BY_ZERO
bool size_changed = size != this->size();
#endif
if(RowsAtCompileTime == 1)
m_storage.resize(size, 1, size);
else
m_storage.resize(size, size, 1);
#ifdef EIGEN_INITIALIZE_MATRICES_BY_ZERO
if(size_changed) EIGEN_INITIALIZE_BY_ZERO_IF_THAT_OPTION_IS_ENABLED
#endif
}
/** Copies the value of the expression \a other into \c *this with automatic resizing.
@@ -290,13 +309,14 @@ class Matrix
EIGEN_STRONG_INLINE explicit Matrix() : m_storage()
{
_check_template_params();
EIGEN_INITIALIZE_BY_ZERO_IF_THAT_OPTION_IS_ENABLED
}
#ifndef EIGEN_PARSED_BY_DOXYGEN
/** \internal */
Matrix(ei_constructor_without_unaligned_array_assert)
: m_storage(ei_constructor_without_unaligned_array_assert())
{}
{EIGEN_INITIALIZE_BY_ZERO_IF_THAT_OPTION_IS_ENABLED}
#endif
/** Constructs a vector or row-vector with given dimension. \only_for_vectors
@@ -312,6 +332,7 @@ class Matrix
EIGEN_STATIC_ASSERT_VECTOR_ONLY(Matrix)
ei_assert(dim > 0);
ei_assert(SizeAtCompileTime == Dynamic || SizeAtCompileTime == dim);
EIGEN_INITIALIZE_BY_ZERO_IF_THAT_OPTION_IS_ENABLED
}
/** This constructor has two very different behaviors, depending on the type of *this.
@@ -337,6 +358,7 @@ class Matrix
{
ei_assert(x > 0 && (RowsAtCompileTime == Dynamic || RowsAtCompileTime == x)
&& y > 0 && (ColsAtCompileTime == Dynamic || ColsAtCompileTime == y));
EIGEN_INITIALIZE_BY_ZERO_IF_THAT_OPTION_IS_ENABLED
}
}
/** constructs an initialized 2D vector with given coefficients */
@@ -402,9 +424,10 @@ class Matrix
void swap(const MatrixBase<OtherDerived>& other);
/** \name Map
* These are convenience functions returning Map objects. The Map() static functions return unaligned Map objects,
* while the AlignedMap() functions return aligned Map objects and thus should be called only with 16-byte-aligned
* \a data pointers.
* These are convenience functions returning Map objects.
*
* \warning Do not use MapAligned in the Eigen 2.0. Mapping aligned arrays will be fully
* supported in Eigen 3.0 (already implemented in the development branch)
*
* \see class Map
*/
@@ -505,7 +528,9 @@ class Matrix
template<typename OtherDerived>
EIGEN_STRONG_INLINE Matrix& _set(const MatrixBase<OtherDerived>& other)
{
_set_selector(other.derived(), typename ei_meta_if<bool(int(OtherDerived::Flags) & EvalBeforeAssigningBit), ei_meta_true, ei_meta_false>::ret());
// this enum introduced to fix compilation with gcc 3.3
enum { cond = int(OtherDerived::Flags) & EvalBeforeAssigningBit };
_set_selector(other.derived(), typename ei_meta_if<bool(cond), ei_meta_true, ei_meta_false>::ret());
return *this;
}
@@ -567,7 +592,8 @@ template<typename _Scalar, int _Rows, int _Cols, int _Options, int _MaxRows, int
template<typename OtherDerived>
inline void Matrix<_Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols>::swap(const MatrixBase<OtherDerived>& other)
{
ei_matrix_swap_impl<Matrix, OtherDerived>::run(*this, *const_cast<MatrixBase<OtherDerived>*>(&other));
// the Eigen:: here is to work around a stupid ICC 11.1 bug.
Eigen::ei_matrix_swap_impl<Matrix, OtherDerived>::run(*this, *const_cast<MatrixBase<OtherDerived>*>(&other));
}

View File

@@ -136,12 +136,6 @@ template<typename Derived> class MatrixBase
*/
};
/** Default constructor. Just checks at compile-time for self-consistency of the flags. */
MatrixBase()
{
ei_assert(ei_are_flags_consistent<Flags>::ret);
}
#ifndef EIGEN_PARSED_BY_DOXYGEN
/** This is the "real scalar" type; if the \a Scalar type is already real numbers
* (e.g. int, float or double) then \a RealScalar is just the same as \a Scalar. If
@@ -165,7 +159,7 @@ template<typename Derived> class MatrixBase
inline int size() const { return rows() * cols(); }
/** \returns the number of nonzero coefficients which is in practice the number
* of stored coefficients. */
inline int nonZeros() const { return derived.nonZeros(); }
inline int nonZeros() const { return size(); }
/** \returns true if either the number of rows or the number of columns is equal to 1.
* In other words, this function returns
* \code rows()==1 || cols()==1 \endcode
@@ -589,7 +583,8 @@ template<typename Derived> class MatrixBase
const LU<PlainMatrixType> lu() const;
const PlainMatrixType inverse() const;
void computeInverse(PlainMatrixType *result) const;
template<typename ResultType>
void computeInverse(MatrixBase<ResultType> *result) const;
Scalar determinant() const;
/////////// Cholesky module ///////////
@@ -627,6 +622,24 @@ template<typename Derived> class MatrixBase
#ifdef EIGEN_MATRIXBASE_PLUGIN
#include EIGEN_MATRIXBASE_PLUGIN
#endif
protected:
/** Default constructor. Do nothing. */
MatrixBase()
{
/* Just checks for self-consistency of the flags.
* Only do it when debugging Eigen, as this borders on paranoiac and could slow compilation down
*/
#ifdef EIGEN_INTERNAL_DEBUGGING
EIGEN_STATIC_ASSERT(ei_are_flags_consistent<Flags>::ret,
INVALID_MATRIXBASE_TEMPLATE_PARAMETERS)
#endif
}
private:
explicit MatrixBase(int);
MatrixBase(int,int);
template<typename OtherDerived> explicit MatrixBase(const MatrixBase<OtherDerived>&);
};
#endif // EIGEN_MATRIXBASE_H

View File

@@ -40,7 +40,7 @@ template <typename T, int Size, int MatrixOptions,
ei_matrix_array()
{
#ifndef EIGEN_DISABLE_UNALIGNED_ARRAY_ASSERT
ei_assert((reinterpret_cast<size_t>(array) & 0xf) == 0
ei_assert((reinterpret_cast<std::size_t>(array) & 0xf) == 0
&& "this assertion is explained here: http://eigen.tuxfamily.org/dox/UnalignedArrayAssert.html **** READ THIS WEB PAGE !!! ****");
#endif
}

View File

@@ -100,6 +100,9 @@ template<typename ExpressionType> class NestByValue
protected:
const ExpressionType m_expression;
private:
NestByValue& operator=(const NestByValue&);
};
/** \returns an expression of the temporary version of *this.

View File

@@ -50,7 +50,7 @@ struct ei_traits<Part<MatrixType, Mode> > : ei_traits<MatrixType>
typedef typename ei_unref<MatrixTypeNested>::type _MatrixTypeNested;
enum {
Flags = (_MatrixTypeNested::Flags & (HereditaryBits) & (~(PacketAccessBit | DirectAccessBit | LinearAccessBit))) | Mode,
CoeffReadCost = _MatrixTypeNested::CoeffReadCost
CoeffReadCost = _MatrixTypeNested::CoeffReadCost + ConditionalJumpCost
};
};
@@ -124,8 +124,10 @@ template<typename MatrixType, unsigned int Mode> class Part
}
protected:
const typename MatrixType::Nested m_matrix;
private:
Part& operator=(const Part&);
};
/** \nonstableyet

View File

@@ -66,11 +66,8 @@ struct ProductReturnType
template<typename Lhs, typename Rhs>
struct ProductReturnType<Lhs,Rhs,CacheFriendlyProduct>
{
typedef typename ei_nested<Lhs,Rhs::ColsAtCompileTime>::type LhsNested;
typedef typename ei_nested<Rhs,Lhs::RowsAtCompileTime,
typename ei_plain_matrix_type_column_major<Rhs>::type
>::type RhsNested;
typedef const Lhs& LhsNested;
typedef const Rhs& RhsNested;
typedef Product<LhsNested, RhsNested, CacheFriendlyProduct> Type;
};
@@ -128,7 +125,7 @@ struct ei_traits<Product<LhsNested, RhsNested, ProductMode> >
RowsAtCompileTime = _LhsNested::RowsAtCompileTime,
ColsAtCompileTime = _RhsNested::ColsAtCompileTime,
InnerSize = EIGEN_ENUM_MIN(_LhsNested::ColsAtCompileTime, _RhsNested::RowsAtCompileTime),
InnerSize = EIGEN_SIZE_MIN(_LhsNested::ColsAtCompileTime, _RhsNested::RowsAtCompileTime),
MaxRowsAtCompileTime = _LhsNested::MaxRowsAtCompileTime,
MaxColsAtCompileTime = _RhsNested::MaxColsAtCompileTime,
@@ -144,7 +141,7 @@ struct ei_traits<Product<LhsNested, RhsNested, ProductMode> >
EvalToRowMajor = RhsRowMajor && (ProductMode==(int)CacheFriendlyProduct ? LhsRowMajor : (!CanVectorizeLhs)),
RemovedBits = ~(EvalToRowMajor ? 0 : RowMajorBit),
RemovedBits = ~((EvalToRowMajor ? 0 : RowMajorBit)|DirectAccessBit),
Flags = ((unsigned int)(LhsFlags | RhsFlags) & HereditaryBits & RemovedBits)
| EvalBeforeAssigningBit
@@ -571,7 +568,7 @@ struct ei_cache_friendly_product_selector<ProductType,LhsRows,ColMajor,HasDirect
else
{
_res = ei_aligned_stack_new(Scalar,res.size());
Map<Matrix<Scalar,DestDerived::RowsAtCompileTime,1> >(_res, res.size()) = res;
Map<Matrix<Scalar,DestDerived::RowsAtCompileTime,1,ColMajor> >(_res, res.size()) = res;
}
ei_cache_friendly_product_colmajor_times_vector(res.size(),
&product.lhs().const_cast_derived().coeffRef(0,0), product.lhs().stride(),
@@ -579,7 +576,7 @@ struct ei_cache_friendly_product_selector<ProductType,LhsRows,ColMajor,HasDirect
if (!EvalToRes)
{
res = Map<Matrix<Scalar,DestDerived::SizeAtCompileTime,1> >(_res, res.size());
res = Map<Matrix<Scalar,DestDerived::SizeAtCompileTime,1,ColMajor> >(_res, res.size());
ei_aligned_stack_delete(Scalar, _res, res.size());
}
}
@@ -617,7 +614,7 @@ struct ei_cache_friendly_product_selector<ProductType,1,LhsOrder,LhsAccess,RhsCo
else
{
_res = ei_aligned_stack_new(Scalar, res.size());
Map<Matrix<Scalar,DestDerived::SizeAtCompileTime,1> >(_res, res.size()) = res;
Map<Matrix<Scalar,DestDerived::SizeAtCompileTime,1,ColMajor> >(_res, res.size()) = res;
}
ei_cache_friendly_product_colmajor_times_vector(res.size(),
&product.rhs().const_cast_derived().coeffRef(0,0), product.rhs().stride(),
@@ -625,7 +622,7 @@ struct ei_cache_friendly_product_selector<ProductType,1,LhsOrder,LhsAccess,RhsCo
if (!EvalToRes)
{
res = Map<Matrix<Scalar,DestDerived::SizeAtCompileTime,1> >(_res, res.size());
res = Map<Matrix<Scalar,DestDerived::SizeAtCompileTime,1,ColMajor> >(_res, res.size());
ei_aligned_stack_delete(Scalar, _res, res.size());
}
}
@@ -650,7 +647,7 @@ struct ei_cache_friendly_product_selector<ProductType,LhsRows,RowMajor,HasDirect
else
{
_rhs = ei_aligned_stack_new(Scalar, product.rhs().size());
Map<Matrix<Scalar,Rhs::SizeAtCompileTime,1> >(_rhs, product.rhs().size()) = product.rhs();
Map<Matrix<Scalar,Rhs::SizeAtCompileTime,1,ColMajor> >(_rhs, product.rhs().size()) = product.rhs();
}
ei_cache_friendly_product_rowmajor_times_vector(&product.lhs().const_cast_derived().coeffRef(0,0), product.lhs().stride(),
_rhs, product.rhs().size(), res);
@@ -678,7 +675,7 @@ struct ei_cache_friendly_product_selector<ProductType,1,LhsOrder,LhsAccess,RhsCo
else
{
_lhs = ei_aligned_stack_new(Scalar, product.lhs().size());
Map<Matrix<Scalar,Lhs::SizeAtCompileTime,1> >(_lhs, product.lhs().size()) = product.lhs();
Map<Matrix<Scalar,Lhs::SizeAtCompileTime,1,ColMajor> >(_lhs, product.lhs().size()) = product.lhs();
}
ei_cache_friendly_product_rowmajor_times_vector(&product.rhs().const_cast_derived().coeffRef(0,0), product.rhs().stride(),
_lhs, product.lhs().size(), res);
@@ -709,7 +706,17 @@ MatrixBase<Derived>::operator+=(const Flagged<Product<Lhs,Rhs,CacheFriendlyProdu
if (other._expression()._useCacheFriendlyProduct())
ei_cache_friendly_product_selector<Product<Lhs,Rhs,CacheFriendlyProduct> >::run(const_cast_derived(), other._expression());
else
lazyAssign(derived() + other._expression());
{
typedef typename ei_cleantype<Lhs>::type _Lhs;
typedef typename ei_cleantype<Rhs>::type _Rhs;
typedef typename ei_nested<_Lhs,_Rhs::ColsAtCompileTime>::type LhsNested;
typedef typename ei_nested<_Rhs,_Lhs::RowsAtCompileTime>::type RhsNested;
Product<LhsNested,RhsNested,NormalProduct> prod(other._expression().lhs(),other._expression().rhs());
lazyAssign(derived() + prod);
}
return derived();
}
@@ -724,12 +731,21 @@ inline Derived& MatrixBase<Derived>::lazyAssign(const Product<Lhs,Rhs,CacheFrien
}
else
{
lazyAssign<Product<Lhs,Rhs,CacheFriendlyProduct> >(product);
typedef typename ei_cleantype<Lhs>::type _Lhs;
typedef typename ei_cleantype<Rhs>::type _Rhs;
typedef typename ei_nested<_Lhs,_Rhs::ColsAtCompileTime>::type LhsNested;
typedef typename ei_nested<_Rhs,_Lhs::RowsAtCompileTime>::type RhsNested;
typedef Product<LhsNested,RhsNested,NormalProduct> NormalProduct;
NormalProduct normal_prod(product.lhs(),product.rhs());
lazyAssign<NormalProduct>(normal_prod);
}
return derived();
}
template<typename T> struct ei_product_copy_rhs
template<typename T,int StorageOrder> struct ei_product_copy_rhs
{
typedef typename ei_meta_if<
(ei_traits<T>::Flags & RowMajorBit)
@@ -739,11 +755,30 @@ template<typename T> struct ei_product_copy_rhs
>::ret type;
};
template<typename T> struct ei_product_copy_lhs
template<typename T> struct ei_product_copy_rhs<T,RowMajorBit>
{
typedef typename ei_meta_if<
(!(ei_traits<T>::Flags & DirectAccessBit)),
typename ei_plain_matrix_type<T>::type,
const T&
>::ret type;
};
template<typename T,int StorageOrder> struct ei_product_copy_lhs
{
typedef typename ei_meta_if<
(!(int(ei_traits<T>::Flags) & DirectAccessBit)),
typename ei_plain_matrix_type<T>::type,
typename ei_plain_matrix_type_row_major<T>::type,
const T&
>::ret type;
};
template<typename T> struct ei_product_copy_lhs<T,RowMajorBit>
{
typedef typename ei_meta_if<
((ei_traits<T>::Flags & RowMajorBit)==0)
|| (!(int(ei_traits<T>::Flags) & DirectAccessBit)),
typename ei_plain_matrix_type_row_major<T>::type,
const T&
>::ret type;
};
@@ -752,9 +787,9 @@ template<typename Lhs, typename Rhs, int ProductMode>
template<typename DestDerived>
inline void Product<Lhs,Rhs,ProductMode>::_cacheFriendlyEvalAndAdd(DestDerived& res) const
{
typedef typename ei_product_copy_lhs<_LhsNested>::type LhsCopy;
typedef typename ei_product_copy_lhs<_LhsNested,DestDerived::Flags&RowMajorBit>::type LhsCopy;
typedef typename ei_unref<LhsCopy>::type _LhsCopy;
typedef typename ei_product_copy_rhs<_RhsNested>::type RhsCopy;
typedef typename ei_product_copy_rhs<_RhsNested,DestDerived::Flags&RowMajorBit>::type RhsCopy;
typedef typename ei_unref<RhsCopy>::type _RhsCopy;
LhsCopy lhs(m_lhs);
RhsCopy rhs(m_rhs);
@@ -764,6 +799,7 @@ inline void Product<Lhs,Rhs,ProductMode>::_cacheFriendlyEvalAndAdd(DestDerived&
_RhsCopy::Flags&RowMajorBit, (const Scalar*)&(rhs.const_cast_derived().coeffRef(0,0)), rhs.stride(),
DestDerived::Flags&RowMajorBit, (Scalar*)&(res.coeffRef(0,0)), res.stride()
);
}
#endif // EIGEN_PRODUCT_H

View File

@@ -117,6 +117,9 @@ template<typename ExpressionType> class SwapWrapper
protected:
ExpressionType& m_expression;
private:
SwapWrapper& operator=(const SwapWrapper&);
};
/** swaps *this with the expression \a other.

View File

@@ -69,7 +69,6 @@ template<typename MatrixType> class Transpose
inline int rows() const { return m_matrix.cols(); }
inline int cols() const { return m_matrix.rows(); }
inline int nonZeros() const { return m_matrix.nonZeros(); }
inline int stride(void) const { return m_matrix.stride(); }
inline Scalar& coeffRef(int row, int col)

View File

@@ -114,9 +114,15 @@ template<> EIGEN_STRONG_INLINE void ei_pstoreu<float>(float* to, const __m128&
template<> EIGEN_STRONG_INLINE void ei_pstoreu<double>(double* to, const __m128d& from) { _mm_storeu_pd(to, from); }
template<> EIGEN_STRONG_INLINE void ei_pstoreu<int>(int* to, const __m128i& from) { _mm_storeu_si128(reinterpret_cast<__m128i*>(to), from); }
#ifdef _MSC_VER
// this fix internal compilation error
template<> EIGEN_STRONG_INLINE float ei_pfirst<__m128>(const __m128& a) { float x = _mm_cvtss_f32(a); return x; }
#if defined(_MSC_VER) && (_MSC_VER <= 1500) && defined(_WIN64)
// The temporary variable fixes an internal compilation error.
// Direct of the struct members fixed bug #62.
template<> EIGEN_STRONG_INLINE float ei_pfirst<__m128>(const __m128& a) { return a.m128_f32[0]; }
template<> EIGEN_STRONG_INLINE double ei_pfirst<__m128d>(const __m128d& a) { return a.m128d_f64[0]; }
template<> EIGEN_STRONG_INLINE int ei_pfirst<__m128i>(const __m128i& a) { int x = _mm_cvtsi128_si32(a); return x; }
#elif defined(_MSC_VER) && (_MSC_VER <= 1500)
// The temporary variable fixes an internal compilation error.
template<> EIGEN_STRONG_INLINE float ei_pfirst<__m128>(const __m128& a) { float x = _mm_cvtss_f32(a); return x; }
template<> EIGEN_STRONG_INLINE double ei_pfirst<__m128d>(const __m128d& a) { double x = _mm_cvtsd_f64(a); return x; }
template<> EIGEN_STRONG_INLINE int ei_pfirst<__m128i>(const __m128i& a) { int x = _mm_cvtsi128_si32(a); return x; }
#else

View File

@@ -30,7 +30,7 @@
#define EIGEN_WORLD_VERSION 2
#define EIGEN_MAJOR_VERSION 0
#define EIGEN_MINOR_VERSION 6
#define EIGEN_MINOR_VERSION 12
#define EIGEN_VERSION_AT_LEAST(x,y,z) (EIGEN_WORLD_VERSION>x || (EIGEN_WORLD_VERSION>=x && \
(EIGEN_MAJOR_VERSION>y || (EIGEN_MAJOR_VERSION>=y && \
@@ -39,7 +39,7 @@
// 16 byte alignment is only useful for vectorization. Since it affects the ABI, we need to enable 16 byte alignment on all
// platforms where vectorization might be enabled. In theory we could always enable alignment, but it can be a cause of problems
// on some platforms, so we just disable it in certain common platform (compiler+architecture combinations) to avoid these problems.
#if defined(__GNUC__) && !(defined(__i386__) || defined(__x86_64__) || defined(__powerpc__) || defined(__ia64__))
#if defined(__GNUC__) && !(defined(__i386__) || defined(__x86_64__) || defined(__powerpc__) || defined(__ppc__) || defined(__ia64__))
#define EIGEN_GCC_AND_ARCH_DOESNT_WANT_ALIGNMENT 1
#else
#define EIGEN_GCC_AND_ARCH_DOESNT_WANT_ALIGNMENT 0
@@ -257,6 +257,9 @@ enum { RowsAtCompileTime = Eigen::ei_traits<Derived>::RowsAtCompileTime, \
_EIGEN_GENERIC_PUBLIC_INTERFACE(Derived, Eigen::MatrixBase<Derived>)
#define EIGEN_ENUM_MIN(a,b) (((int)a <= (int)b) ? (int)a : (int)b)
#define EIGEN_SIZE_MIN(a,b) (((int)a == 1 || (int)b == 1) ? 1 \
: ((int)a == Dynamic || (int)b == Dynamic) ? Dynamic \
: ((int)a <= (int)b) ? (int)a : (int)b)
#define EIGEN_ENUM_MAX(a,b) (((int)a >= (int)b) ? (int)a : (int)b)
// just an empty macro !

View File

@@ -1,5 +1,5 @@
// This file is part of Eigen, a lightweight C++ template library
// for linear algebra. Eigen itself is part of the KDE project.
// for linear algebra.
//
// Copyright (C) 2008 Gael Guennebaud <g.gael@free.fr>
// Copyright (C) 2008-2009 Benoit Jacob <jacob.benoit.1@gmail.com>
@@ -27,7 +27,17 @@
#ifndef EIGEN_MEMORY_H
#define EIGEN_MEMORY_H
#if defined(__APPLE__) || defined(_WIN64)
// FreeBSD 6 seems to have 16-byte aligned malloc
// See http://svn.freebsd.org/viewvc/base/stable/6/lib/libc/stdlib/malloc.c?view=markup
// FreeBSD 7 seems to have 16-byte aligned malloc except on ARM and MIPS architectures
// See http://svn.freebsd.org/viewvc/base/stable/7/lib/libc/stdlib/malloc.c?view=markup
#if defined(__FreeBSD__) && !defined(__arm__) && !defined(__mips__)
#define EIGEN_FREEBSD_MALLOC_ALREADY_ALIGNED 1
#else
#define EIGEN_FREEBSD_MALLOC_ALREADY_ALIGNED 0
#endif
#if defined(__APPLE__) || defined(_WIN64) || EIGEN_FREEBSD_MALLOC_ALREADY_ALIGNED
#define EIGEN_MALLOC_ALREADY_ALIGNED 1
#else
#define EIGEN_MALLOC_ALREADY_ALIGNED 0
@@ -49,10 +59,10 @@
* Fast, but wastes 16 additional bytes of memory.
* Does not throw any exception.
*/
inline void* ei_handmade_aligned_malloc(size_t size)
inline void* ei_handmade_aligned_malloc(std::size_t size)
{
void *original = malloc(size+16);
void *aligned = reinterpret_cast<void*>((reinterpret_cast<size_t>(original) & ~(size_t(15))) + 16);
void *original = std::malloc(size+16);
void *aligned = reinterpret_cast<void*>((reinterpret_cast<std::size_t>(original) & ~(std::size_t(15))) + 16);
*(reinterpret_cast<void**>(aligned) - 1) = original;
return aligned;
}
@@ -61,13 +71,13 @@ inline void* ei_handmade_aligned_malloc(size_t size)
inline void ei_handmade_aligned_free(void *ptr)
{
if(ptr)
free(*(reinterpret_cast<void**>(ptr) - 1));
std::free(*(reinterpret_cast<void**>(ptr) - 1));
}
/** \internal allocates \a size bytes. The returned pointer is guaranteed to have 16 bytes alignment.
* On allocation error, the returned pointer is undefined, but if exceptions are enabled then a std::bad_alloc is thrown.
* On allocation error, the returned pointer is null, and if exceptions are enabled then a std::bad_alloc is thrown.
*/
inline void* ei_aligned_malloc(size_t size)
inline void* ei_aligned_malloc(std::size_t size)
{
#ifdef EIGEN_NO_MALLOC
ei_assert(false && "heap allocation is forbidden (EIGEN_NO_MALLOC is defined)");
@@ -96,40 +106,49 @@ inline void* ei_aligned_malloc(size_t size)
}
/** allocates \a size bytes. If Align is true, then the returned ptr is 16-byte-aligned.
* On allocation error, the returned pointer is undefined, but if exceptions are enabled then a std::bad_alloc is thrown.
* On allocation error, the returned pointer is null, and if exceptions are enabled then a std::bad_alloc is thrown.
*/
template<bool Align> inline void* ei_conditional_aligned_malloc(size_t size)
template<bool Align> inline void* ei_conditional_aligned_malloc(std::size_t size)
{
return ei_aligned_malloc(size);
}
template<> inline void* ei_conditional_aligned_malloc<false>(size_t size)
template<> inline void* ei_conditional_aligned_malloc<false>(std::size_t size)
{
#ifdef EIGEN_NO_MALLOC
ei_assert(false && "heap allocation is forbidden (EIGEN_NO_MALLOC is defined)");
#endif
void *result = malloc(size);
void *result = std::malloc(size);
#ifdef EIGEN_EXCEPTIONS
if(!result) throw std::bad_alloc();
#endif
return result;
}
/** \internal construct the elements of an array.
* The \a size parameter tells on how many objects to call the constructor of T.
*/
template<typename T> inline T* ei_construct_elements_of_array(T *ptr, std::size_t size)
{
for (std::size_t i=0; i < size; ++i) ::new (ptr + i) T;
return ptr;
}
/** allocates \a size objects of type T. The returned pointer is guaranteed to have 16 bytes alignment.
* On allocation error, the returned pointer is undefined, but if exceptions are enabled then a std::bad_alloc is thrown.
* The default constructor of T is called.
*/
template<typename T> inline T* ei_aligned_new(size_t size)
template<typename T> inline T* ei_aligned_new(std::size_t size)
{
void *void_result = ei_aligned_malloc(sizeof(T)*size);
return ::new(void_result) T[size];
T *result = reinterpret_cast<T*>(ei_aligned_malloc(sizeof(T)*size));
return ei_construct_elements_of_array(result, size);
}
template<typename T, bool Align> inline T* ei_conditional_aligned_new(size_t size)
template<typename T, bool Align> inline T* ei_conditional_aligned_new(std::size_t size)
{
void *void_result = ei_conditional_aligned_malloc<Align>(sizeof(T)*size);
return ::new(void_result) T[size];
T *result = reinterpret_cast<T*>(ei_conditional_aligned_malloc<Align>(sizeof(T)*size));
return ei_construct_elements_of_array(result, size);
}
/** \internal free memory allocated with ei_aligned_malloc
@@ -160,13 +179,13 @@ template<bool Align> inline void ei_conditional_aligned_free(void *ptr)
template<> inline void ei_conditional_aligned_free<false>(void *ptr)
{
free(ptr);
std::free(ptr);
}
/** \internal delete the elements of an array.
/** \internal destruct the elements of an array.
* The \a size parameters tells on how many objects to call the destructor of T.
*/
template<typename T> inline void ei_delete_elements_of_array(T *ptr, size_t size)
template<typename T> inline void ei_destruct_elements_of_array(T *ptr, std::size_t size)
{
// always destruct an array starting from the end.
while(size) ptr[--size].~T();
@@ -175,33 +194,62 @@ template<typename T> inline void ei_delete_elements_of_array(T *ptr, size_t size
/** \internal delete objects constructed with ei_aligned_new
* The \a size parameters tells on how many objects to call the destructor of T.
*/
template<typename T> inline void ei_aligned_delete(T *ptr, size_t size)
template<typename T> inline void ei_aligned_delete(T *ptr, std::size_t size)
{
ei_delete_elements_of_array<T>(ptr, size);
ei_destruct_elements_of_array<T>(ptr, size);
ei_aligned_free(ptr);
}
/** \internal delete objects constructed with ei_conditional_aligned_new
* The \a size parameters tells on how many objects to call the destructor of T.
*/
template<typename T, bool Align> inline void ei_conditional_aligned_delete(T *ptr, size_t size)
template<typename T, bool Align> inline void ei_conditional_aligned_delete(T *ptr, std::size_t size)
{
ei_delete_elements_of_array<T>(ptr, size);
ei_destruct_elements_of_array<T>(ptr, size);
ei_conditional_aligned_free<Align>(ptr);
}
/** \internal \returns the number of elements which have to be skipped such that data are 16 bytes aligned */
template<typename Scalar>
inline static int ei_alignmentOffset(const Scalar* ptr, int maxOffset)
/** \internal \returns the index of the first element of the array that is well aligned for vectorization.
*
* \param array the address of the start of the array
* \param size the size of the array
*
* \note If no element of the array is well aligned, the size of the array is returned. Typically,
* for example with SSE, "well aligned" means 16-byte-aligned. If vectorization is disabled or if the
* packet size for the given scalar type is 1, then everything is considered well-aligned.
*
* \note If the scalar type is vectorizable, we rely on the following assumptions: sizeof(Scalar) is a
* power of 2, the packet size in bytes is also a power of 2, and is a multiple of sizeof(Scalar). On the
* other hand, we do not assume that the array address is a multiple of sizeof(Scalar), as that fails for
* example with Scalar=double on certain 32-bit platforms, see bug #79.
*
* There is also the variant ei_first_aligned(const MatrixBase&, Integer) defined in Coeffs.h.
*/
template<typename Scalar, typename Integer>
inline static Integer ei_alignmentOffset(const Scalar* array, Integer size)
{
typedef typename ei_packet_traits<Scalar>::type Packet;
const int PacketSize = ei_packet_traits<Scalar>::size;
const int PacketAlignedMask = PacketSize-1;
const bool Vectorized = PacketSize>1;
return Vectorized
? std::min<int>( (PacketSize - (int((size_t(ptr)/sizeof(Scalar))) & PacketAlignedMask))
& PacketAlignedMask, maxOffset)
: 0;
enum { PacketSize = ei_packet_traits<Scalar>::size,
PacketAlignedMask = PacketSize-1
};
if(PacketSize==1)
{
// Either there is no vectorization, or a packet consists of exactly 1 scalar so that all elements
// of the array have the same aligment.
return 0;
}
else if(size_t(array) & (sizeof(Scalar)-1))
{
// There is vectorization for this scalar type, but the array is not aligned to the size of a single scalar.
// Consequently, no element of the array is well aligned.
return size;
}
else
{
return std::min<Integer>( (PacketSize - (Integer((size_t(array)/sizeof(Scalar))) & PacketAlignedMask))
& PacketAlignedMask, size);
}
}
/** \internal
@@ -225,31 +273,31 @@ inline static int ei_alignmentOffset(const Scalar* ptr, int maxOffset)
#define ei_aligned_stack_free(PTR,SIZE) ei_aligned_free(PTR)
#endif
#define ei_aligned_stack_new(TYPE,SIZE) ::new(ei_aligned_stack_alloc(sizeof(TYPE)*SIZE)) TYPE[SIZE]
#define ei_aligned_stack_delete(TYPE,PTR,SIZE) do {ei_delete_elements_of_array<TYPE>(PTR, SIZE); \
#define ei_aligned_stack_new(TYPE,SIZE) ei_construct_elements_of_array(reinterpret_cast<TYPE*>(ei_aligned_stack_alloc(sizeof(TYPE)*SIZE)), SIZE)
#define ei_aligned_stack_delete(TYPE,PTR,SIZE) do {ei_destruct_elements_of_array<TYPE>(PTR, SIZE); \
ei_aligned_stack_free(PTR,sizeof(TYPE)*SIZE);} while(0)
#if EIGEN_ALIGN
#ifdef EIGEN_EXCEPTIONS
#define EIGEN_MAKE_ALIGNED_OPERATOR_NEW_NOTHROW(NeedsToAlign) \
void* operator new(size_t size, const std::nothrow_t&) throw() { \
void* operator new(std::size_t size, const std::nothrow_t&) throw() { \
try { return Eigen::ei_conditional_aligned_malloc<NeedsToAlign>(size); } \
catch (...) { return 0; } \
return 0; \
}
#else
#define EIGEN_MAKE_ALIGNED_OPERATOR_NEW_NOTHROW(NeedsToAlign) \
void* operator new(size_t size, const std::nothrow_t&) throw() { \
void* operator new(std::size_t size, const std::nothrow_t&) throw() { \
return Eigen::ei_conditional_aligned_malloc<NeedsToAlign>(size); \
}
#endif
#define EIGEN_MAKE_ALIGNED_OPERATOR_NEW_IF(NeedsToAlign) \
void *operator new(size_t size) { \
void *operator new(std::size_t size) { \
return Eigen::ei_conditional_aligned_malloc<NeedsToAlign>(size); \
} \
void *operator new[](size_t size) { \
void *operator new[](std::size_t size) { \
return Eigen::ei_conditional_aligned_malloc<NeedsToAlign>(size); \
} \
void operator delete(void * ptr) throw() { Eigen::ei_conditional_aligned_free<NeedsToAlign>(ptr); } \
@@ -257,7 +305,7 @@ inline static int ei_alignmentOffset(const Scalar* ptr, int maxOffset)
/* in-place new and delete. since (at least afaik) there is no actual */ \
/* memory allocated we can safely let the default implementation handle */ \
/* this particular case. */ \
static void *operator new(size_t size, void *ptr) { return ::operator new(size,ptr); } \
static void *operator new(std::size_t size, void *ptr) { return ::operator new(size,ptr); } \
void operator delete(void * memory, void *ptr) throw() { return ::operator delete(memory,ptr); } \
/* nothrow-new (returns zero instead of std::bad_alloc) */ \
EIGEN_MAKE_ALIGNED_OPERATOR_NEW_NOTHROW(NeedsToAlign) \
@@ -291,8 +339,8 @@ template<class T>
class aligned_allocator
{
public:
typedef size_t size_type;
typedef ptrdiff_t difference_type;
typedef std::size_t size_type;
typedef std::ptrdiff_t difference_type;
typedef T* pointer;
typedef const T* const_pointer;
typedef T& reference;

View File

@@ -74,6 +74,7 @@
THIS_METHOD_IS_ONLY_FOR_COLUMN_MAJOR_MATRICES,
THIS_METHOD_IS_ONLY_FOR_ROW_MAJOR_MATRICES,
INVALID_MATRIX_TEMPLATE_PARAMETERS,
INVALID_MATRIXBASE_TEMPLATE_PARAMETERS,
BOTH_MATRICES_MUST_HAVE_THE_SAME_STORAGE_ORDER,
THIS_METHOD_IS_ONLY_FOR_DIAGONAL_MATRIX
};

View File

@@ -37,6 +37,10 @@
//classes inheriting ei_no_assignment_operator don't generate a default operator=.
class ei_no_assignment_operator
{
#if EIGEN_GCC3_OR_OLDER
protected:
void nevermind_this_is_just_to_work_around_a_stupid_gcc3_warning();
#endif
private:
ei_no_assignment_operator& operator=(const ei_no_assignment_operator&);
};
@@ -157,6 +161,19 @@ template<typename T> struct ei_plain_matrix_type_column_major
> type;
};
/* ei_plain_matrix_type_row_major : same as ei_plain_matrix_type but guaranteed to be row-major
*/
template<typename T> struct ei_plain_matrix_type_row_major
{
typedef Matrix<typename ei_traits<T>::Scalar,
ei_traits<T>::RowsAtCompileTime,
ei_traits<T>::ColsAtCompileTime,
AutoAlign | RowMajor,
ei_traits<T>::MaxRowsAtCompileTime,
ei_traits<T>::MaxColsAtCompileTime
> type;
};
template<typename T> struct ei_must_nest_by_value { enum { ret = false }; };
template<typename T> struct ei_must_nest_by_value<NestByValue<T> > { enum { ret = true }; };

View File

@@ -52,9 +52,9 @@ public:
typedef _Scalar Scalar;
typedef typename NumTraits<Scalar>::Real RealScalar;
typedef Matrix<Scalar,AmbientDimAtCompileTime,1> VectorType;
typedef Matrix<Scalar,AmbientDimAtCompileTime==Dynamic
typedef Matrix<Scalar,int(AmbientDimAtCompileTime)==Dynamic
? Dynamic
: AmbientDimAtCompileTime+1,1> Coefficients;
: int(AmbientDimAtCompileTime)+1,1> Coefficients;
typedef Block<Coefficients,AmbientDimAtCompileTime,1> NormalReturnType;
/** Default constructor without initialization */

View File

@@ -450,22 +450,31 @@ inline Scalar Quaternion<Scalar>::angularDistance(const Quaternion& other) const
template <typename Scalar>
Quaternion<Scalar> Quaternion<Scalar>::slerp(Scalar t, const Quaternion& other) const
{
static const Scalar one = Scalar(1) - precision<Scalar>();
static const Scalar one = Scalar(1) - machine_epsilon<Scalar>();
Scalar d = this->dot(other);
Scalar absD = ei_abs(d);
Scalar scale0;
Scalar scale1;
if (absD>=one)
return *this;
{
scale0 = Scalar(1) - t;
scale1 = t;
}
else
{
// theta is the angle between the 2 quaternions
Scalar theta = std::acos(absD);
Scalar sinTheta = ei_sin(theta);
// theta is the angle between the 2 quaternions
Scalar theta = std::acos(absD);
Scalar sinTheta = ei_sin(theta);
scale0 = ei_sin( ( Scalar(1) - t ) * theta) / sinTheta;
scale1 = ei_sin( ( t * theta) ) / sinTheta;
if (d<0)
scale1 = -scale1;
}
Scalar scale0 = ei_sin( ( Scalar(1) - t ) * theta) / sinTheta;
Scalar scale1 = ei_sin( ( t * theta) ) / sinTheta;
if (d<0)
scale1 = -scale1;
return Quaternion(scale0 * m_coeffs + scale1 * other.m_coeffs);
return Quaternion<Scalar>(scale0 * coeffs() + scale1 * other.coeffs());
}
// set from a rotation matrix

View File

@@ -29,8 +29,8 @@
*** Part 1 : optimized implementations for fixed-size 2,3,4 cases ***
********************************************************************/
template<typename MatrixType>
void ei_compute_inverse_in_size2_case(const MatrixType& matrix, MatrixType* result)
template<typename XprType, typename MatrixType>
void ei_compute_inverse_in_size2_case(const XprType& matrix, MatrixType* result)
{
typedef typename MatrixType::Scalar Scalar;
const Scalar invdet = Scalar(1) / matrix.determinant();
@@ -54,10 +54,10 @@ bool ei_compute_inverse_in_size2_case_with_check(const XprType& matrix, MatrixTy
return true;
}
template<typename MatrixType>
void ei_compute_inverse_in_size3_case(const MatrixType& matrix, MatrixType* result)
template<typename Derived, typename OtherDerived>
void ei_compute_inverse_in_size3_case(const Derived& matrix, OtherDerived* result)
{
typedef typename MatrixType::Scalar Scalar;
typedef typename Derived::Scalar Scalar;
const Scalar det_minor00 = matrix.minor(0,0).determinant();
const Scalar det_minor10 = matrix.minor(1,0).determinant();
const Scalar det_minor20 = matrix.minor(2,0).determinant();
@@ -75,140 +75,204 @@ void ei_compute_inverse_in_size3_case(const MatrixType& matrix, MatrixType* resu
result->coeffRef(2, 2) = matrix.minor(2,2).determinant() * invdet;
}
template<typename MatrixType>
bool ei_compute_inverse_in_size4_case_helper(const MatrixType& matrix, MatrixType* result)
template<typename Derived, typename OtherDerived, typename Scalar = typename Derived::Scalar>
struct ei_compute_inverse_in_size4_case
{
/* Let's split M into four 2x2 blocks:
* (P Q)
* (R S)
* If P is invertible, with inverse denoted by P_inverse, and if
* (S - R*P_inverse*Q) is also invertible, then the inverse of M is
* (P' Q')
* (R' S')
* where
* S' = (S - R*P_inverse*Q)^(-1)
* P' = P1 + (P1*Q) * S' *(R*P_inverse)
* Q' = -(P_inverse*Q) * S'
* R' = -S' * (R*P_inverse)
*/
typedef Block<MatrixType,2,2> XprBlock22;
typedef typename MatrixBase<XprBlock22>::PlainMatrixType Block22;
Block22 P_inverse;
if(ei_compute_inverse_in_size2_case_with_check(matrix.template block<2,2>(0,0), &P_inverse))
static void run(const Derived& matrix, OtherDerived& result)
{
const Block22 Q = matrix.template block<2,2>(0,2);
const Block22 P_inverse_times_Q = P_inverse * Q;
const XprBlock22 R = matrix.template block<2,2>(2,0);
const Block22 R_times_P_inverse = R * P_inverse;
const Block22 R_times_P_inverse_times_Q = R_times_P_inverse * Q;
const XprBlock22 S = matrix.template block<2,2>(2,2);
const Block22 X = S - R_times_P_inverse_times_Q;
Block22 Y;
ei_compute_inverse_in_size2_case(X, &Y);
result->template block<2,2>(2,2) = Y;
result->template block<2,2>(2,0) = - Y * R_times_P_inverse;
const Block22 Z = P_inverse_times_Q * Y;
result->template block<2,2>(0,2) = - Z;
result->template block<2,2>(0,0) = P_inverse + Z * R_times_P_inverse;
return true;
result.coeffRef(0,0) = matrix.minor(0,0).determinant();
result.coeffRef(1,0) = -matrix.minor(0,1).determinant();
result.coeffRef(2,0) = matrix.minor(0,2).determinant();
result.coeffRef(3,0) = -matrix.minor(0,3).determinant();
result.coeffRef(0,2) = matrix.minor(2,0).determinant();
result.coeffRef(1,2) = -matrix.minor(2,1).determinant();
result.coeffRef(2,2) = matrix.minor(2,2).determinant();
result.coeffRef(3,2) = -matrix.minor(2,3).determinant();
result.coeffRef(0,1) = -matrix.minor(1,0).determinant();
result.coeffRef(1,1) = matrix.minor(1,1).determinant();
result.coeffRef(2,1) = -matrix.minor(1,2).determinant();
result.coeffRef(3,1) = matrix.minor(1,3).determinant();
result.coeffRef(0,3) = -matrix.minor(3,0).determinant();
result.coeffRef(1,3) = matrix.minor(3,1).determinant();
result.coeffRef(2,3) = -matrix.minor(3,2).determinant();
result.coeffRef(3,3) = matrix.minor(3,3).determinant();
result /= (matrix.col(0).cwise()*result.row(0).transpose()).sum();
}
else
{
return false;
}
}
};
template<typename MatrixType>
void ei_compute_inverse_in_size4_case(const MatrixType& matrix, MatrixType* result)
#ifdef EIGEN_VECTORIZE_SSE
// The SSE code for the 4x4 float matrix inverse in this file comes from the file
// ftp://download.intel.com/design/PentiumIII/sml/24504301.pdf
// its copyright information is:
// Copyright (C) 1999 Intel Corporation
// See page ii of that document for legal stuff. Not being lawyers, we just assume
// here that if Intel makes this document publically available, with source code
// and detailed explanations, it's because they want their CPUs to be fed with
// good code, and therefore they presumably don't mind us using it in Eigen.
template<typename Derived, typename OtherDerived>
struct ei_compute_inverse_in_size4_case<Derived, OtherDerived, float>
{
if(ei_compute_inverse_in_size4_case_helper(matrix, result))
static void run(const Derived& matrix, OtherDerived& result)
{
// good ! The topleft 2x2 block was invertible, so the 2x2 blocks approach is successful.
return;
// Variables (Streaming SIMD Extensions registers) which will contain cofactors and, later, the
// lines of the inverted matrix.
__m128 minor0, minor1, minor2, minor3;
// Variables which will contain the lines of the reference matrix and, later (after the transposition),
// the columns of the original matrix.
__m128 row0, row1, row2, row3;
// Temporary variables and the variable that will contain the matrix determinant.
__m128 det, tmp1;
// Matrix transposition
const float *src = matrix.data();
tmp1 = _mm_loadh_pi(_mm_castpd_ps(_mm_load_sd((double*)src)), (__m64*)(src+ 4));
row1 = _mm_loadh_pi(_mm_castpd_ps(_mm_load_sd((double*)(src+8))), (__m64*)(src+12));
row0 = _mm_shuffle_ps(tmp1, row1, 0x88);
row1 = _mm_shuffle_ps(row1, tmp1, 0xDD);
tmp1 = _mm_loadh_pi(_mm_castpd_ps(_mm_load_sd((double*)(src+ 2))), (__m64*)(src+ 6));
row3 = _mm_loadh_pi(_mm_castpd_ps(_mm_load_sd((double*)(src+10))), (__m64*)(src+14));
row2 = _mm_shuffle_ps(tmp1, row3, 0x88);
row3 = _mm_shuffle_ps(row3, tmp1, 0xDD);
// Cofactors calculation. Because in the process of cofactor computation some pairs in three-
// element products are repeated, it is not reasonable to load these pairs anew every time. The
// values in the registers with these pairs are formed using shuffle instruction. Cofactors are
// calculated row by row (4 elements are placed in 1 SP FP SIMD floating point register).
tmp1 = _mm_mul_ps(row2, row3);
tmp1 = _mm_shuffle_ps(tmp1, tmp1, 0xB1);
minor0 = _mm_mul_ps(row1, tmp1);
minor1 = _mm_mul_ps(row0, tmp1);
tmp1 = _mm_shuffle_ps(tmp1, tmp1, 0x4E);
minor0 = _mm_sub_ps(_mm_mul_ps(row1, tmp1), minor0);
minor1 = _mm_sub_ps(_mm_mul_ps(row0, tmp1), minor1);
minor1 = _mm_shuffle_ps(minor1, minor1, 0x4E);
// -----------------------------------------------
tmp1 = _mm_mul_ps(row1, row2);
tmp1 = _mm_shuffle_ps(tmp1, tmp1, 0xB1);
minor0 = _mm_add_ps(_mm_mul_ps(row3, tmp1), minor0);
minor3 = _mm_mul_ps(row0, tmp1);
tmp1 = _mm_shuffle_ps(tmp1, tmp1, 0x4E);
minor0 = _mm_sub_ps(minor0, _mm_mul_ps(row3, tmp1));
minor3 = _mm_sub_ps(_mm_mul_ps(row0, tmp1), minor3);
minor3 = _mm_shuffle_ps(minor3, minor3, 0x4E);
// -----------------------------------------------
tmp1 = _mm_mul_ps(_mm_shuffle_ps(row1, row1, 0x4E), row3);
tmp1 = _mm_shuffle_ps(tmp1, tmp1, 0xB1);
row2 = _mm_shuffle_ps(row2, row2, 0x4E);
minor0 = _mm_add_ps(_mm_mul_ps(row2, tmp1), minor0);
minor2 = _mm_mul_ps(row0, tmp1);
tmp1 = _mm_shuffle_ps(tmp1, tmp1, 0x4E);
minor0 = _mm_sub_ps(minor0, _mm_mul_ps(row2, tmp1));
minor2 = _mm_sub_ps(_mm_mul_ps(row0, tmp1), minor2);
minor2 = _mm_shuffle_ps(minor2, minor2, 0x4E);
// -----------------------------------------------
tmp1 = _mm_mul_ps(row0, row1);
tmp1 = _mm_shuffle_ps(tmp1, tmp1, 0xB1);
minor2 = _mm_add_ps(_mm_mul_ps(row3, tmp1), minor2);
minor3 = _mm_sub_ps(_mm_mul_ps(row2, tmp1), minor3);
tmp1 = _mm_shuffle_ps(tmp1, tmp1, 0x4E);
minor2 = _mm_sub_ps(_mm_mul_ps(row3, tmp1), minor2);
minor3 = _mm_sub_ps(minor3, _mm_mul_ps(row2, tmp1));
// -----------------------------------------------
tmp1 = _mm_mul_ps(row0, row3);
tmp1 = _mm_shuffle_ps(tmp1, tmp1, 0xB1);
minor1 = _mm_sub_ps(minor1, _mm_mul_ps(row2, tmp1));
minor2 = _mm_add_ps(_mm_mul_ps(row1, tmp1), minor2);
tmp1 = _mm_shuffle_ps(tmp1, tmp1, 0x4E);
minor1 = _mm_add_ps(_mm_mul_ps(row2, tmp1), minor1);
minor2 = _mm_sub_ps(minor2, _mm_mul_ps(row1, tmp1));
// -----------------------------------------------
tmp1 = _mm_mul_ps(row0, row2);
tmp1 = _mm_shuffle_ps(tmp1, tmp1, 0xB1);
minor1 = _mm_add_ps(_mm_mul_ps(row3, tmp1), minor1);
minor3 = _mm_sub_ps(minor3, _mm_mul_ps(row1, tmp1));
tmp1 = _mm_shuffle_ps(tmp1, tmp1, 0x4E);
minor1 = _mm_sub_ps(minor1, _mm_mul_ps(row3, tmp1));
minor3 = _mm_add_ps(_mm_mul_ps(row1, tmp1), minor3);
// Evaluation of determinant and its reciprocal value. In the original Intel document,
// 1/det was evaluated using a fast rcpps command with subsequent approximation using
// the Newton-Raphson algorithm. Here, we go for a IEEE-compliant division instead,
// so as to not compromise precision at all.
det = _mm_mul_ps(row0, minor0);
det = _mm_add_ps(_mm_shuffle_ps(det, det, 0x4E), det);
det = _mm_add_ss(_mm_shuffle_ps(det, det, 0xB1), det);
// tmp1= _mm_rcp_ss(det);
// det= _mm_sub_ss(_mm_add_ss(tmp1, tmp1), _mm_mul_ss(det, _mm_mul_ss(tmp1, tmp1)));
det = _mm_div_ss(_mm_set_ss(1.0f), det); // <--- yay, one original line not copied from Intel
det = _mm_shuffle_ps(det, det, 0x00);
// warning, Intel's variable naming is very confusing: now 'det' is 1/det !
// Multiplication of cofactors by 1/det. Storing the inverse matrix to the address in pointer src.
minor0 = _mm_mul_ps(det, minor0);
float *dst = result.data();
_mm_storel_pi((__m64*)(dst), minor0);
_mm_storeh_pi((__m64*)(dst+2), minor0);
minor1 = _mm_mul_ps(det, minor1);
_mm_storel_pi((__m64*)(dst+4), minor1);
_mm_storeh_pi((__m64*)(dst+6), minor1);
minor2 = _mm_mul_ps(det, minor2);
_mm_storel_pi((__m64*)(dst+ 8), minor2);
_mm_storeh_pi((__m64*)(dst+10), minor2);
minor3 = _mm_mul_ps(det, minor3);
_mm_storel_pi((__m64*)(dst+12), minor3);
_mm_storeh_pi((__m64*)(dst+14), minor3);
}
else
{
// rare case: the topleft 2x2 block is not invertible (but the matrix itself is assumed to be).
// since this is a rare case, we don't need to optimize it. We just want to handle it with little
// additional code.
MatrixType m(matrix);
m.row(0).swap(m.row(2));
m.row(1).swap(m.row(3));
if(ei_compute_inverse_in_size4_case_helper(m, result))
{
// good, the topleft 2x2 block of m is invertible. Since m is different from matrix in that some
// rows were permuted, the actual inverse of matrix is derived from the inverse of m by permuting
// the corresponding columns.
result->col(0).swap(result->col(2));
result->col(1).swap(result->col(3));
}
else
{
// last possible case. Since matrix is assumed to be invertible, this last case has to work.
// first, undo the swaps previously made
m.row(0).swap(m.row(2));
m.row(1).swap(m.row(3));
// swap row 0 with the the row among 0 and 1 that has the biggest 2 first coeffs
int swap0with = ei_abs(m.coeff(0,0))+ei_abs(m.coeff(0,1))>ei_abs(m.coeff(1,0))+ei_abs(m.coeff(1,1)) ? 0 : 1;
m.row(0).swap(m.row(swap0with));
// swap row 1 with the the row among 2 and 3 that has the biggest 2 first coeffs
int swap1with = ei_abs(m.coeff(2,0))+ei_abs(m.coeff(2,1))>ei_abs(m.coeff(3,0))+ei_abs(m.coeff(3,1)) ? 2 : 3;
m.row(1).swap(m.row(swap1with));
ei_compute_inverse_in_size4_case_helper(m, result);
result->col(1).swap(result->col(swap1with));
result->col(0).swap(result->col(swap0with));
}
}
}
};
#endif
/***********************************************
*** Part 2 : selector and MatrixBase methods ***
***********************************************/
template<typename MatrixType, int Size = MatrixType::RowsAtCompileTime>
template<typename Derived, typename OtherDerived, int Size = Derived::RowsAtCompileTime>
struct ei_compute_inverse
{
static inline void run(const MatrixType& matrix, MatrixType* result)
static inline void run(const Derived& matrix, OtherDerived* result)
{
LU<MatrixType> lu(matrix);
LU<Derived> lu(matrix);
lu.computeInverse(result);
}
};
template<typename MatrixType>
struct ei_compute_inverse<MatrixType, 1>
template<typename Derived, typename OtherDerived>
struct ei_compute_inverse<Derived, OtherDerived, 1>
{
static inline void run(const MatrixType& matrix, MatrixType* result)
static inline void run(const Derived& matrix, OtherDerived* result)
{
typedef typename MatrixType::Scalar Scalar;
typedef typename Derived::Scalar Scalar;
result->coeffRef(0,0) = Scalar(1) / matrix.coeff(0,0);
}
};
template<typename MatrixType>
struct ei_compute_inverse<MatrixType, 2>
template<typename Derived, typename OtherDerived>
struct ei_compute_inverse<Derived, OtherDerived, 2>
{
static inline void run(const MatrixType& matrix, MatrixType* result)
static inline void run(const Derived& matrix, OtherDerived* result)
{
ei_compute_inverse_in_size2_case(matrix, result);
}
};
template<typename MatrixType>
struct ei_compute_inverse<MatrixType, 3>
template<typename Derived, typename OtherDerived>
struct ei_compute_inverse<Derived, OtherDerived, 3>
{
static inline void run(const MatrixType& matrix, MatrixType* result)
static inline void run(const Derived& matrix, OtherDerived* result)
{
ei_compute_inverse_in_size3_case(matrix, result);
}
};
template<typename MatrixType>
struct ei_compute_inverse<MatrixType, 4>
template<typename Derived, typename OtherDerived>
struct ei_compute_inverse<Derived, OtherDerived, 4>
{
static inline void run(const MatrixType& matrix, MatrixType* result)
static inline void run(const Derived& matrix, OtherDerived* result)
{
ei_compute_inverse_in_size4_case(matrix, result);
ei_compute_inverse_in_size4_case<Derived, OtherDerived>::run(matrix, *result);
}
};
@@ -226,11 +290,12 @@ struct ei_compute_inverse<MatrixType, 4>
* \sa inverse()
*/
template<typename Derived>
inline void MatrixBase<Derived>::computeInverse(PlainMatrixType *result) const
template<typename OtherDerived>
inline void MatrixBase<Derived>::computeInverse(MatrixBase<OtherDerived> *result) const
{
ei_assert(rows() == cols());
EIGEN_STATIC_ASSERT(NumTraits<Scalar>::HasFloatingPoint,NUMERIC_TYPE_MUST_BE_FLOATING_POINT)
ei_compute_inverse<PlainMatrixType>::run(eval(), result);
ei_compute_inverse<PlainMatrixType, OtherDerived>::run(eval(), static_cast<OtherDerived*>(result));
}
/** \lu_module

View File

@@ -68,7 +68,7 @@ template<typename MatrixType> class LU
typedef Matrix<Scalar, 1, MatrixType::ColsAtCompileTime> RowVectorType;
typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> ColVectorType;
enum { MaxSmallDimAtCompileTime = EIGEN_ENUM_MIN(
enum { MaxSmallDimAtCompileTime = EIGEN_SIZE_MIN(
MatrixType::MaxColsAtCompileTime,
MatrixType::MaxRowsAtCompileTime)
};
@@ -297,7 +297,8 @@ template<typename MatrixType> class LU
*
* \sa MatrixBase::computeInverse(), inverse()
*/
inline void computeInverse(MatrixType *result) const
template<typename ResultType>
inline void computeInverse(ResultType *result) const
{
solve(MatrixType::Identity(m_lu.rows(), m_lu.cols()), result);
}
@@ -508,7 +509,7 @@ bool LU<MatrixType>::solve(
if(!isSurjective())
{
// is c is in the image of U ?
RealScalar biggest_in_c = m_rank>0 ? c.corner(TopLeft, m_rank, c.cols()).cwise().abs().maxCoeff() : 0;
RealScalar biggest_in_c = m_rank>0 ? c.corner(TopLeft, m_rank, c.cols()).cwise().abs().maxCoeff() : RealScalar(0);
for(int col = 0; col < c.cols(); ++col)
for(int row = m_rank; row < c.rows(); ++row)
if(!ei_isMuchSmallerThan(c.coeff(row,col), biggest_in_c, m_precision))

View File

@@ -293,7 +293,7 @@ void Tridiagonalization<MatrixType>::_compute(MatrixType& matA, CoeffVectorType&
{
int starti = i+1;
int alignedEnd = starti;
if (PacketSize>1)
if (PacketSize>1 && (int(MatrixType::Flags)&RowMajor) == 0)
{
int alignedStart = (starti) + ei_alignmentOffset(&matA.coeffRef(starti,j1), n-starti);
alignedEnd = alignedStart + ((n-alignedStart)/PacketSize)*PacketSize;

View File

@@ -49,7 +49,7 @@ template<typename MatrixType> class SVD
enum {
PacketSize = ei_packet_traits<Scalar>::size,
AlignmentMask = int(PacketSize)-1,
MinSize = EIGEN_ENUM_MIN(MatrixType::RowsAtCompileTime, MatrixType::ColsAtCompileTime)
MinSize = EIGEN_SIZE_MIN(MatrixType::RowsAtCompileTime, MatrixType::ColsAtCompileTime)
};
typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, 1> ColVector;
@@ -61,6 +61,8 @@ template<typename MatrixType> class SVD
public:
SVD() {} // a user who relied on compiler-generated default compiler reported problems with MSVC in 2.0.7
SVD(const MatrixType& matrix)
: m_matU(matrix.rows(), std::min(matrix.rows(), matrix.cols())),
m_matV(matrix.cols(),matrix.cols()),
@@ -107,6 +109,8 @@ void SVD<MatrixType>::compute(const MatrixType& matrix)
const int m = matrix.rows();
const int n = matrix.cols();
const int nu = std::min(m,n);
ei_assert(m>=n && "In Eigen 2.0, SVD only works for MxN matrices with M>=N. Sorry!");
ei_assert(m>1 && "In Eigen 2.0, SVD doesn't work on 1x1 matrices");
m_matU.resize(m, nu);
m_matU.setZero();

View File

@@ -99,6 +99,8 @@ template<typename _Scalar> class AmbiVector
allocSize = allocSize/sizeof(Scalar) + (allocSize%sizeof(Scalar)>0?1:0);
Scalar* newBuffer = new Scalar[allocSize];
memcpy(newBuffer, m_buffer, copyElements * sizeof(ListEl));
delete[] m_buffer;
m_buffer = newBuffer;
}
protected:
@@ -238,8 +240,11 @@ Scalar& AmbiVector<Scalar>::coeffRef(int i)
else
{
if (m_llSize>=m_allocatedElements)
{
reallocateSparse();
ei_internal_assert(m_llSize<m_size && "internal error: overflow in sparse mode");
llElements = reinterpret_cast<ListEl*>(m_buffer);
}
ei_internal_assert(m_llSize<m_allocatedElements && "internal error: overflow in sparse mode");
// let's insert a new coefficient
ListEl& el = llElements[m_llSize];
el.value = Scalar(0);
@@ -365,6 +370,9 @@ class AmbiVector<_Scalar>::Iterator
int m_cachedIndex; // current coordinate
Scalar m_cachedValue; // current value
bool m_isDense; // mode of the vector
private:
Iterator& operator=(const Iterator&);
};

View File

@@ -289,9 +289,11 @@ class DynamicSparseMatrix<Scalar,_Flags>::InnerIterator : public SparseVector<Sc
inline int row() const { return IsRowMajor ? m_outer : Base::index(); }
inline int col() const { return IsRowMajor ? Base::index() : m_outer; }
protected:
const int m_outer;
private:
InnerIterator& operator=(const InnerIterator&);
};
#endif // EIGEN_DYNAMIC_SPARSEMATRIX_H

View File

@@ -53,6 +53,9 @@ class SparseInnerVectorSet : ei_no_assignment_operator,
inline InnerIterator(const SparseInnerVectorSet& xpr, int outer)
: MatrixType::InnerIterator(xpr.m_matrix, xpr.m_outerStart + outer)
{}
private:
InnerIterator& operator=(const InnerIterator&);
};
inline SparseInnerVectorSet(const MatrixType& matrix, int outerStart, int outerSize)
@@ -110,6 +113,8 @@ class SparseInnerVectorSet<DynamicSparseMatrix<_Scalar, _Options>, Size>
inline InnerIterator(const SparseInnerVectorSet& xpr, int outer)
: MatrixType::InnerIterator(xpr.m_matrix, xpr.m_outerStart + outer)
{}
private:
InnerIterator& operator=(const InnerIterator&);
};
inline SparseInnerVectorSet(const MatrixType& matrix, int outerStart, int outerSize)

View File

@@ -156,6 +156,9 @@ template<typename ExpressionType> class SparseCwise
protected:
ExpressionTypeNested m_matrix;
private:
SparseCwise& operator=(const SparseCwise&);
};
template<typename Derived>

View File

@@ -126,6 +126,8 @@ class SparseCwiseBinaryOp<BinaryOp,Lhs,Rhs>::InnerIterator
EIGEN_STRONG_INLINE InnerIterator(const SparseCwiseBinaryOp& binOp, int outer)
: Base(binOp,outer)
{}
private:
InnerIterator& operator=(const InnerIterator&);
};
/***************************************************************************
@@ -197,6 +199,9 @@ class ei_sparse_cwise_binary_op_inner_iterator_selector<BinaryOp, Lhs, Rhs, Deri
const BinaryOp& m_functor;
Scalar m_value;
int m_id;
private:
ei_sparse_cwise_binary_op_inner_iterator_selector& operator=(const ei_sparse_cwise_binary_op_inner_iterator_selector&);
};
// sparse - sparse (product)
@@ -250,6 +255,9 @@ class ei_sparse_cwise_binary_op_inner_iterator_selector<ei_scalar_product_op<T>,
LhsIterator m_lhsIter;
RhsIterator m_rhsIter;
const BinaryFunc& m_functor;
private:
ei_sparse_cwise_binary_op_inner_iterator_selector& operator=(const ei_sparse_cwise_binary_op_inner_iterator_selector&);
};
// sparse - dense (product)
@@ -290,6 +298,9 @@ class ei_sparse_cwise_binary_op_inner_iterator_selector<ei_scalar_product_op<T>,
LhsIterator m_lhsIter;
const BinaryFunc m_functor;
const int m_outer;
private:
ei_sparse_cwise_binary_op_inner_iterator_selector& operator=(const ei_sparse_cwise_binary_op_inner_iterator_selector&);
};
// sparse - dense (product)

View File

@@ -90,6 +90,9 @@ class SparseCwiseUnaryOp<UnaryOp,MatrixType>::InnerIterator
protected:
MatrixTypeIterator m_iter;
const UnaryOp m_functor;
private:
InnerIterator& operator=(const InnerIterator&);
};
template<typename Derived>

View File

@@ -120,6 +120,8 @@ class ei_sparse_diagonal_product_inner_iterator_selector
const SparseDiagonalProductType& expr, int outer)
: Base(expr.rhs().innerVector(outer) .cwise()* expr.lhs().diagonal(), 0)
{}
private:
ei_sparse_diagonal_product_inner_iterator_selector& operator=(const ei_sparse_diagonal_product_inner_iterator_selector&);
};
template<typename Lhs, typename Rhs, typename SparseDiagonalProductType>

View File

@@ -64,16 +64,21 @@ template<typename ExpressionType, unsigned int Added, unsigned int Removed> clas
protected:
ExpressionTypeNested m_matrix;
private:
SparseFlagged& operator=(const SparseFlagged&);
};
template<typename ExpressionType, unsigned int Added, unsigned int Removed>
class SparseFlagged<ExpressionType,Added,Removed>::InnerIterator : public ExpressionType::InnerIterator
{
public:
EIGEN_STRONG_INLINE InnerIterator(const SparseFlagged& xpr, int outer)
: ExpressionType::InnerIterator(xpr.m_matrix, outer)
{}
private:
InnerIterator& operator=(const InnerIterator&);
};
template<typename ExpressionType, unsigned int Added, unsigned int Removed>

View File

@@ -259,19 +259,21 @@ class SparseMatrix
m_data.resize(k,0);
}
/** Resizes the matrix to a \a rows x \a cols matrix and initializes it to zero
* \sa resizeNonZeros(int), reserve(), setZero()
*/
void resize(int rows, int cols)
{
// std::cerr << this << " resize " << rows << "x" << cols << "\n";
const int outerSize = IsRowMajor ? rows : cols;
m_innerSize = IsRowMajor ? cols : rows;
m_data.clear();
if (m_outerSize != outerSize)
if (m_outerSize != outerSize || m_outerSize==0)
{
delete[] m_outerIndex;
m_outerIndex = new int [outerSize+1];
m_outerSize = outerSize;
memset(m_outerIndex, 0, (m_outerSize+1)*sizeof(int));
}
memset(m_outerIndex, 0, (m_outerSize+1)*sizeof(int));
}
void resizeNonZeros(int size)
{
@@ -442,6 +444,9 @@ class SparseMatrix<Scalar,_Flags>::InnerIterator
int m_id;
const int m_start;
const int m_end;
private:
InnerIterator& operator=(const InnerIterator&);
};
#endif // EIGEN_SPARSEMATRIX_H

View File

@@ -97,7 +97,7 @@ struct ei_traits<SparseProduct<LhsNested, RhsNested, ProductMode> >
RowsAtCompileTime = _LhsNested::RowsAtCompileTime,
ColsAtCompileTime = _RhsNested::ColsAtCompileTime,
InnerSize = EIGEN_ENUM_MIN(_LhsNested::ColsAtCompileTime, _RhsNested::RowsAtCompileTime),
InnerSize = EIGEN_SIZE_MIN(_LhsNested::ColsAtCompileTime, _RhsNested::RowsAtCompileTime),
MaxRowsAtCompileTime = _LhsNested::MaxRowsAtCompileTime,
MaxColsAtCompileTime = _RhsNested::MaxColsAtCompileTime,

View File

@@ -62,15 +62,20 @@ template<typename MatrixType> class SparseTranspose
protected:
const typename MatrixType::Nested m_matrix;
private:
SparseTranspose& operator=(const SparseTranspose&);
};
template<typename MatrixType> class SparseTranspose<MatrixType>::InnerIterator : public MatrixType::InnerIterator
{
public:
EIGEN_STRONG_INLINE InnerIterator(const SparseTranspose& trans, int outer)
: MatrixType::InnerIterator(trans.m_matrix, outer)
{}
private:
InnerIterator& operator=(const InnerIterator&);
};
template<typename MatrixType> class SparseTranspose<MatrixType>::ReverseInnerIterator : public MatrixType::ReverseInnerIterator

View File

@@ -360,6 +360,9 @@ class SparseVector<Scalar,_Flags>::InnerIterator
const CompressedStorage<Scalar>& m_data;
int m_id;
const int m_end;
private:
InnerIterator& operator=(const InnerIterator&);
};
#endif // EIGEN_SPARSEVECTOR_H

View File

@@ -43,8 +43,11 @@ struct ei_solve_triangular_selector<Lhs,Rhs,LowerTriangular,RowMajor|IsSparse>
{
lastVal = it.value();
lastIndex = it.index();
if(lastIndex == i)
break;
tmp -= lastVal * other.coeff(lastIndex,col);
}
if (Lhs::Flags & UnitDiagBit)
other.coeffRef(i,col) = tmp;
else

View File

@@ -1,8 +1,8 @@
// This file is part of Eigen, a lightweight C++ template library
// for linear algebra. Eigen itself is part of the KDE project.
// for linear algebra.
//
// Copyright (C) 2008 Gael Guennebaud <g.gael@free.fr>
// Copyright (C) 2006-2008 Benoit Jacob <jacob.benoit.1@gmail.com>
// Copyright (C) 2009 Benoit Jacob <jacob.benoit.1@gmail.com>
//
// Eigen is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
@@ -26,8 +26,15 @@
#ifndef EIGEN_BENCH_TIMER_H
#define EIGEN_BENCH_TIMER_H
#include <sys/time.h>
#if defined(_WIN32) || defined(__CYGWIN__)
#define NOMINMAX
#define WIN32_LEAN_AND_MEAN
#include <windows.h>
#else
#include <time.h>
#include <unistd.h>
#endif
#include <cstdlib>
#include <numeric>
@@ -35,12 +42,25 @@ namespace Eigen
{
/** Elapsed time timer keeping the best try.
*
* On POSIX platforms we use clock_gettime with CLOCK_PROCESS_CPUTIME_ID.
* On Windows we use QueryPerformanceCounter
*
* Important: on linux, you must link with -lrt
*/
class BenchTimer
{
public:
BenchTimer() { reset(); }
BenchTimer()
{
#if defined(_WIN32) || defined(__CYGWIN__)
LARGE_INTEGER freq;
QueryPerformanceFrequency(&freq);
m_frequency = (double)freq.QuadPart;
#endif
reset();
}
~BenchTimer() {}
@@ -51,23 +71,34 @@ public:
m_best = std::min(m_best, getTime() - m_start);
}
/** Return the best elapsed time.
/** Return the best elapsed time in seconds.
*/
inline double value(void)
{
return m_best;
return m_best;
}
#if defined(_WIN32) || defined(__CYGWIN__)
inline double getTime(void)
#else
static inline double getTime(void)
#endif
{
struct timeval tv;
struct timezone tz;
gettimeofday(&tv, &tz);
return (double)tv.tv_sec + 1.e-6 * (double)tv.tv_usec;
#ifdef WIN32
LARGE_INTEGER query_ticks;
QueryPerformanceCounter(&query_ticks);
return query_ticks.QuadPart/m_frequency;
#else
timespec ts;
clock_gettime(CLOCK_PROCESS_CPUTIME_ID, &ts);
return double(ts.tv_sec) + 1e-9 * double(ts.tv_nsec);
#endif
}
protected:
#if defined(_WIN32) || defined(__CYGWIN__)
double m_frequency;
#endif
double m_best, m_start;
};

View File

@@ -57,10 +57,10 @@ void makeFloor(const MatrixBase<OtherDerived>& other) { derived() = derived().cw
template<typename OtherDerived>
void makeCeil(const MatrixBase<OtherDerived>& other) { derived() = derived().cwise().max(other.derived()); }
const typename Cwise<Derived>::ScalarAddReturnType
operator+(const Scalar& scalar) const { return cwise() + scalar }
const const CwiseUnaryOp<ei_scalar_add_op<Scalar>, Derived>
operator+(const Scalar& scalar) const { return cwise() + scalar; }
friend const typename Cwise<Derived>::ScalarAddReturnType
friend const CwiseUnaryOp<ei_scalar_add_op<Scalar>, Derived>
operator+(const Scalar& scalar, const MatrixBase<Derived>& mat) { return mat + scalar; }
\endcode
@@ -81,7 +81,7 @@ In order to add support for a custom type \c T you need:
3 - define a couple of math functions for your type such as: ei_sqrt, ei_abs, etc...
(see the file Eigen/src/Core/MathFunctions.h)
Here is a concrete example adding support for the Adolc's \c adouble type. <a href="http://www.math.tu-dresden.de/~adol-c/">Adolc</a> is an automatic differentiation library. The type \c adouble is basically a real value tracking the values of any number of partial derivatives.
Here is a concrete example adding support for the Adolc's \c adouble type. <a href="https://projects.coin-or.org/ADOL-C">Adolc</a> is an automatic differentiation library. The type \c adouble is basically a real value tracking the values of any number of partial derivatives.
\code
#ifndef ADLOCSUPPORT_H

View File

@@ -229,17 +229,24 @@ Of course, fixed-size matrices can't be resized.
\subsection TutorialMap Map
Any memory buffer can be mapped as an Eigen expression:
<table class="tutorial_code"><tr><td>
Any memory buffer can be mapped as an Eigen expression using the Map() static method:
\code
std::vector<float> stlarray(10);
Map<VectorXf>(&stlarray[0], stlarray.size()).setOnes();
int data[4] = 1, 2, 3, 4;
Matrix2i mat2x2(data);
MatrixXi mat2x2 = Map<Matrix2i>(data);
MatrixXi mat2x2 = Map<MatrixXi>(data,2,2);
VectorXf::Map(&stlarray[0], stlarray.size()).squaredNorm();
\endcode
Here VectorXf::Map returns an object of class Map<VectorXf>, which behaves like a VectorXf except that it uses the existing array. You can write to this object, that will write to the existing array. You can also construct a named obtect to reuse it:
\code
float array[rows*cols];
Map<MatrixXf> m(array,rows,cols);
m = othermatrix1 * othermatrix2;
m.eigenvalues();
\endcode
In the fixed-size case, no need to pass sizes:
\code
float array[9];
Map<Matrix3d> m(array);
Matrix3d::Map(array).setIdentity();
\endcode
</td></tr></table>

View File

@@ -11,7 +11,7 @@ namespace Eigen {
\section summary Executive summary
Using STL containers on \ref FixedSizeVectorizable "fixed-size vectorizable Eigen types" requires taking the following two steps:
Using STL containers on \ref FixedSizeVectorizable "fixed-size vectorizable Eigen types", or classes having members of such types, requires taking the following two steps:
\li A 16-byte-aligned allocator must be used. Eigen does provide one ready for use: aligned_allocator.
\li If you want to use the std::vector container, you need to \#include <Eigen/StdVector>.

View File

@@ -55,16 +55,17 @@ Note that here, Eigen::Vector2d is only used as an example, more generally the i
\section c2 Cause 2: STL Containers
If you use STL Containers such as std::vector, std::map, ..., with Eigen objects, like this,
If you use STL Containers such as std::vector, std::map, ..., with Eigen objects, or with classes containing Eigen objects, like this,
\code
std::vector<Eigen::Matrix2f> my_vector;
std::map<int, Eigen::Matrix2f> my_map;
struct my_class { ... Eigen::Matrix2f m; ... };
std::map<int, my_class> my_map;
\endcode
then you need to read this separate page: \ref StlContainers "Using STL Containers with Eigen".
Note that here, Eigen::Matrix2f is only used as an example, more generally the issue arises for all \ref FixedSizeVectorizable "fixed-size vectorizable Eigen types".
Note that here, Eigen::Matrix2f is only used as an example, more generally the issue arises for all \ref FixedSizeVectorizable "fixed-size vectorizable Eigen types" and \ref StructHavingEigenMembers "structures having such Eigen objects as member".
\section c3 Cause 3: Passing Eigen objects by value

7
eigen2.pc.in Normal file
View File

@@ -0,0 +1,7 @@
Name: Eigen2
Description: A C++ template library for linear algebra: vectors, matrices, and related algorithms
Requires:
Version: ${EIGEN_VERSION_NUMBER}
Libs:
Cflags: -I${INCLUDE_INSTALL_DIR}

View File

@@ -1,27 +1,19 @@
#!/bin/sh
# TODO : actually exit on exit, currently it only exit from the ()
# TODO : display error msg on stderr instead of stdout
# configuration
USER='orzel'
# You should call this script with USER set as you want, else some default
# will be used
USER=${USER:-'orzel'}
# step 1 : update
hg pull -u || (echo "update failed"; exit 1)
# step 2 : build
# step 1 : build
# todo if 'build is not there, create one:
#mkdir build
(cd build && cmake .. && make -j3 doc) || (echo "make failed"; exit 1)
(cd build && cmake .. && make -j3 doc) || { echo "make failed"; exit 1; }
#todo: n+1 where n = number of cpus
#step 3 : upload
BRANCH=`hg branch`
if [ $BRANCH == "default" ]
then
BRANCH='devel'
fi
#step 2 : upload
# (the '/' at the end of path are very important, see rsync documentation)
rsync -az build/doc/html/ $USER@ssh.tuxfamily.org:eigen/eigen.tuxfamily.org-web/htdocs/dox-$BRANCH/ || (echo "upload failed"; exit 1)
rsync -az build/doc/html/ $USER@ssh.tuxfamily.org:eigen/eigen.tuxfamily.org-web/htdocs/dox-2.0/ || { echo "upload failed"; exit 1; }
echo "Uploaded successfully"

View File

@@ -1,3 +1,7 @@
option(EIGEN_DEFAULT_TO_ROW_MAJOR "Use row-major as default matrix storage order" OFF)
if(EIGEN_DEFAULT_TO_ROW_MAJOR)
add_definitions("-DEIGEN_DEFAULT_TO_ROW_MAJOR")
endif()
find_package(GSL)
if(GSL_FOUND AND GSL_VERSION_MINOR LESS 9)
@@ -93,12 +97,20 @@ else(CMAKE_COMPILER_IS_GNUCXX)
endif(CMAKE_COMPILER_IS_GNUCXX)
option(EIGEN_NO_ASSERTION_CHECKING "Disable checking of assertions" OFF)
if(EIGEN_NO_ASSERTION_CHECKING)
add_definitions("-DEIGEN_NO_ASSERTION_CHECKING=1")
endif()
# similar to set_target_properties but append the property instead of overwriting it
macro(ei_add_target_property target prop value)
get_target_property(previous ${target} ${prop})
set_target_properties(${target} PROPERTIES ${prop} "${previous} ${value}")
if(previous MATCHES "NOTFOUND")
set_target_properties(${target} PROPERTIES ${prop} "${value}")
else()
set_target_properties(${target} PROPERTIES ${prop} "${previous} ${value}")
endif()
endmacro(ei_add_target_property)
@@ -134,13 +146,9 @@ macro(ei_add_test testname)
option(EIGEN_DEBUG_ASSERTS "Enable debuging of assertions" OFF)
if(EIGEN_DEBUG_ASSERTS)
set_target_properties(${targetname} PROPERTIES COMPILE_DEFINITIONS "-DEIGEN_DEBUG_ASSERTS=1")
set_target_properties(${targetname} PROPERTIES COMPILE_DEFINITIONS "EIGEN_DEBUG_ASSERTS=1")
endif(EIGEN_DEBUG_ASSERTS)
else(NOT EIGEN_NO_ASSERTION_CHECKING)
set_target_properties(${targetname} PROPERTIES COMPILE_DEFINITIONS "-DEIGEN_NO_ASSERTION_CHECKING=1")
endif(NOT EIGEN_NO_ASSERTION_CHECKING)
if(${ARGC} GREATER 1)
@@ -180,6 +188,7 @@ ei_add_test(meta)
ei_add_test(sizeof)
ei_add_test(dynalloc)
ei_add_test(nomalloc)
ei_add_test(first_aligned)
ei_add_test(mixingtypes)
ei_add_test(packetmath)
ei_add_test(unalignedassert)
@@ -200,7 +209,7 @@ ei_add_test(array)
ei_add_test(triangular)
ei_add_test(cholesky " " "${GSL_LIBRARIES}")
ei_add_test(lu ${EI_OFLAG})
ei_add_test(determinant)
ei_add_test(determinant ${EI_OFLAG})
ei_add_test(inverse)
ei_add_test(qr)
ei_add_test(eigensolver " " "${GSL_LIBRARIES}")
@@ -215,13 +224,17 @@ ei_add_test(newstdvector)
if(QT4_FOUND)
ei_add_test(qtvector " " "${QT_QTCORE_LIBRARY}")
endif(QT4_FOUND)
ei_add_test(sparse_vector)
ei_add_test(sparse_basic)
ei_add_test(sparse_solvers " " "${SPARSE_LIBS}")
ei_add_test(sparse_product)
if(NOT EIGEN_DEFAULT_TO_ROW_MAJOR)
ei_add_test(sparse_vector)
ei_add_test(sparse_basic)
ei_add_test(sparse_solvers " " "${SPARSE_LIBS}")
ei_add_test(sparse_product)
endif()
ei_add_test(swap)
ei_add_test(visitor)
ei_add_test(prec_inverse_4x4 ${EI_OFLAG})
# print a summary of the different options
message("************************************************************")
message("*** Eigen's unit tests configuration summary ***")
@@ -265,6 +278,12 @@ else(EIGEN_TEST_NO_EXPLICIT_VECTORIZATION)
message("Explicit vec: AUTO")
endif(EIGEN_TEST_NO_EXPLICIT_VECTORIZATION)
if(EIGEN_DEFAULT_TO_ROW_MAJOR)
message("Default order: Row-major")
else()
message("Default order: Column-major")
endif()
message("CXX: ${CMAKE_CXX_COMPILER}")
if(CMAKE_COMPILER_IS_GNUCXX)
execute_process(COMMAND ${CMAKE_CXX_COMPILER} --version COMMAND head -n 1 OUTPUT_VARIABLE EIGEN_CXX_VERSION_STRING OUTPUT_STRIP_TRAILING_WHITESPACE)

64
test/first_aligned.cpp Normal file
View File

@@ -0,0 +1,64 @@
// This file is part of Eigen, a lightweight C++ template library
// for linear algebra.
//
// Copyright (C) 2009 Benoit Jacob <jacob.benoit.1@gmail.com>
//
// Eigen is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 3 of the License, or (at your option) any later version.
//
// Alternatively, you can redistribute it and/or
// modify it under the terms of the GNU General Public License as
// published by the Free Software Foundation; either version 2 of
// the License, or (at your option) any later version.
//
// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU Lesser General Public
// License and a copy of the GNU General Public License along with
// Eigen. If not, see <http://www.gnu.org/licenses/>.
#include "main.h"
template<typename Scalar>
void test_first_aligned_helper(Scalar *array, int size)
{
const int packet_size = sizeof(Scalar) * ei_packet_traits<Scalar>::size;
VERIFY(((size_t(array) + sizeof(Scalar) * ei_alignmentOffset(array, size)) % packet_size) == 0);
}
template<typename Scalar>
void test_none_aligned_helper(Scalar *array, int size)
{
VERIFY(ei_packet_traits<Scalar>::size == 1 || ei_alignmentOffset(array, size) == size);
}
struct some_non_vectorizable_type { float x; };
void test_first_aligned()
{
EIGEN_ALIGN_128 float array_float[100];
test_first_aligned_helper(array_float, 50);
test_first_aligned_helper(array_float+1, 50);
test_first_aligned_helper(array_float+2, 50);
test_first_aligned_helper(array_float+3, 50);
test_first_aligned_helper(array_float+4, 50);
test_first_aligned_helper(array_float+5, 50);
EIGEN_ALIGN_128 double array_double[100];
test_first_aligned_helper(array_double, 50);
test_first_aligned_helper(array_double+1, 50);
test_first_aligned_helper(array_double+2, 50);
double *array_double_plus_4_bytes = (double*)(size_t(array_double)+4);
test_none_aligned_helper(array_double_plus_4_bytes, 50);
test_none_aligned_helper(array_double_plus_4_bytes+1, 50);
some_non_vectorizable_type array_nonvec[100];
test_first_aligned_helper(array_nonvec, 100);
test_none_aligned_helper(array_nonvec, 100);
}

View File

@@ -121,7 +121,8 @@ template<typename Scalar> void lines()
VERIFY_IS_APPROX(result, center);
// check conversions between two types of lines
CoeffsType converted_coeffs(HLine(PLine(line_u)).coeffs());
PLine pl(line_u); // gcc 3.3 will commit suicide if we don't name this variable
CoeffsType converted_coeffs(HLine(pl).coeffs());
converted_coeffs *= line_u.coeffs()(0)/converted_coeffs(0);
VERIFY(line_u.coeffs().isApprox(converted_coeffs));
}

View File

@@ -43,11 +43,9 @@ template<typename MatrixType> void inverse(const MatrixType& m)
mzero = MatrixType::Zero(rows, cols),
identity = MatrixType::Identity(rows, rows);
if (ei_is_same_type<RealScalar,float>::ret)
while(ei_abs(m1.determinant()) < RealScalar(0.1) && rows <= 8)
{
// let's build a more stable to inverse matrix
MatrixType a = MatrixType::Random(rows,cols);
m1 += m1 * m1.adjoint() + a * a.adjoint();
m1 = MatrixType::Random(rows, cols);
}
m2 = m1.inverse();

View File

@@ -1,7 +1,7 @@
// This file is part of Eigen, a lightweight C++ template library
// for linear algebra. Eigen itself is part of the KDE project.
//
// Copyright (C) 2006-2008 Benoit Jacob <jacob.benoit.1@gmail.com>
// Copyright (C) 2007-2010 Benoit Jacob <jacob.benoit.1@gmail.com>
//
// Eigen is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
@@ -24,7 +24,7 @@
#include "main.h"
template<typename VectorType> void map_class(const VectorType& m)
template<typename VectorType> void map_class_vector(const VectorType& m)
{
typedef typename VectorType::Scalar Scalar;
@@ -50,6 +50,34 @@ template<typename VectorType> void map_class(const VectorType& m)
delete[] array3;
}
template<typename MatrixType> void map_class_matrix(const MatrixType& m)
{
typedef typename MatrixType::Scalar Scalar;
int rows = m.rows(), cols = m.cols(), size = rows*cols;
// test Map.h
Scalar* array1 = ei_aligned_new<Scalar>(size);
for(int i = 0; i < size; i++) array1[i] = Scalar(1);
Scalar* array2 = ei_aligned_new<Scalar>(size);
for(int i = 0; i < size; i++) array2[i] = Scalar(1);
Scalar* array3 = new Scalar[size+1];
for(int i = 0; i < size+1; i++) array3[i] = Scalar(1);
Scalar* array3unaligned = size_t(array3)%16 == 0 ? array3+1 : array3;
Map<MatrixType, Aligned>(array1, rows, cols) = MatrixType::Ones(rows,cols);
Map<MatrixType>(array2, rows, cols) = Map<MatrixType>(array1, rows, cols);
Map<MatrixType>(array3unaligned, rows, cols) = Map<MatrixType>(array1, rows, cols);
MatrixType ma1 = Map<MatrixType>(array1, rows, cols);
MatrixType ma2 = Map<MatrixType, Aligned>(array2, rows, cols);
VERIFY_IS_APPROX(ma1, ma2);
MatrixType ma3 = Map<MatrixType>(array3unaligned, rows, cols);
VERIFY_IS_APPROX(ma1, ma3);
ei_aligned_delete(array1, size);
ei_aligned_delete(array2, size);
delete[] array3;
}
template<typename VectorType> void map_static_methods(const VectorType& m)
{
typedef typename VectorType::Scalar Scalar;
@@ -80,11 +108,17 @@ template<typename VectorType> void map_static_methods(const VectorType& m)
void test_map()
{
for(int i = 0; i < g_repeat; i++) {
CALL_SUBTEST( map_class(Matrix<float, 1, 1>()) );
CALL_SUBTEST( map_class(Vector4d()) );
CALL_SUBTEST( map_class(RowVector4f()) );
CALL_SUBTEST( map_class(VectorXcf(8)) );
CALL_SUBTEST( map_class(VectorXi(12)) );
CALL_SUBTEST( map_class_vector(Matrix<float, 1, 1>()) );
CALL_SUBTEST( map_class_vector(Vector4d()) );
CALL_SUBTEST( map_class_vector(RowVector4f()) );
CALL_SUBTEST( map_class_vector(VectorXcf(8)) );
CALL_SUBTEST( map_class_vector(VectorXi(12)) );
CALL_SUBTEST( map_class_matrix(Matrix<float, 1, 1>()) );
CALL_SUBTEST( map_class_matrix(Matrix4d()) );
CALL_SUBTEST( map_class_matrix(Matrix<float,3,5>()) );
CALL_SUBTEST( map_class_matrix(MatrixXcf(ei_random<int>(1,10),ei_random<int>(1,10))) );
CALL_SUBTEST( map_class_matrix(MatrixXi(ei_random<int>(1,10),ei_random<int>(1,10))) );
CALL_SUBTEST( map_static_methods(Matrix<double, 1, 1>()) );
CALL_SUBTEST( map_static_methods(Vector3f()) );

99
test/prec_inverse_4x4.cpp Normal file
View File

@@ -0,0 +1,99 @@
// This file is part of Eigen, a lightweight C++ template library
// for linear algebra.
//
// Copyright (C) 2009 Benoit Jacob <jacob.benoit.1@gmail.com>
//
// Eigen is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 3 of the License, or (at your option) any later version.
//
// Alternatively, you can redistribute it and/or
// modify it under the terms of the GNU General Public License as
// published by the Free Software Foundation; either version 2 of
// the License, or (at your option) any later version.
//
// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU Lesser General Public
// License and a copy of the GNU General Public License along with
// Eigen. If not, see <http://www.gnu.org/licenses/>.
#include "main.h"
#include <Eigen/LU>
#include <algorithm>
template<typename T> std::string type_name() { return "other"; }
template<> std::string type_name<float>() { return "float"; }
template<> std::string type_name<double>() { return "double"; }
template<> std::string type_name<int>() { return "int"; }
template<> std::string type_name<std::complex<float> >() { return "complex<float>"; }
template<> std::string type_name<std::complex<double> >() { return "complex<double>"; }
template<> std::string type_name<std::complex<int> >() { return "complex<int>"; }
#define EIGEN_DEBUG_VAR(x) std::cerr << #x << " = " << x << std::endl;
template<typename T> inline typename NumTraits<T>::Real epsilon()
{
return std::numeric_limits<typename NumTraits<T>::Real>::epsilon();
}
template<typename MatrixType> void inverse_permutation_4x4()
{
typedef typename MatrixType::Scalar Scalar;
typedef typename MatrixType::RealScalar RealScalar;
Vector4i indices(0,1,2,3);
for(int i = 0; i < 24; ++i)
{
MatrixType m = MatrixType::Zero();
m(indices(0),0) = 1;
m(indices(1),1) = 1;
m(indices(2),2) = 1;
m(indices(3),3) = 1;
MatrixType inv = m.inverse();
double error = double( (m*inv-MatrixType::Identity()).norm() / epsilon<Scalar>() );
VERIFY(error == 0.0);
std::next_permutation(indices.data(),indices.data()+4);
}
}
template<typename MatrixType> void inverse_general_4x4(int repeat)
{
typedef typename MatrixType::Scalar Scalar;
typedef typename MatrixType::RealScalar RealScalar;
double error_sum = 0., error_max = 0.;
for(int i = 0; i < repeat; ++i)
{
MatrixType m;
RealScalar absdet;
do {
m = MatrixType::Random();
absdet = ei_abs(m.determinant());
} while(absdet < 10 * epsilon<Scalar>());
MatrixType inv = m.inverse();
double error = double( (m*inv-MatrixType::Identity()).norm() * absdet / epsilon<Scalar>() );
error_sum += error;
error_max = std::max(error_max, error);
}
std::cerr << "inverse_general_4x4, Scalar = " << type_name<Scalar>() << std::endl;
double error_avg = error_sum / repeat;
EIGEN_DEBUG_VAR(error_avg);
EIGEN_DEBUG_VAR(error_max);
VERIFY(error_avg < (NumTraits<Scalar>::IsComplex ? 8.0 : 1.0));
VERIFY(error_max < (NumTraits<Scalar>::IsComplex ? 64.0 : 20.0));
}
void test_prec_inverse_4x4()
{
CALL_SUBTEST((inverse_permutation_4x4<Matrix4f>()));
CALL_SUBTEST(( inverse_general_4x4<Matrix4f>(200000 * g_repeat) ));
CALL_SUBTEST((inverse_permutation_4x4<Matrix<double,4,4,RowMajor> >()));
CALL_SUBTEST(( inverse_general_4x4<Matrix<double,4,4,RowMajor> >(200000 * g_repeat) ));
CALL_SUBTEST((inverse_permutation_4x4<Matrix4cf>()));
CALL_SUBTEST((inverse_general_4x4<Matrix4cf>(50000 * g_repeat)));
}

View File

@@ -46,7 +46,7 @@ template<typename MatrixType> void product(const MatrixType& m)
typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, MatrixType::RowsAtCompileTime> RowSquareMatrixType;
typedef Matrix<Scalar, MatrixType::ColsAtCompileTime, MatrixType::ColsAtCompileTime> ColSquareMatrixType;
typedef Matrix<Scalar, MatrixType::RowsAtCompileTime, MatrixType::ColsAtCompileTime,
MatrixType::Flags&RowMajorBit> OtherMajorMatrixType;
MatrixType::Options^RowMajor> OtherMajorMatrixType;
int rows = m.rows();
int cols = m.cols();
@@ -77,6 +77,7 @@ template<typename MatrixType> void product(const MatrixType& m)
// begin testing Product.h: only associativity for now
// (we use Transpose.h but this doesn't count as a test for it)
VERIFY_IS_APPROX((m1*m1.transpose())*m2, m1*(m1.transpose()*m2));
m3 = m1;
m3 *= m1.transpose() * m2;
@@ -137,6 +138,7 @@ template<typename MatrixType> void product(const MatrixType& m)
res2 = square2;
res2 += (m1.transpose() * m2).lazy();
VERIFY_IS_APPROX(res2, square2 + m1.transpose() * m2);
if (NumTraits<Scalar>::HasFloatingPoint && std::min(rows,cols)>1)
{
VERIFY(areNotApprox(res2,square2 + m2.transpose() * m1));

View File

@@ -26,10 +26,12 @@
#define EIGEN_WORK_AROUND_QT_BUG_CALLING_WRONG_OPERATOR_NEW_FIXED_IN_QT_4_5
#include "main.h"
#include <QtCore/QVector>
#include <Eigen/Geometry>
#include <Eigen/QtAlignedMalloc>
#include <QtCore/QVector>
template<typename MatrixType>
void check_qtvector_matrix(const MatrixType& m)
{

View File

@@ -68,12 +68,20 @@ template<typename Scalar> void sparse_solvers(int rows, int cols)
VERIFY_IS_APPROX(refMat2.template marked<LowerTriangular>().solveTriangular(vec2),
m2.template marked<LowerTriangular>().solveTriangular(vec3));
// lower - transpose
initSparse<Scalar>(density, refMat2, m2, ForceNonZeroDiag|MakeLowerTriangular, &zeroCoords, &nonzeroCoords);
VERIFY_IS_APPROX(refMat2.template marked<LowerTriangular>().transpose().solveTriangular(vec2),
m2.template marked<LowerTriangular>().transpose().solveTriangular(vec3));
// upper
initSparse<Scalar>(density, refMat2, m2, ForceNonZeroDiag|MakeUpperTriangular, &zeroCoords, &nonzeroCoords);
VERIFY_IS_APPROX(refMat2.template marked<UpperTriangular>().solveTriangular(vec2),
m2.template marked<UpperTriangular>().solveTriangular(vec3));
// TODO test row major
// upper - transpose
initSparse<Scalar>(density, refMat2, m2, ForceNonZeroDiag|MakeUpperTriangular, &zeroCoords, &nonzeroCoords);
VERIFY_IS_APPROX(refMat2.template marked<UpperTriangular>().transpose().solveTriangular(vec2),
m2.template marked<UpperTriangular>().transpose().solveTriangular(vec3));
}
// test LLT

View File

@@ -95,5 +95,8 @@ void test_svd()
// complex are not implemented yet
// CALL_SUBTEST( svd(MatrixXcd(6,6)) );
// CALL_SUBTEST( svd(MatrixXcf(3,3)) );
SVD<MatrixXf> s;
MatrixXf m = MatrixXf::Random(10,1);
s.compute(m);
}
}

View File

@@ -39,7 +39,7 @@
# VERSION=opensuse-11.1
# WORK_DIR=/home/gael/Coding/eigen2/cdash
# # get the last version of the script
# wget http://bitbucket.org/eigen/eigen2/raw/tip/test/testsuite.cmake -o $WORK_DIR/testsuite.cmake
# wget http://bitbucket.org/eigen/eigen/raw/tip/test/testsuite.cmake -o $WORK_DIR/testsuite.cmake
# COMMON="ctest -S $WORK_DIR/testsuite.cmake,EIGEN_WORK_DIR=$WORK_DIR,EIGEN_SITE=$SITE,EIGEN_MODE=$1,EIGEN_BUILD_STRING=$OS_VERSION-$ARCH"
# $COMMON-gcc-3.4.6,EIGEN_CXX=g++-3.4
# $COMMON-gcc-4.0.1,EIGEN_CXX=g++-4.0.1
@@ -133,7 +133,7 @@ endif(NOT EIGEN_MODE)
## mandatory variables (the default should be ok in most cases):
SET (CTEST_CVS_COMMAND "hg")
SET (CTEST_CVS_CHECKOUT "${CTEST_CVS_COMMAND} clone -r 2.0 http://bitbucket.org/eigen/eigen2 \"${CTEST_SOURCE_DIRECTORY}\"")
SET (CTEST_CVS_CHECKOUT "${CTEST_CVS_COMMAND} clone -r 2.0 http://bitbucket.org/eigen/eigen \"${CTEST_SOURCE_DIRECTORY}\"")
# which ctest command to use for running the dashboard
SET (CTEST_COMMAND "${EIGEN_CMAKE_DIR}ctest -D ${EIGEN_MODE}")

View File

@@ -44,12 +44,21 @@ void test_vectorization_logic()
#ifdef EIGEN_VECTORIZE
#ifdef EIGEN_DEFAULT_TO_ROW_MAJOR
VERIFY(test_assign(Vector4f(),Vector4f(),
LinearVectorization,CompleteUnrolling));
VERIFY(test_assign(Vector4f(),Vector4f()+Vector4f(),
LinearVectorization,CompleteUnrolling));
VERIFY(test_assign(Vector4f(),Vector4f().cwise() * Vector4f(),
LinearVectorization,CompleteUnrolling));
#else
VERIFY(test_assign(Vector4f(),Vector4f(),
InnerVectorization,CompleteUnrolling));
VERIFY(test_assign(Vector4f(),Vector4f()+Vector4f(),
InnerVectorization,CompleteUnrolling));
VERIFY(test_assign(Vector4f(),Vector4f().cwise() * Vector4f(),
InnerVectorization,CompleteUnrolling));
#endif
VERIFY(test_assign(Matrix4f(),Matrix4f(),
InnerVectorization,CompleteUnrolling));
@@ -92,8 +101,10 @@ void test_vectorization_logic()
VERIFY(test_sum(Matrix<float,16,16>().block<4,4>(1,2),
NoVectorization,CompleteUnrolling));
#ifndef EIGEN_DEFAULT_TO_ROW_MAJOR
VERIFY(test_sum(Matrix<float,16,16>().block<8,1>(1,2),
LinearVectorization,CompleteUnrolling));
#endif
VERIFY(test_sum(Matrix<double,7,3>(),
NoVectorization,CompleteUnrolling));

View File

@@ -40,7 +40,7 @@ template<typename MatrixType> void matrixVisitor(const MatrixType& p)
m(i) = ei_random<Scalar>();
Scalar minc = Scalar(1000), maxc = Scalar(-1000);
int minrow,mincol,maxrow,maxcol;
int minrow=0,mincol=0,maxrow=0,maxcol=0;
for(int j = 0; j < cols; j++)
for(int i = 0; i < rows; i++)
{
@@ -86,7 +86,7 @@ template<typename VectorType> void vectorVisitor(const VectorType& w)
v(i) = ei_random<Scalar>();
Scalar minc = Scalar(1000), maxc = Scalar(-1000);
int minidx,maxidx;
int minidx=0,maxidx=0;
for(int i = 0; i < size; i++)
{
if(v(i) < minc)