GitHub-Proxy/eigen
1
0
Fork 0
mirror of https://gitlab.com/libeigen/eigen.git synced 2025-08-08 17:59:00 +08:00
Code Issues Packages Projects Releases Wiki Activity

some cleaning and add copyrights

Browse Source
This commit is contained in:
Gael Guennebaud 2012-02-10 19:38:31 +01:00
parent 16da7299dd
commit ef7f1371b2
2 changed files with 16 additions and 15 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

12
Eigen/src/IterativeLinearSolvers/BiCGSTAB.h
View File

@ -2,6 +2,7 @@
// for linear algebra. // for linear algebra.
// //
// Copyright (C) 2011 Gael Guennebaud <gael.guennebaud@inria.fr> // Copyright (C) 2011 Gael Guennebaud <gael.guennebaud@inria.fr>
// Copyright (C) 2012 Désiré Nuentsa-Wakam <desire.nuentsa_wakam@inria.fr>
// //
// Eigen is free software; you can redistribute it and/or // Eigen is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public // modify it under the terms of the GNU Lesser General Public
@ -70,7 +71,6 @@ bool bicgstab(const MatrixType& mat, const Rhs& rhs, Dest& x,
while ( r.squaredNorm()/r0_sqnorm > tol2 && i<maxIters ) while ( r.squaredNorm()/r0_sqnorm > tol2 && i<maxIters )
{ {
// std::cout<<i<<" : Relative residual norm " << sqrt(r.squaredNorm()/r0_sqnorm)<<"\n";
Scalar rho_old = rho; Scalar rho_old = rho;
rho = r0.dot(r); rho = r0.dot(r);
@ -216,17 +216,19 @@ public:
template<typename Rhs,typename Dest> template<typename Rhs,typename Dest>
void _solveWithGuess(const Rhs& b, Dest& x) const void _solveWithGuess(const Rhs& b, Dest& x) const
{ {
bool ok; bool failed = false;
for(int j=0; j<b.cols(); ++j) for(int j=0; j<b.cols(); ++j)
{ {
m_iterations = Base::maxIterations(); m_iterations = Base::maxIterations();
m_error = Base::m_tolerance; m_error = Base::m_tolerance;
typename Dest::ColXpr xj(x,j); typename Dest::ColXpr xj(x,j);
ok = internal::bicgstab(*mp_matrix, b.col(j), xj, Base::m_preconditioner, m_iterations, m_error); if(!internal::bicgstab(*mp_matrix, b.col(j), xj, Base::m_preconditioner, m_iterations, m_error))
failed = true;
} }
if (ok == false) m_info = NumericalIssue; m_info = failed ? NumericalIssue
else m_info = m_error <= Base::m_tolerance ? Success : NoConvergence; : m_error <= Base::m_tolerance ? Success
: NoConvergence;
m_isInitialized = true; m_isInitialized = true;
} }

17
Eigen/src/IterativeLinearSolvers/IncompleteLUT.h
View File

@ -1,7 +1,7 @@
// This file is part of Eigen, a lightweight C++ template library // This file is part of Eigen, a lightweight C++ template library
// for linear algebra. // for linear algebra.
// //
// Copyright (C) 2011 Gael Guennebaud <gael.guennebaud@inria.fr> // Copyright (C) 2012 Désiré Nuentsa-Wakam <desire.nuentsa_wakam@inria.fr>
// //
// Eigen is free software; you can redistribute it and/or // Eigen is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public // modify it under the terms of the GNU Lesser General Public
@ -63,11 +63,11 @@ class IncompleteLUT
public: public:
typedef Matrix<Scalar,Dynamic,Dynamic> MatrixType; typedef Matrix<Scalar,Dynamic,Dynamic> MatrixType;
IncompleteLUT() : m_droptol(NumTraits<Scalar>::dummy_precision()),m_fillfactor(10),m_isInitialized(false),m_analysisIsOk(false),m_factorizationIsOk(false) {}; IncompleteLUT() : m_droptol(NumTraits<Scalar>::dummy_precision()),m_fillfactor(10),m_analysisIsOk(false),m_factorizationIsOk(false),m_isInitialized(false) {}
template<typename MatrixType> template<typename MatrixType>
IncompleteLUT(const MatrixType& mat, RealScalar droptol, int fillfactor) IncompleteLUT(const MatrixType& mat, RealScalar droptol, int fillfactor)
: m_droptol(droptol),m_fillfactor(fillfactor),m_isInitialized(false),m_analysisIsOk(false),m_factorizationIsOk(false) : m_droptol(droptol),m_fillfactor(fillfactor),m_analysisIsOk(false),m_factorizationIsOk(false),m_isInitialized(false)
{ {
eigen_assert(fillfactor != 0); eigen_assert(fillfactor != 0);
compute(mat); compute(mat);
@ -105,7 +105,7 @@ class IncompleteLUT
Vector u(n) ; /* real values of the row -- maximum size is n -- */ Vector u(n) ; /* real values of the row -- maximum size is n -- */
VectorXi ju(n); /*column position of the values in u -- maximum size is n*/ VectorXi ju(n); /*column position of the values in u -- maximum size is n*/
VectorXi jr(n); /* Indicate the position of the nonzero elements in the vector u -- A zero location is indicated by -1*/ VectorXi jr(n); /* Indicate the position of the nonzero elements in the vector u -- A zero location is indicated by -1*/
int j, k, ii, jj, jpos, minrow, len; int j, k, jj, jpos, minrow, len;
Scalar fact, prod; Scalar fact, prod;
RealScalar rownorm; RealScalar rownorm;
@ -320,8 +320,8 @@ protected:
FactorType m_lu; FactorType m_lu;
RealScalar m_droptol; RealScalar m_droptol;
int m_fillfactor; int m_fillfactor;
bool m_factorizationIsOk;
bool m_analysisIsOk; bool m_analysisIsOk;
bool m_factorizationIsOk;
bool m_isInitialized; bool m_isInitialized;
template <typename VectorV, typename VectorI> template <typename VectorV, typename VectorI>
int QuickSplit(VectorV &row, VectorI &ind, int ncut); int QuickSplit(VectorV &row, VectorI &ind, int ncut);
@ -369,10 +369,9 @@ void IncompleteLUT<Scalar>::setFillfactor(int fillfactor)
**/ **/
template <typename Scalar> template <typename Scalar>
template <typename VectorV, typename VectorI> template <typename VectorV, typename VectorI>
int IncompleteLUT<Scalar>::QuickSplit(VectorV &row, VectorI &ind, int ncut) int IncompleteLUT<Scalar>::QuickSplit(VectorV &row, VectorI &ind, int ncut)
{ {
int i,j,mid; int mid;
Scalar d;
int n = row.size(); /* lenght of the vector */ int n = row.size(); /* lenght of the vector */
int first, last ; int first, last ;
@ -384,7 +383,7 @@ int IncompleteLUT<Scalar>::QuickSplit(VectorV &row, VectorI &ind, int ncut)
do { do {
mid = first; mid = first;
RealScalar abskey = std::abs(row(mid)); RealScalar abskey = std::abs(row(mid));
for (j = first + 1; j <= last; j++) { for (int j = first + 1; j <= last; j++) {
if ( std::abs(row(j)) > abskey) { if ( std::abs(row(j)) > abskey) {
++mid; ++mid;
std::swap(row(mid), row(j)); std::swap(row(mid), row(j));