aboutsummaryrefslogtreecommitdiffhomepage
diff options
context:
space:
mode:
authorGravatar Gael Guennebaud <g.gael@free.fr>2009-03-11 08:45:53 +0000
committerGravatar Gael Guennebaud <g.gael@free.fr>2009-03-11 08:45:53 +0000
commitf697ea6d306f05690abcb7dbfa2d6f8856f16ddb (patch)
treeb5414765429d0932cb0225edee14465f72d1ac46
parent14691d68363b4903621400d670d0b3649cb5b02a (diff)
fix a few compilation errors and warnings (ICC)
-rw-r--r--Eigen/src/Array/PartialRedux.h18
-rw-r--r--Eigen/src/Array/Replicate.h24
-rw-r--r--Eigen/src/Core/SolveTriangular.h4
-rw-r--r--Eigen/src/Core/products/GeneralMatrixMatrix.h50
-rw-r--r--Eigen/src/Geometry/Transform.h2
5 files changed, 49 insertions, 49 deletions
diff --git a/Eigen/src/Array/PartialRedux.h b/Eigen/src/Array/PartialRedux.h
index 0f868677c..ee956bb1f 100644
--- a/Eigen/src/Array/PartialRedux.h
+++ b/Eigen/src/Array/PartialRedux.h
@@ -280,22 +280,22 @@ template<typename ExpressionType, int Direction> class PartialRedux
{
return Reverse<ExpressionType, Direction>( _expression() );
}
-
+
+ const Replicate<ExpressionType,Direction==Vertical?Dynamic:1,Direction==Horizontal?Dynamic:1>
+ replicate(int factor) const;
+
template<int Factor>
- const Replicate<ExpressionType,Direction==Vertical?Factor:1,Direction==Horizontal?Factor:1>
+ const Replicate<ExpressionType,(Direction==Vertical?Factor:1),(Direction==Horizontal?Factor:1)>
replicate(int factor = Factor) const;
-
+
/////////// Geometry module ///////////
const Homogeneous<ExpressionType,Direction> homogeneous() const;
-
- const Replicate<ExpressionType,Direction==Vertical?Dynamic:1,Direction==Horizontal?Dynamic:1>
- replicate(int factor) const;
typedef typename ExpressionType::PlainMatrixType 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,
@@ -311,13 +311,13 @@ template<typename ExpressionType, int Direction> class PartialRedux
Direction==Vertical ? 1 : int(ei_traits<ExpressionType>::RowsAtCompileTime),
Direction==Horizontal ? 1 : int(ei_traits<ExpressionType>::ColsAtCompileTime)>
HNormalized_Factors;
- typedef CwiseBinaryOp<ei_scalar_quotient_op<typename ei_traits<ExpressionType>::Scalar>,
+ typedef CwiseBinaryOp<ei_scalar_quotient_op<typename ei_traits<ExpressionType>::Scalar>,
NestByValue<HNormalized_Block>,
NestByValue<Replicate<NestByValue<HNormalized_Factors>,
Direction==Vertical ? HNormalized_SizeMinusOne : 1,
Direction==Horizontal ? HNormalized_SizeMinusOne : 1> > >
HNormalizedReturnType;
-
+
const HNormalizedReturnType hnormalized() const;
protected:
diff --git a/Eigen/src/Array/Replicate.h b/Eigen/src/Array/Replicate.h
index d3220c154..87ff5eaa8 100644
--- a/Eigen/src/Array/Replicate.h
+++ b/Eigen/src/Array/Replicate.h
@@ -25,7 +25,7 @@
#ifndef EIGEN_REPLICATE_H
#define EIGEN_REPLICATE_H
-/** \nonstableyet
+/** \nonstableyet
* \class Replicate
*
* \brief Expression of the multiple replication of a matrix or vector
@@ -70,11 +70,11 @@ template<typename MatrixType,int RowFactor,int ColFactor> class Replicate
EIGEN_GENERIC_PUBLIC_INTERFACE(Replicate)
inline Replicate(const MatrixType& matrix)
- : m_matrix(matrix)
+ : m_matrix(matrix), m_rowFactor(RowFactor), m_colFactor(ColFactor)
{
ei_assert(RowFactor!=Dynamic && ColFactor!=Dynamic);
}
-
+
inline Replicate(const MatrixType& matrix, int rowFactor, int colFactor)
: m_matrix(matrix), m_rowFactor(rowFactor), m_colFactor(colFactor)
{}
@@ -93,9 +93,9 @@ template<typename MatrixType,int RowFactor,int ColFactor> class Replicate
const ei_int_if_dynamic<ColFactor> m_colFactor;
};
-/** \nonstableyet
+/** \nonstableyet
* \return an expression of the replication of \c *this
- *
+ *
* Example: \include MatrixBase_replicate.cpp
* Output: \verbinclude MatrixBase_replicate.out
*
@@ -109,9 +109,9 @@ MatrixBase<Derived>::replicate() const
return derived();
}
-/** \nonstableyet
+/** \nonstableyet
* \return an expression of the replication of \c *this
- *
+ *
* Example: \include MatrixBase_replicate_int_int.cpp
* Output: \verbinclude MatrixBase_replicate_int_int.out
*
@@ -124,9 +124,9 @@ MatrixBase<Derived>::replicate(int rowFactor,int colFactor) const
return Replicate<Derived,Dynamic,Dynamic>(derived(),rowFactor,colFactor);
}
-/** \nonstableyet
+/** \nonstableyet
* \return an expression of the replication of each column (or row) of \c *this
- *
+ *
* Example: \include DirectionWise_replicate_int.cpp
* Output: \verbinclude DirectionWise_replicate_int.out
*
@@ -140,9 +140,9 @@ PartialRedux<ExpressionType,Direction>::replicate(int factor) const
(_expression(),Direction==Vertical?factor:1,Direction==Horizontal?factor:1);
}
-/** \nonstableyet
+/** \nonstableyet
* \return an expression of the replication of each column (or row) of \c *this
- *
+ *
* Example: \include DirectionWise_replicate.cpp
* Output: \verbinclude DirectionWise_replicate.out
*
@@ -156,5 +156,5 @@ PartialRedux<ExpressionType,Direction>::replicate(int factor) const
return Replicate<ExpressionType,Direction==Vertical?Factor:1,Direction==Horizontal?Factor:1>
(_expression(),Direction==Vertical?factor:1,Direction==Horizontal?factor:1);
}
-
+
#endif // EIGEN_REPLICATE_H
diff --git a/Eigen/src/Core/SolveTriangular.h b/Eigen/src/Core/SolveTriangular.h
index 12fb0e1d1..f8fe74310 100644
--- a/Eigen/src/Core/SolveTriangular.h
+++ b/Eigen/src/Core/SolveTriangular.h
@@ -33,8 +33,8 @@ template<typename Lhs, typename Rhs,
? LowerTriangular
: (int(Lhs::Flags) & UpperTriangularBit)
? UpperTriangular
- : -1,
- int StorageOrder = ei_is_part<Lhs>::value ? -1 // this is to solve ambiguous specializations
+ : 0xffffff,
+ int StorageOrder = ei_is_part<Lhs>::value ? 0xffffff // this is to solve ambiguous specializations
: int(Lhs::Flags) & (RowMajorBit|SparseBit)
>
struct ei_solve_triangular_selector;
diff --git a/Eigen/src/Core/products/GeneralMatrixMatrix.h b/Eigen/src/Core/products/GeneralMatrixMatrix.h
index 3543de9b8..9dedb7f14 100644
--- a/Eigen/src/Core/products/GeneralMatrixMatrix.h
+++ b/Eigen/src/Core/products/GeneralMatrixMatrix.h
@@ -69,7 +69,7 @@ static void ei_cache_friendly_product(
typedef typename ei_packet_traits<Scalar>::type PacketType;
-
+
#ifndef EIGEN_USE_ALT_PRODUCT
@@ -83,31 +83,31 @@ static void ei_cache_friendly_product(
// register block size along the N direction
nr = HalfRegisterCount/2,
-
+
// register block size along the M direction
mr = 2 * PacketSize,
-
+
// max cache block size along the K direction
Max_kc = ei_L2_block_traits<EIGEN_TUNE_FOR_CPU_CACHE_SIZE,Scalar>::width,
-
+
// max cache block size along the M direction
Max_mc = 2*Max_kc
};
int kc = std::min<int>(Max_kc,depth); // cache block size along the K direction
int mc = std::min<int>(Max_mc,rows); // cache block size along the M direction
-
+
Scalar* blockA = ei_aligned_stack_new(Scalar, kc*mc);
Scalar* blockB = ei_aligned_stack_new(Scalar, kc*cols*PacketSize);
-
+
// number of columns which can be processed by packet of nr columns
int packet_cols = (cols/nr)*nr;
-
+
// GEMM_VAR1
for(int k2=0; k2<depth; k2+=kc)
{
const int actual_kc = std::min(k2+kc,depth)-k2;
-
+
// we have selected one row panel of rhs and one column panel of lhs
// pack rhs's panel into a sequential chunk of memory
// and expand each coeff to a constant packet for further reuse
@@ -132,12 +132,12 @@ static void ei_cache_friendly_product(
}
}
}
-
+
// => GEPP_VAR1
for(int i2=0; i2<rows; i2+=mc)
{
const int actual_mc = std::min(i2+mc,rows)-i2;
-
+
// We have selected a mc x kc block of lhs
// Let's pack it in a clever order for further purely sequential access
int count = 0;
@@ -176,11 +176,11 @@ static void ei_cache_friendly_product(
{
const Scalar* blA = &blockA[i*actual_kc];
#ifdef EIGEN_VECTORIZE_SSE
- _mm_prefetch(&blA[0], _MM_HINT_T0);
+ _mm_prefetch((const char*)(&blA[0]), _MM_HINT_T0);
#endif
-
+
// TODO move the res loads to the stores
-
+
// gets res block as register
PacketType C0, C1, C2, C3, C4, C5, C6, C7;
C0 = ei_ploadu(&res[(j2+0)*resStride + i2 + i]);
@@ -191,7 +191,7 @@ static void ei_cache_friendly_product(
C5 = ei_ploadu(&res[(j2+1)*resStride + i2 + i + PacketSize]);
if(nr==4) C6 = ei_ploadu(&res[(j2+2)*resStride + i2 + i + PacketSize]);
if(nr==4) C7 = ei_ploadu(&res[(j2+3)*resStride + i2 + i + PacketSize]);
-
+
// performs "inner" product
// TODO let's check wether the flowing peeled loop could not be
// optimized via optimal prefetching from one loop to the other
@@ -258,7 +258,7 @@ static void ei_cache_friendly_product(
if(nr==4) C6 = ei_pmadd(B2, A1, C6);
if(nr==4) C3 = ei_pmadd(B3, A0, C3);
if(nr==4) C7 = ei_pmadd(B3, A1, C7);
-
+
blB += 4*nr*PacketSize;
blA += 4*mr;
}
@@ -281,7 +281,7 @@ static void ei_cache_friendly_product(
if(nr==4) C6 = ei_pmadd(B2, A1, C6);
if(nr==4) C3 = ei_pmadd(B3, A0, C3);
if(nr==4) C7 = ei_pmadd(B3, A1, C7);
-
+
blB += nr*PacketSize;
blA += mr;
}
@@ -299,9 +299,9 @@ static void ei_cache_friendly_product(
{
const Scalar* blA = &blockA[i*actual_kc];
#ifdef EIGEN_VECTORIZE_SSE
- _mm_prefetch(&blA[0], _MM_HINT_T0);
+ _mm_prefetch((const char*)(&blA[0]), _MM_HINT_T0);
#endif
-
+
// gets a 1 x nr res block as registers
Scalar C0(0), C1(0), C2(0), C3(0);
const Scalar* blB = &blockB[j2*actual_kc*PacketSize];
@@ -318,7 +318,7 @@ static void ei_cache_friendly_product(
C1 += B1 * A0;
if(nr==4) C2 += B2 * A0;
if(nr==4) C3 += B3 * A0;
-
+
blB += nr*PacketSize;
}
res[(j2+0)*resStride + i2 + i] += C0;
@@ -344,10 +344,10 @@ static void ei_cache_friendly_product(
}
}
}
-
+
ei_aligned_stack_delete(Scalar, blockA, kc*mc);
ei_aligned_stack_delete(Scalar, blockB, kc*cols*PacketSize);
-
+
#else // alternate product from cylmor
enum {
@@ -364,18 +364,18 @@ static void ei_cache_friendly_product(
// maximal size of the blocks fitted in L2 cache
MaxL2BlockSize = ei_L2_block_traits<EIGEN_TUNE_FOR_CPU_CACHE_SIZE,Scalar>::width
};
-
+
const bool resIsAligned = (PacketSize==1) || (((resStride%PacketSize) == 0) && (size_t(res)%16==0));
const int remainingSize = depth % PacketSize;
const int size = depth - remainingSize; // third dimension of the product clamped to packet boundaries
-
+
const int l2BlockRows = MaxL2BlockSize > rows ? rows : 512;
const int l2BlockCols = MaxL2BlockSize > cols ? cols : 128;
const int l2BlockSize = MaxL2BlockSize > size ? size : 256;
const int l2BlockSizeAligned = (1 + std::max(l2BlockSize,l2BlockCols)/PacketSize)*PacketSize;
const bool needRhsCopy = (PacketSize>1) && ((rhsStride%PacketSize!=0) || (size_t(rhs)%16!=0));
-
+
Scalar* EIGEN_RESTRICT block = new Scalar[l2BlockRows*size];
// for(int i=0; i<l2BlockRows*l2BlockSize; ++i)
// block[i] = 0;
@@ -491,7 +491,7 @@ static void ei_cache_friendly_product(
delete[] block;
#endif
-
+
}
#endif // EIGEN_EXTERN_INSTANTIATIONS
diff --git a/Eigen/src/Geometry/Transform.h b/Eigen/src/Geometry/Transform.h
index a8dc164e4..47e4b2925 100644
--- a/Eigen/src/Geometry/Transform.h
+++ b/Eigen/src/Geometry/Transform.h
@@ -386,7 +386,7 @@ public:
Transform& fromPositionOrientationScale(const MatrixBase<PositionDerived> &position,
const OrientationType& orientation, const MatrixBase<ScaleDerived> &scale);
- inline const MatrixType inverse(TransformTraits traits = Mode) const;
+ inline const MatrixType inverse(TransformTraits traits = (TransformTraits)Mode) const;
/** \returns a const pointer to the column major internal matrix */
const Scalar* data() const { return m_matrix.data(); }