Back to home page

EIC code displayed by LXR
 
 

    

Warning, file /include/eigen3/Eigen/src/PardisoSupport/PardisoSupport.h was not indexed or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).

0001 /*
0002  Copyright (c) 2011, Intel Corporation. All rights reserved.
0003 
0004  Redistribution and use in source and binary forms, with or without modification,
0005  are permitted provided that the following conditions are met:
0006 
0007  * Redistributions of source code must retain the above copyright notice, this
0008    list of conditions and the following disclaimer.
0009  * Redistributions in binary form must reproduce the above copyright notice,
0010    this list of conditions and the following disclaimer in the documentation
0011    and/or other materials provided with the distribution.
0012  * Neither the name of Intel Corporation nor the names of its contributors may
0013    be used to endorse or promote products derived from this software without
0014    specific prior written permission.
0015 
0016  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
0017  ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
0018  WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
0019  DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
0020  ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
0021  (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
0022  LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
0023  ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
0024  (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
0025  SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
0026 
0027  ********************************************************************************
0028  *   Content : Eigen bindings to Intel(R) MKL PARDISO
0029  ********************************************************************************
0030 */
0031 
0032 #ifndef EIGEN_PARDISOSUPPORT_H
0033 #define EIGEN_PARDISOSUPPORT_H
0034 
0035 namespace Eigen { 
0036 
0037 template<typename _MatrixType> class PardisoLU;
0038 template<typename _MatrixType, int Options=Upper> class PardisoLLT;
0039 template<typename _MatrixType, int Options=Upper> class PardisoLDLT;
0040 
0041 namespace internal
0042 {
0043   template<typename IndexType>
0044   struct pardiso_run_selector
0045   {
0046     static IndexType run( _MKL_DSS_HANDLE_t pt, IndexType maxfct, IndexType mnum, IndexType type, IndexType phase, IndexType n, void *a,
0047                       IndexType *ia, IndexType *ja, IndexType *perm, IndexType nrhs, IndexType *iparm, IndexType msglvl, void *b, void *x)
0048     {
0049       IndexType error = 0;
0050       ::pardiso(pt, &maxfct, &mnum, &type, &phase, &n, a, ia, ja, perm, &nrhs, iparm, &msglvl, b, x, &error);
0051       return error;
0052     }
0053   };
0054   template<>
0055   struct pardiso_run_selector<long long int>
0056   {
0057     typedef long long int IndexType;
0058     static IndexType run( _MKL_DSS_HANDLE_t pt, IndexType maxfct, IndexType mnum, IndexType type, IndexType phase, IndexType n, void *a,
0059                       IndexType *ia, IndexType *ja, IndexType *perm, IndexType nrhs, IndexType *iparm, IndexType msglvl, void *b, void *x)
0060     {
0061       IndexType error = 0;
0062       ::pardiso_64(pt, &maxfct, &mnum, &type, &phase, &n, a, ia, ja, perm, &nrhs, iparm, &msglvl, b, x, &error);
0063       return error;
0064     }
0065   };
0066 
0067   template<class Pardiso> struct pardiso_traits;
0068 
0069   template<typename _MatrixType>
0070   struct pardiso_traits< PardisoLU<_MatrixType> >
0071   {
0072     typedef _MatrixType MatrixType;
0073     typedef typename _MatrixType::Scalar Scalar;
0074     typedef typename _MatrixType::RealScalar RealScalar;
0075     typedef typename _MatrixType::StorageIndex StorageIndex;
0076   };
0077 
0078   template<typename _MatrixType, int Options>
0079   struct pardiso_traits< PardisoLLT<_MatrixType, Options> >
0080   {
0081     typedef _MatrixType MatrixType;
0082     typedef typename _MatrixType::Scalar Scalar;
0083     typedef typename _MatrixType::RealScalar RealScalar;
0084     typedef typename _MatrixType::StorageIndex StorageIndex;
0085   };
0086 
0087   template<typename _MatrixType, int Options>
0088   struct pardiso_traits< PardisoLDLT<_MatrixType, Options> >
0089   {
0090     typedef _MatrixType MatrixType;
0091     typedef typename _MatrixType::Scalar Scalar;
0092     typedef typename _MatrixType::RealScalar RealScalar;
0093     typedef typename _MatrixType::StorageIndex StorageIndex;    
0094   };
0095 
0096 } // end namespace internal
0097 
0098 template<class Derived>
0099 class PardisoImpl : public SparseSolverBase<Derived>
0100 {
0101   protected:
0102     typedef SparseSolverBase<Derived> Base;
0103     using Base::derived;
0104     using Base::m_isInitialized;
0105     
0106     typedef internal::pardiso_traits<Derived> Traits;
0107   public:
0108     using Base::_solve_impl;
0109     
0110     typedef typename Traits::MatrixType MatrixType;
0111     typedef typename Traits::Scalar Scalar;
0112     typedef typename Traits::RealScalar RealScalar;
0113     typedef typename Traits::StorageIndex StorageIndex;
0114     typedef SparseMatrix<Scalar,RowMajor,StorageIndex> SparseMatrixType;
0115     typedef Matrix<Scalar,Dynamic,1> VectorType;
0116     typedef Matrix<StorageIndex, 1, MatrixType::ColsAtCompileTime> IntRowVectorType;
0117     typedef Matrix<StorageIndex, MatrixType::RowsAtCompileTime, 1> IntColVectorType;
0118     typedef Array<StorageIndex,64,1,DontAlign> ParameterType;
0119     enum {
0120       ScalarIsComplex = NumTraits<Scalar>::IsComplex,
0121       ColsAtCompileTime = Dynamic,
0122       MaxColsAtCompileTime = Dynamic
0123     };
0124 
0125     PardisoImpl()
0126       : m_analysisIsOk(false), m_factorizationIsOk(false)
0127     {
0128       eigen_assert((sizeof(StorageIndex) >= sizeof(_INTEGER_t) && sizeof(StorageIndex) <= 8) && "Non-supported index type");
0129       m_iparm.setZero();
0130       m_msglvl = 0; // No output
0131       m_isInitialized = false;
0132     }
0133 
0134     ~PardisoImpl()
0135     {
0136       pardisoRelease();
0137     }
0138 
0139     inline Index cols() const { return m_size; }
0140     inline Index rows() const { return m_size; }
0141   
0142     /** \brief Reports whether previous computation was successful.
0143       *
0144       * \returns \c Success if computation was successful,
0145       *          \c NumericalIssue if the matrix appears to be negative.
0146       */
0147     ComputationInfo info() const
0148     {
0149       eigen_assert(m_isInitialized && "Decomposition is not initialized.");
0150       return m_info;
0151     }
0152 
0153     /** \warning for advanced usage only.
0154       * \returns a reference to the parameter array controlling PARDISO.
0155       * See the PARDISO manual to know how to use it. */
0156     ParameterType& pardisoParameterArray()
0157     {
0158       return m_iparm;
0159     }
0160     
0161     /** Performs a symbolic decomposition on the sparcity of \a matrix.
0162       *
0163       * This function is particularly useful when solving for several problems having the same structure.
0164       * 
0165       * \sa factorize()
0166       */
0167     Derived& analyzePattern(const MatrixType& matrix);
0168     
0169     /** Performs a numeric decomposition of \a matrix
0170       *
0171       * The given matrix must has the same sparcity than the matrix on which the symbolic decomposition has been performed.
0172       *
0173       * \sa analyzePattern()
0174       */
0175     Derived& factorize(const MatrixType& matrix);
0176 
0177     Derived& compute(const MatrixType& matrix);
0178 
0179     template<typename Rhs,typename Dest>
0180     void _solve_impl(const MatrixBase<Rhs> &b, MatrixBase<Dest> &dest) const;
0181 
0182   protected:
0183     void pardisoRelease()
0184     {
0185       if(m_isInitialized) // Factorization ran at least once
0186       {
0187         internal::pardiso_run_selector<StorageIndex>::run(m_pt, 1, 1, m_type, -1, internal::convert_index<StorageIndex>(m_size),0, 0, 0, m_perm.data(), 0,
0188                                                           m_iparm.data(), m_msglvl, NULL, NULL);
0189         m_isInitialized = false;
0190       }
0191     }
0192 
0193     void pardisoInit(int type)
0194     {
0195       m_type = type;
0196       bool symmetric = std::abs(m_type) < 10;
0197       m_iparm[0] = 1;   // No solver default
0198       m_iparm[1] = 2;   // use Metis for the ordering
0199       m_iparm[2] = 0;   // Reserved. Set to zero. (??Numbers of processors, value of OMP_NUM_THREADS??)
0200       m_iparm[3] = 0;   // No iterative-direct algorithm
0201       m_iparm[4] = 0;   // No user fill-in reducing permutation
0202       m_iparm[5] = 0;   // Write solution into x, b is left unchanged
0203       m_iparm[6] = 0;   // Not in use
0204       m_iparm[7] = 2;   // Max numbers of iterative refinement steps
0205       m_iparm[8] = 0;   // Not in use
0206       m_iparm[9] = 13;  // Perturb the pivot elements with 1E-13
0207       m_iparm[10] = symmetric ? 0 : 1; // Use nonsymmetric permutation and scaling MPS
0208       m_iparm[11] = 0;  // Not in use
0209       m_iparm[12] = symmetric ? 0 : 1;  // Maximum weighted matching algorithm is switched-off (default for symmetric).
0210                                         // Try m_iparm[12] = 1 in case of inappropriate accuracy
0211       m_iparm[13] = 0;  // Output: Number of perturbed pivots
0212       m_iparm[14] = 0;  // Not in use
0213       m_iparm[15] = 0;  // Not in use
0214       m_iparm[16] = 0;  // Not in use
0215       m_iparm[17] = -1; // Output: Number of nonzeros in the factor LU
0216       m_iparm[18] = -1; // Output: Mflops for LU factorization
0217       m_iparm[19] = 0;  // Output: Numbers of CG Iterations
0218       
0219       m_iparm[20] = 0;  // 1x1 pivoting
0220       m_iparm[26] = 0;  // No matrix checker
0221       m_iparm[27] = (sizeof(RealScalar) == 4) ? 1 : 0;
0222       m_iparm[34] = 1;  // C indexing
0223       m_iparm[36] = 0;  // CSR
0224       m_iparm[59] = 0;  // 0 - In-Core ; 1 - Automatic switch between In-Core and Out-of-Core modes ; 2 - Out-of-Core
0225       
0226       memset(m_pt, 0, sizeof(m_pt));
0227     }
0228 
0229   protected:
0230     // cached data to reduce reallocation, etc.
0231     
0232     void manageErrorCode(Index error) const
0233     {
0234       switch(error)
0235       {
0236         case 0:
0237           m_info = Success;
0238           break;
0239         case -4:
0240         case -7:
0241           m_info = NumericalIssue;
0242           break;
0243         default:
0244           m_info = InvalidInput;
0245       }
0246     }
0247 
0248     mutable SparseMatrixType m_matrix;
0249     mutable ComputationInfo m_info;
0250     bool m_analysisIsOk, m_factorizationIsOk;
0251     StorageIndex m_type, m_msglvl;
0252     mutable void *m_pt[64];
0253     mutable ParameterType m_iparm;
0254     mutable IntColVectorType m_perm;
0255     Index m_size;
0256     
0257 };
0258 
0259 template<class Derived>
0260 Derived& PardisoImpl<Derived>::compute(const MatrixType& a)
0261 {
0262   m_size = a.rows();
0263   eigen_assert(a.rows() == a.cols());
0264 
0265   pardisoRelease();
0266   m_perm.setZero(m_size);
0267   derived().getMatrix(a);
0268   
0269   Index error;
0270   error = internal::pardiso_run_selector<StorageIndex>::run(m_pt, 1, 1, m_type, 12, internal::convert_index<StorageIndex>(m_size),
0271                                                             m_matrix.valuePtr(), m_matrix.outerIndexPtr(), m_matrix.innerIndexPtr(),
0272                                                             m_perm.data(), 0, m_iparm.data(), m_msglvl, NULL, NULL);
0273   manageErrorCode(error);
0274   m_analysisIsOk = true;
0275   m_factorizationIsOk = true;
0276   m_isInitialized = true;
0277   return derived();
0278 }
0279 
0280 template<class Derived>
0281 Derived& PardisoImpl<Derived>::analyzePattern(const MatrixType& a)
0282 {
0283   m_size = a.rows();
0284   eigen_assert(m_size == a.cols());
0285 
0286   pardisoRelease();
0287   m_perm.setZero(m_size);
0288   derived().getMatrix(a);
0289   
0290   Index error;
0291   error = internal::pardiso_run_selector<StorageIndex>::run(m_pt, 1, 1, m_type, 11, internal::convert_index<StorageIndex>(m_size),
0292                                                             m_matrix.valuePtr(), m_matrix.outerIndexPtr(), m_matrix.innerIndexPtr(),
0293                                                             m_perm.data(), 0, m_iparm.data(), m_msglvl, NULL, NULL);
0294   
0295   manageErrorCode(error);
0296   m_analysisIsOk = true;
0297   m_factorizationIsOk = false;
0298   m_isInitialized = true;
0299   return derived();
0300 }
0301 
0302 template<class Derived>
0303 Derived& PardisoImpl<Derived>::factorize(const MatrixType& a)
0304 {
0305   eigen_assert(m_analysisIsOk && "You must first call analyzePattern()");
0306   eigen_assert(m_size == a.rows() && m_size == a.cols());
0307   
0308   derived().getMatrix(a);
0309 
0310   Index error;
0311   error = internal::pardiso_run_selector<StorageIndex>::run(m_pt, 1, 1, m_type, 22, internal::convert_index<StorageIndex>(m_size),
0312                                                             m_matrix.valuePtr(), m_matrix.outerIndexPtr(), m_matrix.innerIndexPtr(),
0313                                                             m_perm.data(), 0, m_iparm.data(), m_msglvl, NULL, NULL);
0314   
0315   manageErrorCode(error);
0316   m_factorizationIsOk = true;
0317   return derived();
0318 }
0319 
0320 template<class Derived>
0321 template<typename BDerived,typename XDerived>
0322 void PardisoImpl<Derived>::_solve_impl(const MatrixBase<BDerived> &b, MatrixBase<XDerived>& x) const
0323 {
0324   if(m_iparm[0] == 0) // Factorization was not computed
0325   {
0326     m_info = InvalidInput;
0327     return;
0328   }
0329 
0330   //Index n = m_matrix.rows();
0331   Index nrhs = Index(b.cols());
0332   eigen_assert(m_size==b.rows());
0333   eigen_assert(((MatrixBase<BDerived>::Flags & RowMajorBit) == 0 || nrhs == 1) && "Row-major right hand sides are not supported");
0334   eigen_assert(((MatrixBase<XDerived>::Flags & RowMajorBit) == 0 || nrhs == 1) && "Row-major matrices of unknowns are not supported");
0335   eigen_assert(((nrhs == 1) || b.outerStride() == b.rows()));
0336 
0337 
0338 //  switch (transposed) {
0339 //    case SvNoTrans    : m_iparm[11] = 0 ; break;
0340 //    case SvTranspose  : m_iparm[11] = 2 ; break;
0341 //    case SvAdjoint    : m_iparm[11] = 1 ; break;
0342 //    default:
0343 //      //std::cerr << "Eigen: transposition  option \"" << transposed << "\" not supported by the PARDISO backend\n";
0344 //      m_iparm[11] = 0;
0345 //  }
0346 
0347   Scalar* rhs_ptr = const_cast<Scalar*>(b.derived().data());
0348   Matrix<Scalar,Dynamic,Dynamic,ColMajor> tmp;
0349   
0350   // Pardiso cannot solve in-place
0351   if(rhs_ptr == x.derived().data())
0352   {
0353     tmp = b;
0354     rhs_ptr = tmp.data();
0355   }
0356   
0357   Index error;
0358   error = internal::pardiso_run_selector<StorageIndex>::run(m_pt, 1, 1, m_type, 33, internal::convert_index<StorageIndex>(m_size),
0359                                                             m_matrix.valuePtr(), m_matrix.outerIndexPtr(), m_matrix.innerIndexPtr(),
0360                                                             m_perm.data(), internal::convert_index<StorageIndex>(nrhs), m_iparm.data(), m_msglvl,
0361                                                             rhs_ptr, x.derived().data());
0362 
0363   manageErrorCode(error);
0364 }
0365 
0366 
0367 /** \ingroup PardisoSupport_Module
0368   * \class PardisoLU
0369   * \brief A sparse direct LU factorization and solver based on the PARDISO library
0370   *
0371   * This class allows to solve for A.X = B sparse linear problems via a direct LU factorization
0372   * using the Intel MKL PARDISO library. The sparse matrix A must be squared and invertible.
0373   * The vectors or matrices X and B can be either dense or sparse.
0374   *
0375   * By default, it runs in in-core mode. To enable PARDISO's out-of-core feature, set:
0376   * \code solver.pardisoParameterArray()[59] = 1; \endcode
0377   *
0378   * \tparam _MatrixType the type of the sparse matrix A, it must be a SparseMatrix<>
0379   *
0380   * \implsparsesolverconcept
0381   *
0382   * \sa \ref TutorialSparseSolverConcept, class SparseLU
0383   */
0384 template<typename MatrixType>
0385 class PardisoLU : public PardisoImpl< PardisoLU<MatrixType> >
0386 {
0387   protected:
0388     typedef PardisoImpl<PardisoLU> Base;
0389     using Base::pardisoInit;
0390     using Base::m_matrix;
0391     friend class PardisoImpl< PardisoLU<MatrixType> >;
0392 
0393   public:
0394 
0395     typedef typename Base::Scalar Scalar;
0396     typedef typename Base::RealScalar RealScalar;
0397 
0398     using Base::compute;
0399     using Base::solve;
0400 
0401     PardisoLU()
0402       : Base()
0403     {
0404       pardisoInit(Base::ScalarIsComplex ? 13 : 11);
0405     }
0406 
0407     explicit PardisoLU(const MatrixType& matrix)
0408       : Base()
0409     {
0410       pardisoInit(Base::ScalarIsComplex ? 13 : 11);
0411       compute(matrix);
0412     }
0413   protected:
0414     void getMatrix(const MatrixType& matrix)
0415     {
0416       m_matrix = matrix;
0417       m_matrix.makeCompressed();
0418     }
0419 };
0420 
0421 /** \ingroup PardisoSupport_Module
0422   * \class PardisoLLT
0423   * \brief A sparse direct Cholesky (LLT) factorization and solver based on the PARDISO library
0424   *
0425   * This class allows to solve for A.X = B sparse linear problems via a LL^T Cholesky factorization
0426   * using the Intel MKL PARDISO library. The sparse matrix A must be selfajoint and positive definite.
0427   * The vectors or matrices X and B can be either dense or sparse.
0428   *
0429   * By default, it runs in in-core mode. To enable PARDISO's out-of-core feature, set:
0430   * \code solver.pardisoParameterArray()[59] = 1; \endcode
0431   *
0432   * \tparam MatrixType the type of the sparse matrix A, it must be a SparseMatrix<>
0433   * \tparam UpLo can be any bitwise combination of Upper, Lower. The default is Upper, meaning only the upper triangular part has to be used.
0434   *         Upper|Lower can be used to tell both triangular parts can be used as input.
0435   *
0436   * \implsparsesolverconcept
0437   *
0438   * \sa \ref TutorialSparseSolverConcept, class SimplicialLLT
0439   */
0440 template<typename MatrixType, int _UpLo>
0441 class PardisoLLT : public PardisoImpl< PardisoLLT<MatrixType,_UpLo> >
0442 {
0443   protected:
0444     typedef PardisoImpl< PardisoLLT<MatrixType,_UpLo> > Base;
0445     using Base::pardisoInit;
0446     using Base::m_matrix;
0447     friend class PardisoImpl< PardisoLLT<MatrixType,_UpLo> >;
0448 
0449   public:
0450 
0451     typedef typename Base::Scalar Scalar;
0452     typedef typename Base::RealScalar RealScalar;
0453     typedef typename Base::StorageIndex StorageIndex;
0454     enum { UpLo = _UpLo };
0455     using Base::compute;
0456 
0457     PardisoLLT()
0458       : Base()
0459     {
0460       pardisoInit(Base::ScalarIsComplex ? 4 : 2);
0461     }
0462 
0463     explicit PardisoLLT(const MatrixType& matrix)
0464       : Base()
0465     {
0466       pardisoInit(Base::ScalarIsComplex ? 4 : 2);
0467       compute(matrix);
0468     }
0469     
0470   protected:
0471     
0472     void getMatrix(const MatrixType& matrix)
0473     {
0474       // PARDISO supports only upper, row-major matrices
0475       PermutationMatrix<Dynamic,Dynamic,StorageIndex> p_null;
0476       m_matrix.resize(matrix.rows(), matrix.cols());
0477       m_matrix.template selfadjointView<Upper>() = matrix.template selfadjointView<UpLo>().twistedBy(p_null);
0478       m_matrix.makeCompressed();
0479     }
0480 };
0481 
0482 /** \ingroup PardisoSupport_Module
0483   * \class PardisoLDLT
0484   * \brief A sparse direct Cholesky (LDLT) factorization and solver based on the PARDISO library
0485   *
0486   * This class allows to solve for A.X = B sparse linear problems via a LDL^T Cholesky factorization
0487   * using the Intel MKL PARDISO library. The sparse matrix A is assumed to be selfajoint and positive definite.
0488   * For complex matrices, A can also be symmetric only, see the \a Options template parameter.
0489   * The vectors or matrices X and B can be either dense or sparse.
0490   *
0491   * By default, it runs in in-core mode. To enable PARDISO's out-of-core feature, set:
0492   * \code solver.pardisoParameterArray()[59] = 1; \endcode
0493   *
0494   * \tparam MatrixType the type of the sparse matrix A, it must be a SparseMatrix<>
0495   * \tparam Options can be any bitwise combination of Upper, Lower, and Symmetric. The default is Upper, meaning only the upper triangular part has to be used.
0496   *         Symmetric can be used for symmetric, non-selfadjoint complex matrices, the default being to assume a selfadjoint matrix.
0497   *         Upper|Lower can be used to tell both triangular parts can be used as input.
0498   *
0499   * \implsparsesolverconcept
0500   *
0501   * \sa \ref TutorialSparseSolverConcept, class SimplicialLDLT
0502   */
0503 template<typename MatrixType, int Options>
0504 class PardisoLDLT : public PardisoImpl< PardisoLDLT<MatrixType,Options> >
0505 {
0506   protected:
0507     typedef PardisoImpl< PardisoLDLT<MatrixType,Options> > Base;
0508     using Base::pardisoInit;
0509     using Base::m_matrix;
0510     friend class PardisoImpl< PardisoLDLT<MatrixType,Options> >;
0511 
0512   public:
0513 
0514     typedef typename Base::Scalar Scalar;
0515     typedef typename Base::RealScalar RealScalar;
0516     typedef typename Base::StorageIndex StorageIndex;
0517     using Base::compute;
0518     enum { UpLo = Options&(Upper|Lower) };
0519 
0520     PardisoLDLT()
0521       : Base()
0522     {
0523       pardisoInit(Base::ScalarIsComplex ? ( bool(Options&Symmetric) ? 6 : -4 ) : -2);
0524     }
0525 
0526     explicit PardisoLDLT(const MatrixType& matrix)
0527       : Base()
0528     {
0529       pardisoInit(Base::ScalarIsComplex ? ( bool(Options&Symmetric) ? 6 : -4 ) : -2);
0530       compute(matrix);
0531     }
0532     
0533     void getMatrix(const MatrixType& matrix)
0534     {
0535       // PARDISO supports only upper, row-major matrices
0536       PermutationMatrix<Dynamic,Dynamic,StorageIndex> p_null;
0537       m_matrix.resize(matrix.rows(), matrix.cols());
0538       m_matrix.template selfadjointView<Upper>() = matrix.template selfadjointView<UpLo>().twistedBy(p_null);
0539       m_matrix.makeCompressed();
0540     }
0541 };
0542 
0543 } // end namespace Eigen
0544 
0545 #endif // EIGEN_PARDISOSUPPORT_H
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!