mirror of
https://gitlab.com/libeigen/eigen.git
synced 2025-08-15 13:15:57 +08:00
bug #1330: Cholmod supports double precision only, so let's trigger a static assertion if the scalar type does not match this requirement.
This commit is contained in:
parent
03fd417f66
commit
18038fc829
@ -14,34 +14,40 @@ namespace Eigen {
|
|||||||
|
|
||||||
namespace internal {
|
namespace internal {
|
||||||
|
|
||||||
template<typename Scalar, typename CholmodType>
|
template<typename Scalar> struct cholmod_configure_matrix;
|
||||||
void cholmod_configure_matrix(CholmodType& mat)
|
|
||||||
{
|
template<> struct cholmod_configure_matrix<double> {
|
||||||
if (internal::is_same<Scalar,float>::value)
|
template<typename CholmodType>
|
||||||
{
|
static void run(CholmodType& mat) {
|
||||||
mat.xtype = CHOLMOD_REAL;
|
|
||||||
mat.dtype = CHOLMOD_SINGLE;
|
|
||||||
}
|
|
||||||
else if (internal::is_same<Scalar,double>::value)
|
|
||||||
{
|
|
||||||
mat.xtype = CHOLMOD_REAL;
|
mat.xtype = CHOLMOD_REAL;
|
||||||
mat.dtype = CHOLMOD_DOUBLE;
|
mat.dtype = CHOLMOD_DOUBLE;
|
||||||
}
|
}
|
||||||
else if (internal::is_same<Scalar,std::complex<float> >::value)
|
};
|
||||||
{
|
|
||||||
mat.xtype = CHOLMOD_COMPLEX;
|
template<> struct cholmod_configure_matrix<std::complex<double> > {
|
||||||
mat.dtype = CHOLMOD_SINGLE;
|
template<typename CholmodType>
|
||||||
}
|
static void run(CholmodType& mat) {
|
||||||
else if (internal::is_same<Scalar,std::complex<double> >::value)
|
|
||||||
{
|
|
||||||
mat.xtype = CHOLMOD_COMPLEX;
|
mat.xtype = CHOLMOD_COMPLEX;
|
||||||
mat.dtype = CHOLMOD_DOUBLE;
|
mat.dtype = CHOLMOD_DOUBLE;
|
||||||
}
|
}
|
||||||
else
|
};
|
||||||
{
|
|
||||||
eigen_assert(false && "Scalar type not supported by CHOLMOD");
|
// Other scalar types are not yet suppotred by Cholmod
|
||||||
}
|
// template<> struct cholmod_configure_matrix<float> {
|
||||||
}
|
// template<typename CholmodType>
|
||||||
|
// static void run(CholmodType& mat) {
|
||||||
|
// mat.xtype = CHOLMOD_REAL;
|
||||||
|
// mat.dtype = CHOLMOD_SINGLE;
|
||||||
|
// }
|
||||||
|
// };
|
||||||
|
//
|
||||||
|
// template<> struct cholmod_configure_matrix<std::complex<float> > {
|
||||||
|
// template<typename CholmodType>
|
||||||
|
// static void run(CholmodType& mat) {
|
||||||
|
// mat.xtype = CHOLMOD_COMPLEX;
|
||||||
|
// mat.dtype = CHOLMOD_SINGLE;
|
||||||
|
// }
|
||||||
|
// };
|
||||||
|
|
||||||
} // namespace internal
|
} // namespace internal
|
||||||
|
|
||||||
@ -88,7 +94,7 @@ cholmod_sparse viewAsCholmod(SparseMatrix<_Scalar,_Options,_Index>& mat)
|
|||||||
}
|
}
|
||||||
|
|
||||||
// setup res.xtype
|
// setup res.xtype
|
||||||
internal::cholmod_configure_matrix<_Scalar>(res);
|
internal::cholmod_configure_matrix<_Scalar>::run(res);
|
||||||
|
|
||||||
res.stype = 0;
|
res.stype = 0;
|
||||||
|
|
||||||
@ -131,7 +137,7 @@ cholmod_dense viewAsCholmod(MatrixBase<Derived>& mat)
|
|||||||
res.x = (void*)(mat.derived().data());
|
res.x = (void*)(mat.derived().data());
|
||||||
res.z = 0;
|
res.z = 0;
|
||||||
|
|
||||||
internal::cholmod_configure_matrix<Scalar>(res);
|
internal::cholmod_configure_matrix<Scalar>::run(res);
|
||||||
|
|
||||||
return res;
|
return res;
|
||||||
}
|
}
|
||||||
@ -172,14 +178,16 @@ class CholmodBase : internal::noncopyable
|
|||||||
CholmodBase()
|
CholmodBase()
|
||||||
: m_cholmodFactor(0), m_info(Success), m_isInitialized(false)
|
: m_cholmodFactor(0), m_info(Success), m_isInitialized(false)
|
||||||
{
|
{
|
||||||
m_shiftOffset[0] = m_shiftOffset[1] = RealScalar(0.0);
|
EIGEN_STATIC_ASSERT((internal::is_same<double,RealScalar>::value), CHOLMOD_SUPPORTS_DOUBLE_PRECISION_ONLY);
|
||||||
|
m_shiftOffset[0] = m_shiftOffset[1] = 0.0;
|
||||||
cholmod_start(&m_cholmod);
|
cholmod_start(&m_cholmod);
|
||||||
}
|
}
|
||||||
|
|
||||||
CholmodBase(const MatrixType& matrix)
|
CholmodBase(const MatrixType& matrix)
|
||||||
: m_cholmodFactor(0), m_info(Success), m_isInitialized(false)
|
: m_cholmodFactor(0), m_info(Success), m_isInitialized(false)
|
||||||
{
|
{
|
||||||
m_shiftOffset[0] = m_shiftOffset[1] = RealScalar(0.0);
|
EIGEN_STATIC_ASSERT((internal::is_same<double,RealScalar>::value), CHOLMOD_SUPPORTS_DOUBLE_PRECISION_ONLY);
|
||||||
|
m_shiftOffset[0] = m_shiftOffset[1] = 0.0;
|
||||||
cholmod_start(&m_cholmod);
|
cholmod_start(&m_cholmod);
|
||||||
compute(matrix);
|
compute(matrix);
|
||||||
}
|
}
|
||||||
@ -277,7 +285,7 @@ class CholmodBase : internal::noncopyable
|
|||||||
eigen_assert(m_analysisIsOk && "You must first call analyzePattern()");
|
eigen_assert(m_analysisIsOk && "You must first call analyzePattern()");
|
||||||
cholmod_sparse A = viewAsCholmod(matrix.template selfadjointView<UpLo>());
|
cholmod_sparse A = viewAsCholmod(matrix.template selfadjointView<UpLo>());
|
||||||
cholmod_factorize_p(&A, m_shiftOffset, 0, 0, m_cholmodFactor, &m_cholmod);
|
cholmod_factorize_p(&A, m_shiftOffset, 0, 0, m_cholmodFactor, &m_cholmod);
|
||||||
|
|
||||||
// If the factorization failed, minor is the column at which it did. On success minor == n.
|
// If the factorization failed, minor is the column at which it did. On success minor == n.
|
||||||
this->m_info = (m_cholmodFactor->minor == m_cholmodFactor->n ? Success : NumericalIssue);
|
this->m_info = (m_cholmodFactor->minor == m_cholmodFactor->n ? Success : NumericalIssue);
|
||||||
m_factorizationIsOk = true;
|
m_factorizationIsOk = true;
|
||||||
@ -344,7 +352,7 @@ class CholmodBase : internal::noncopyable
|
|||||||
*/
|
*/
|
||||||
Derived& setShift(const RealScalar& offset)
|
Derived& setShift(const RealScalar& offset)
|
||||||
{
|
{
|
||||||
m_shiftOffset[0] = offset;
|
m_shiftOffset[0] = double(offset);
|
||||||
return derived();
|
return derived();
|
||||||
}
|
}
|
||||||
|
|
||||||
@ -355,7 +363,7 @@ class CholmodBase : internal::noncopyable
|
|||||||
protected:
|
protected:
|
||||||
mutable cholmod_common m_cholmod;
|
mutable cholmod_common m_cholmod;
|
||||||
cholmod_factor* m_cholmodFactor;
|
cholmod_factor* m_cholmodFactor;
|
||||||
RealScalar m_shiftOffset[2];
|
double m_shiftOffset[2];
|
||||||
mutable ComputationInfo m_info;
|
mutable ComputationInfo m_info;
|
||||||
bool m_isInitialized;
|
bool m_isInitialized;
|
||||||
int m_factorizationIsOk;
|
int m_factorizationIsOk;
|
||||||
@ -378,6 +386,8 @@ class CholmodBase : internal::noncopyable
|
|||||||
*
|
*
|
||||||
* This class supports all kind of SparseMatrix<>: row or column major; upper, lower, or both; compressed or non compressed.
|
* This class supports all kind of SparseMatrix<>: row or column major; upper, lower, or both; compressed or non compressed.
|
||||||
*
|
*
|
||||||
|
* \warning Only double precision real and complex scalar types are supported by Cholmod.
|
||||||
|
*
|
||||||
* \sa \ref TutorialSparseDirectSolvers, class CholmodSupernodalLLT, class SimplicialLLT
|
* \sa \ref TutorialSparseDirectSolvers, class CholmodSupernodalLLT, class SimplicialLLT
|
||||||
*/
|
*/
|
||||||
template<typename _MatrixType, int _UpLo = Lower>
|
template<typename _MatrixType, int _UpLo = Lower>
|
||||||
@ -425,6 +435,8 @@ class CholmodSimplicialLLT : public CholmodBase<_MatrixType, _UpLo, CholmodSimpl
|
|||||||
*
|
*
|
||||||
* This class supports all kind of SparseMatrix<>: row or column major; upper, lower, or both; compressed or non compressed.
|
* This class supports all kind of SparseMatrix<>: row or column major; upper, lower, or both; compressed or non compressed.
|
||||||
*
|
*
|
||||||
|
* \warning Only double precision real and complex scalar types are supported by Cholmod.
|
||||||
|
*
|
||||||
* \sa \ref TutorialSparseDirectSolvers, class CholmodSupernodalLLT, class SimplicialLDLT
|
* \sa \ref TutorialSparseDirectSolvers, class CholmodSupernodalLLT, class SimplicialLDLT
|
||||||
*/
|
*/
|
||||||
template<typename _MatrixType, int _UpLo = Lower>
|
template<typename _MatrixType, int _UpLo = Lower>
|
||||||
@ -470,6 +482,8 @@ class CholmodSimplicialLDLT : public CholmodBase<_MatrixType, _UpLo, CholmodSimp
|
|||||||
*
|
*
|
||||||
* This class supports all kind of SparseMatrix<>: row or column major; upper, lower, or both; compressed or non compressed.
|
* This class supports all kind of SparseMatrix<>: row or column major; upper, lower, or both; compressed or non compressed.
|
||||||
*
|
*
|
||||||
|
* \warning Only double precision real and complex scalar types are supported by Cholmod.
|
||||||
|
*
|
||||||
* \sa \ref TutorialSparseDirectSolvers
|
* \sa \ref TutorialSparseDirectSolvers
|
||||||
*/
|
*/
|
||||||
template<typename _MatrixType, int _UpLo = Lower>
|
template<typename _MatrixType, int _UpLo = Lower>
|
||||||
@ -517,6 +531,8 @@ class CholmodSupernodalLLT : public CholmodBase<_MatrixType, _UpLo, CholmodSuper
|
|||||||
*
|
*
|
||||||
* This class supports all kind of SparseMatrix<>: row or column major; upper, lower, or both; compressed or non compressed.
|
* This class supports all kind of SparseMatrix<>: row or column major; upper, lower, or both; compressed or non compressed.
|
||||||
*
|
*
|
||||||
|
* \warning Only double precision real and complex scalar types are supported by Cholmod.
|
||||||
|
*
|
||||||
* \sa \ref TutorialSparseDirectSolvers
|
* \sa \ref TutorialSparseDirectSolvers
|
||||||
*/
|
*/
|
||||||
template<typename _MatrixType, int _UpLo = Lower>
|
template<typename _MatrixType, int _UpLo = Lower>
|
||||||
|
@ -92,7 +92,8 @@
|
|||||||
THE_STORAGE_ORDER_OF_BOTH_SIDES_MUST_MATCH,
|
THE_STORAGE_ORDER_OF_BOTH_SIDES_MUST_MATCH,
|
||||||
OBJECT_ALLOCATED_ON_STACK_IS_TOO_BIG,
|
OBJECT_ALLOCATED_ON_STACK_IS_TOO_BIG,
|
||||||
IMPLICIT_CONVERSION_TO_SCALAR_IS_FOR_INNER_PRODUCT_ONLY,
|
IMPLICIT_CONVERSION_TO_SCALAR_IS_FOR_INNER_PRODUCT_ONLY,
|
||||||
STORAGE_LAYOUT_DOES_NOT_MATCH
|
STORAGE_LAYOUT_DOES_NOT_MATCH,
|
||||||
|
CHOLMOD_SUPPORTS_DOUBLE_PRECISION_ONLY
|
||||||
};
|
};
|
||||||
};
|
};
|
||||||
|
|
||||||
|
Loading…
x
Reference in New Issue
Block a user