[eigen] patches for eigen3.0-beta2 |
[ Thread Index |
Date Index
| More lists.tuxfamily.org/eigen Archives
]
In an attempt to clean up the ROS eigen3 meta package, I was wondering if we can integrate the following 3 patches into
Eigen trunk, such that when beta2 comes out we don't have to apply them in ROS.
These are:
* tridiag_3x3.patch: patch provided by Gael that fixes the float eigen decomposition precision (without it MatrixXd
matrices work well but MatrixXf ones do not)
* gcc43.patch: patch provided by Benoit that makes GCC 4.3 on Ubuntu 9.04 happy (without it Eigen3 fails to compile)
* warning.patch: patch that silences a compilation warning in Ubuntu 10.04
Thanks a lot in advance!
Cheers,
Radu.
--- Eigen/src/Eigenvalues/Tridiagonalization.h 2010-08-31 09:02:15.000000000 -0700
+++ Eigen/src/Eigenvalues/Tridiagonalization.h 2010-08-31 09:02:30.000000000 -0700
@@ -433,8 +433,8 @@
void ei_tridiagonalization_inplace(MatrixType& mat, DiagonalType& diag, SubDiagonalType& subdiag, bool extractQ)
{
typedef typename MatrixType::Index Index;
- Index n = mat.rows();
- ei_assert(mat.cols()==n && diag.size()==n && subdiag.size()==n-1);
+ //Index n = mat.rows();
+ ei_assert(mat.cols()==mat.rows() && diag.size()==mat.rows() && subdiag.size()==mat.rows()-1);
ei_tridiagonalization_inplace_selector<MatrixType>::run(mat, diag, subdiag, extractQ);
}
--- Eigen/src/Eigenvalues/Tridiagonalization.h Fri Jul 16 14:02:20 2010 +0200
+++ Eigen/src/Eigenvalues/Tridiagonalization.h Fri Jul 16 16:38:58 2010 +0200
@@ -458,7 +458,7 @@
};
/** \internal
- * Specialization for 3x3 matrices.
+ * Specialization for 3x3 real matrices.
* Especially useful for plane fitting.
*/
template<typename MatrixType>
@@ -472,7 +472,7 @@
{
diag[0] = ei_real(mat(0,0));
RealScalar v1norm2 = ei_abs2(mat(2,0));
- if (ei_isMuchSmallerThan(v1norm2, RealScalar(1)))
+ if(v1norm2 == RealScalar(0) && ei_imag(mat(1,0))==RealScalar(0))
{
diag[1] = ei_real(mat(1,1));
diag[2] = ei_real(mat(2,2));
@@ -495,8 +495,8 @@
if (extractQ)
{
mat << 1, 0, 0,
- 0, m01, m02,
- 0, m02, -m01;
+ 0, m01, m02,
+ 0, m02, -m01;
}
}
}
--- Eigen/src/Eigenvalues/Tridiagonalization.h Fri Jul 16 16:38:58 2010 +0200
+++ Eigen/src/Eigenvalues/Tridiagonalization.h Fri Jul 16 18:22:00 2010 +0200
@@ -384,7 +384,9 @@
}
// forward declaration, implementation at the end of this file
-template<typename MatrixType, int Size=MatrixType::ColsAtCompileTime>
+template<typename MatrixType,
+ int Size=MatrixType::ColsAtCompileTime,
+ bool IsComplex=NumTraits<typename MatrixType::Scalar>::IsComplex>
struct ei_tridiagonalization_inplace_selector;
/** \brief Performs a full tridiagonalization in place
@@ -439,7 +441,7 @@
/** \internal
* General full tridiagonalization
*/
-template<typename MatrixType, int Size>
+template<typename MatrixType, int Size, bool IsComplex>
struct ei_tridiagonalization_inplace_selector
{
typedef typename Tridiagonalization<MatrixType>::CoeffVectorType CoeffVectorType;
@@ -462,7 +464,7 @@
* Especially useful for plane fitting.
*/
template<typename MatrixType>
-struct ei_tridiagonalization_inplace_selector<MatrixType,3>
+struct ei_tridiagonalization_inplace_selector<MatrixType,3,false>
{
typedef typename MatrixType::Scalar Scalar;
typedef typename MatrixType::RealScalar RealScalar;
@@ -470,14 +472,14 @@
template<typename DiagonalType, typename SubDiagonalType>
static void run(MatrixType& mat, DiagonalType& diag, SubDiagonalType& subdiag, bool extractQ)
{
- diag[0] = ei_real(mat(0,0));
+ diag[0] = mat(0,0);
RealScalar v1norm2 = ei_abs2(mat(2,0));
- if(v1norm2 == RealScalar(0) && ei_imag(mat(1,0))==RealScalar(0))
+ if(v1norm2 == RealScalar(0))
{
- diag[1] = ei_real(mat(1,1));
- diag[2] = ei_real(mat(2,2));
- subdiag[0] = ei_real(mat(1,0));
- subdiag[1] = ei_real(mat(2,1));
+ diag[1] = mat(1,1);
+ diag[2] = mat(2,2);
+ subdiag[0] = mat(1,0);
+ subdiag[1] = mat(2,1);
if (extractQ)
mat.setIdentity();
}
@@ -485,13 +487,13 @@
{
RealScalar beta = ei_sqrt(ei_abs2(mat(1,0)) + v1norm2);
RealScalar invBeta = RealScalar(1)/beta;
- Scalar m01 = ei_conj(mat(1,0)) * invBeta;
- Scalar m02 = ei_conj(mat(2,0)) * invBeta;
- Scalar q = RealScalar(2)*m01*ei_conj(mat(2,1)) + m02*(mat(2,2) - mat(1,1));
- diag[1] = ei_real(mat(1,1) + m02*q);
- diag[2] = ei_real(mat(2,2) - m02*q);
+ Scalar m01 = mat(1,0) * invBeta;
+ Scalar m02 = mat(2,0) * invBeta;
+ Scalar q = RealScalar(2)*m01*mat(2,1) + m02*(mat(2,2) - mat(1,1));
+ diag[1] = mat(1,1) + m02*q;
+ diag[2] = mat(2,2) - m02*q;
subdiag[0] = beta;
- subdiag[1] = ei_real(ei_conj(mat(2,1)) - m01 * q);
+ subdiag[1] = mat(2,1) - m01 * q;
if (extractQ)
{
mat << 1, 0, 0,
--- Eigen/src/Core/MatrixBase.h
+++ Eigen/src/Core/MatrixBase.h
@@ -333,29 +333,33 @@ template<typename Derived> class MatrixB
/////////// Geometry module ///////////
template<typename OtherDerived>
PlainObject cross(const MatrixBase<OtherDerived>& other) const;
template<typename OtherDerived>
PlainObject cross3(const MatrixBase<OtherDerived>& other) const;
PlainObject unitOrthogonal(void) const;
Matrix<Scalar,3,1> eulerAngles(Index a0, Index a1, Index a2) const;
- const ScalarMultipleReturnType operator*(const UniformScaling<Scalar>& s) const;
+ ScalarMultipleReturnType operator*(const UniformScaling<Scalar>& s) const;
enum {
SizeMinusOne = SizeAtCompileTime==Dynamic ? Dynamic : SizeAtCompileTime-1
};
typedef Block<Derived,
ei_traits<Derived>::ColsAtCompileTime==1 ? SizeMinusOne : 1,
ei_traits<Derived>::ColsAtCompileTime==1 ? 1 : SizeMinusOne> StartMinusOne;
typedef CwiseUnaryOp<ei_scalar_quotient1_op<typename ei_traits<Derived>::Scalar>,
StartMinusOne > HNormalizedReturnType;
- const HNormalizedReturnType hnormalized() const;
- typedef Homogeneous<Derived,MatrixBase<Derived>::ColsAtCompileTime==1?Vertical:Horizontal> HomogeneousReturnType;
- const HomogeneousReturnType homogeneous() const;
+ HNormalizedReturnType hnormalized() const;
+
+ // put this as separate enum value to work around possible GCC 4.3 bug (?)
+ enum { HomogeneousReturnTypeDirection = ColsAtCompileTime==1?Vertical:Horizontal };
+ typedef Homogeneous<Derived, HomogeneousReturnTypeDirection> HomogeneousReturnType;
+
+ HomogeneousReturnType homogeneous() const;
////////// Householder module ///////////
void makeHouseholderInPlace(Scalar& tau, RealScalar& beta);
template<typename EssentialPart>
void makeHouseholder(EssentialPart& essential,
Scalar& tau, RealScalar& beta) const;
template<typename EssentialPart>
--- Eigen/src/Core/VectorwiseOp.h
+++ Eigen/src/Core/VectorwiseOp.h
@@ -456,17 +456,17 @@ template<typename ExpressionType, int Di
operator-(const DenseBase<OtherDerived>& other) const
{
EIGEN_STATIC_ASSERT_VECTOR_ONLY(OtherDerived);
return m_matrix - extendedTo(other.derived());
}
/////////// Geometry module ///////////
- const Homogeneous<ExpressionType,Direction> homogeneous() const;
+ Homogeneous<ExpressionType,Direction> homogeneous() const;
typedef typename ExpressionType::PlainObject CrossReturnType;
template<typename OtherDerived>
const CrossReturnType cross(const MatrixBase<OtherDerived>& other) const;
enum {
HNormalized_Size = Direction==Vertical ? ei_traits<ExpressionType>::RowsAtCompileTime
: ei_traits<ExpressionType>::ColsAtCompileTime,
@@ -484,17 +484,17 @@ template<typename ExpressionType, int Di
HNormalized_Factors;
typedef CwiseBinaryOp<ei_scalar_quotient_op<typename ei_traits<ExpressionType>::Scalar>,
HNormalized_Block,
Replicate<HNormalized_Factors,
Direction==Vertical ? HNormalized_SizeMinusOne : 1,
Direction==Horizontal ? HNormalized_SizeMinusOne : 1> >
HNormalizedReturnType;
- const HNormalizedReturnType hnormalized() const;
+ HNormalizedReturnType hnormalized() const;
protected:
ExpressionTypeNested m_matrix;
};
/** \returns a VectorwiseOp wrapper of *this providing additional partial reduction operations
*
* Example: \include MatrixBase_colwise.cpp
--- Eigen/src/Core/util/StaticAssert.h
+++ Eigen/src/Core/util/StaticAssert.h
@@ -98,22 +98,22 @@
};
// Specialized implementation for MSVC to avoid "conditional
// expression is constant" warnings. This implementation doesn't
// appear to work under GCC, hence the multiple implementations.
#ifdef _MSC_VER
#define EIGEN_STATIC_ASSERT(CONDITION,MSG) \
- {Eigen3::ei_static_assert<CONDITION ? true : false>::MSG;}
+ {Eigen3::ei_static_assert<(CONDITION)>::MSG;}
#else
#define EIGEN_STATIC_ASSERT(CONDITION,MSG) \
- if (Eigen3::ei_static_assert<CONDITION ? true : false>::MSG) {}
+ if (Eigen3::ei_static_assert<(CONDITION)>::MSG) {}
#endif
#endif // not CXX0X
#else // EIGEN_NO_STATIC_ASSERT
#define EIGEN_STATIC_ASSERT(CONDITION,MSG) ei_assert((CONDITION) && #MSG);
--- Eigen/src/Geometry/Homogeneous.h
+++ Eigen/src/Geometry/Homogeneous.h
@@ -130,48 +130,48 @@ template<typename MatrixType,int _Direct
* \only_for_vectors
*
* Example: \include MatrixBase_homogeneous.cpp
* Output: \verbinclude MatrixBase_homogeneous.out
*
* \sa class Homogeneous
*/
template<typename Derived>
-inline const typename MatrixBase<Derived>::HomogeneousReturnType
+inline typename MatrixBase<Derived>::HomogeneousReturnType
MatrixBase<Derived>::homogeneous() const
{
EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived);
return derived();
}
/** \geometry_module
*
* \returns a matrix expression of homogeneous column (or row) vectors
*
* Example: \include VectorwiseOp_homogeneous.cpp
* Output: \verbinclude VectorwiseOp_homogeneous.out
*
* \sa MatrixBase::homogeneous() */
template<typename ExpressionType, int Direction>
-inline const Homogeneous<ExpressionType,Direction>
+inline Homogeneous<ExpressionType,Direction>
VectorwiseOp<ExpressionType,Direction>::homogeneous() const
{
return _expression();
}
/** \geometry_module
*
* \returns an expression of the homogeneous normalized vector of \c *this
*
* Example: \include MatrixBase_hnormalized.cpp
* Output: \verbinclude MatrixBase_hnormalized.out
*
* \sa VectorwiseOp::hnormalized() */
template<typename Derived>
-inline const typename MatrixBase<Derived>::HNormalizedReturnType
+inline typename MatrixBase<Derived>::HNormalizedReturnType
MatrixBase<Derived>::hnormalized() const
{
EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived);
return StartMinusOne(derived(),0,0,
ColsAtCompileTime==1?size()-1:1,
ColsAtCompileTime==1?1:size()-1) / coeff(size()-1);
}
@@ -187,17 +187,17 @@ MatrixBase<Derived>::hnormalized() const
*
* \returns an expression of the homogeneous normalized vector of \c *this
*
* Example: \include DirectionWise_hnormalized.cpp
* Output: \verbinclude DirectionWise_hnormalized.out
*
* \sa MatrixBase::hnormalized() */
template<typename ExpressionType, int Direction>
-inline const typename VectorwiseOp<ExpressionType,Direction>::HNormalizedReturnType
+inline typename VectorwiseOp<ExpressionType,Direction>::HNormalizedReturnType
VectorwiseOp<ExpressionType,Direction>::hnormalized() const
{
return HNormalized_Block(_expression(),0,0,
Direction==Vertical ? _expression().rows()-1 : _expression().rows(),
Direction==Horizontal ? _expression().cols()-1 : _expression().cols()).cwiseQuotient(
Replicate<HNormalized_Factors,
Direction==Vertical ? HNormalized_SizeMinusOne : 1,
Direction==Horizontal ? HNormalized_SizeMinusOne : 1>
--- Eigen/src/Geometry/Scaling.h
+++ Eigen/src/Geometry/Scaling.h
@@ -110,17 +110,17 @@ public:
bool isApprox(const UniformScaling& other, typename NumTraits<Scalar>::Real prec = NumTraits<Scalar>::dummy_precision()) const
{ return ei_isApprox(m_factor, other.factor(), prec); }
};
/** Concatenates a linear transformation matrix and a uniform scaling */
// NOTE this operator is defiend in MatrixBase and not as a friend function
// of UniformScaling to fix an internal crash of Intel's ICC
-template<typename Derived> const typename MatrixBase<Derived>::ScalarMultipleReturnType
+template<typename Derived> typename MatrixBase<Derived>::ScalarMultipleReturnType
MatrixBase<Derived>::operator*(const UniformScaling<Scalar>& s) const
{ return derived() * s.factor(); }
/** Constructs a uniform scaling from scale factor \a s */
static inline UniformScaling<float> Scaling(float s) { return UniformScaling<float>(s); }
/** Constructs a uniform scaling from scale factor \a s */
static inline UniformScaling<double> Scaling(double s) { return UniformScaling<double>(s); }
/** Constructs a uniform scaling from scale factor \a s */
--- Eigen/Core 2010-09-02 22:11:52.000000000 -0700
+++ Eigen/Core 2010-09-03 13:49:57.000000000 -0700
@@ -26,6 +26,8 @@
#ifndef EIGEN_CORE_H
#define EIGEN_CORE_H
+#define EIGEN_NO_STATIC_ASSERT
+
// first thing Eigen does: prevent MSVC from committing suicide
#include "src/Core/util/DisableMSVCWarnings.h"