Merged eigen/eigen into default

This commit is contained in:
Benoit Steiner
2017-02-26 10:02:14 -08:00
19 changed files with 198 additions and 77 deletions

View File

@@ -54,7 +54,7 @@ struct is_input_scalar<Sizes<> > {
static const bool value = true;
};
#ifndef EIGEN_EMULATE_CXX11_META_H
template <typename std::size_t... Indices>
template <typename std::ptrdiff_t... Indices>
struct is_input_scalar<Sizes<Indices...> > {
static const bool value = (Sizes<Indices...>::total_size == 1);
};

View File

@@ -126,7 +126,7 @@ class TensorStorage<T, DSizes<IndexType, NumIndices_>, Options_>
}
else
m_data = 0;
EIGEN_INTERNAL_DENSE_STORAGE_CTOR_PLUGIN
EIGEN_INTERNAL_DENSE_STORAGE_CTOR_PLUGIN({})
}
m_dimensions = nbDimensions;
}

View File

@@ -61,10 +61,11 @@ struct MatrixExponentialScalingOp
* After exit, \f$ (V+U)(V-U)^{-1} \f$ is the Pad&eacute;
* approximant of \f$ \exp(A) \f$ around \f$ A = 0 \f$.
*/
template <typename MatrixType>
void matrix_exp_pade3(const MatrixType &A, MatrixType &U, MatrixType &V)
template <typename MatA, typename MatU, typename MatV>
void matrix_exp_pade3(const MatA& A, MatU& U, MatV& V)
{
typedef typename NumTraits<typename traits<MatrixType>::Scalar>::Real RealScalar;
typedef typename MatA::PlainObject MatrixType;
typedef typename NumTraits<typename traits<MatA>::Scalar>::Real RealScalar;
const RealScalar b[] = {120.L, 60.L, 12.L, 1.L};
const MatrixType A2 = A * A;
const MatrixType tmp = b[3] * A2 + b[1] * MatrixType::Identity(A.rows(), A.cols());
@@ -77,9 +78,10 @@ void matrix_exp_pade3(const MatrixType &A, MatrixType &U, MatrixType &V)
* After exit, \f$ (V+U)(V-U)^{-1} \f$ is the Pad&eacute;
* approximant of \f$ \exp(A) \f$ around \f$ A = 0 \f$.
*/
template <typename MatrixType>
void matrix_exp_pade5(const MatrixType &A, MatrixType &U, MatrixType &V)
template <typename MatA, typename MatU, typename MatV>
void matrix_exp_pade5(const MatA& A, MatU& U, MatV& V)
{
typedef typename MatA::PlainObject MatrixType;
typedef typename NumTraits<typename traits<MatrixType>::Scalar>::Real RealScalar;
const RealScalar b[] = {30240.L, 15120.L, 3360.L, 420.L, 30.L, 1.L};
const MatrixType A2 = A * A;
@@ -94,9 +96,10 @@ void matrix_exp_pade5(const MatrixType &A, MatrixType &U, MatrixType &V)
* After exit, \f$ (V+U)(V-U)^{-1} \f$ is the Pad&eacute;
* approximant of \f$ \exp(A) \f$ around \f$ A = 0 \f$.
*/
template <typename MatrixType>
void matrix_exp_pade7(const MatrixType &A, MatrixType &U, MatrixType &V)
template <typename MatA, typename MatU, typename MatV>
void matrix_exp_pade7(const MatA& A, MatU& U, MatV& V)
{
typedef typename MatA::PlainObject MatrixType;
typedef typename NumTraits<typename traits<MatrixType>::Scalar>::Real RealScalar;
const RealScalar b[] = {17297280.L, 8648640.L, 1995840.L, 277200.L, 25200.L, 1512.L, 56.L, 1.L};
const MatrixType A2 = A * A;
@@ -114,9 +117,10 @@ void matrix_exp_pade7(const MatrixType &A, MatrixType &U, MatrixType &V)
* After exit, \f$ (V+U)(V-U)^{-1} \f$ is the Pad&eacute;
* approximant of \f$ \exp(A) \f$ around \f$ A = 0 \f$.
*/
template <typename MatrixType>
void matrix_exp_pade9(const MatrixType &A, MatrixType &U, MatrixType &V)
template <typename MatA, typename MatU, typename MatV>
void matrix_exp_pade9(const MatA& A, MatU& U, MatV& V)
{
typedef typename MatA::PlainObject MatrixType;
typedef typename NumTraits<typename traits<MatrixType>::Scalar>::Real RealScalar;
const RealScalar b[] = {17643225600.L, 8821612800.L, 2075673600.L, 302702400.L, 30270240.L,
2162160.L, 110880.L, 3960.L, 90.L, 1.L};
@@ -135,9 +139,10 @@ void matrix_exp_pade9(const MatrixType &A, MatrixType &U, MatrixType &V)
* After exit, \f$ (V+U)(V-U)^{-1} \f$ is the Pad&eacute;
* approximant of \f$ \exp(A) \f$ around \f$ A = 0 \f$.
*/
template <typename MatrixType>
void matrix_exp_pade13(const MatrixType &A, MatrixType &U, MatrixType &V)
template <typename MatA, typename MatU, typename MatV>
void matrix_exp_pade13(const MatA& A, MatU& U, MatV& V)
{
typedef typename MatA::PlainObject MatrixType;
typedef typename NumTraits<typename traits<MatrixType>::Scalar>::Real RealScalar;
const RealScalar b[] = {64764752532480000.L, 32382376266240000.L, 7771770303897600.L,
1187353796428800.L, 129060195264000.L, 10559470521600.L, 670442572800.L,
@@ -162,9 +167,10 @@ void matrix_exp_pade13(const MatrixType &A, MatrixType &U, MatrixType &V)
* This function activates only if your long double is double-double or quadruple.
*/
#if LDBL_MANT_DIG > 64
template <typename MatrixType>
void matrix_exp_pade17(const MatrixType &A, MatrixType &U, MatrixType &V)
template <typename MatA, typename MatU, typename MatV>
void matrix_exp_pade17(const MatA& A, MatU& U, MatV& V)
{
typedef typename MatA::PlainObject MatrixType;
typedef typename NumTraits<typename traits<MatrixType>::Scalar>::Real RealScalar;
const RealScalar b[] = {830034394580628357120000.L, 415017197290314178560000.L,
100610229646136770560000.L, 15720348382208870400000.L,
@@ -342,9 +348,10 @@ struct matrix_exp_computeUV<MatrixType, long double>
* \param arg argument of matrix exponential (should be plain object)
* \param result variable in which result will be stored
*/
template <typename MatrixType, typename ResultType>
void matrix_exp_compute(const MatrixType& arg, ResultType &result)
template <typename ArgType, typename ResultType>
void matrix_exp_compute(const ArgType& arg, ResultType &result)
{
typedef typename ArgType::PlainObject MatrixType;
#if LDBL_MANT_DIG > 112 // rarely happens
typedef typename traits<MatrixType>::Scalar Scalar;
typedef typename NumTraits<Scalar>::Real RealScalar;
@@ -354,11 +361,11 @@ void matrix_exp_compute(const MatrixType& arg, ResultType &result)
return;
}
#endif
typename MatrixType::PlainObject U, V;
MatrixType U, V;
int squarings;
matrix_exp_computeUV<MatrixType>::run(arg, U, V, squarings); // Pade approximant is (U+V) / (-U+V)
typename MatrixType::PlainObject numer = U + V;
typename MatrixType::PlainObject denom = -U + V;
MatrixType numer = U + V;
MatrixType denom = -U + V;
result = denom.partialPivLu().solve(numer);
for (int i=0; i<squarings; i++)
result *= result; // undo scaling by repeated squaring