Eigen::IterScaling< MatrixType_ > Class Template Reference

iterative scaling algorithm to equilibrate rows and column norms in matrices More...

Public Types

typedef MatrixType::Index Index
 
typedef MatrixType_ MatrixType
 
typedef MatrixType::Scalar Scalar
 

Public Member Functions

void compute (const MatrixType &mat)
 
void computeRef (MatrixType &mat)
 
 IterScaling ()
 
 IterScaling (const MatrixType &matrix)
 
VectorXdLeftScaling ()
 
VectorXdRightScaling ()
 
void setTolerance (double tol)
 
 ~IterScaling ()
 

Protected Member Functions

void init ()
 

Protected Attributes

ComputationInfo m_info
 
bool m_isInitialized
 
VectorXd m_left
 
MatrixType m_matrix
 
int m_maxits
 
VectorXd m_right
 
double m_tol
 

Detailed Description

template<typename MatrixType_>
class Eigen::IterScaling< MatrixType_ >

iterative scaling algorithm to equilibrate rows and column norms in matrices

This class can be used as a preprocessing tool to accelerate the convergence of iterative methods

This feature is useful to limit the pivoting amount during LU/ILU factorization The scaling strategy as presented here preserves the symmetry of the problem NOTE It is assumed that the matrix does not have empty row or column,

Example with key steps

VectorXd x(n), b(n);
// fill A and b;
IterScaling<SparseMatrix<double> > scal;
// Compute the left and right scaling vectors. The matrix is equilibrated at output
scal.computeRef(A);
// Scale the right hand side
b = scal.LeftScaling().cwiseProduct(b);
// Now, solve the equilibrated linear system with any available solver
// Scale back the computed solution
x = scal.RightScaling().cwiseProduct(x);
SparseMatrix< double > A(n, n)
Matrix< double, Dynamic, 1 > VectorXd
Template Parameters
MatrixType_the type of the matrix. It should be a real square sparsematrix

References : D. Ruiz and B. Ucar, A Symmetry Preserving Algorithm for Matrix Scaling, INRIA Research report RR-7552

See also
IncompleteLUT

Definition at line 50 of file Scaling.h.

Member Typedef Documentation

◆ Index

template<typename MatrixType_ >
typedef MatrixType::Index Eigen::IterScaling< MatrixType_ >::Index

Definition at line 55 of file Scaling.h.

◆ MatrixType

template<typename MatrixType_ >
typedef MatrixType_ Eigen::IterScaling< MatrixType_ >::MatrixType

Definition at line 53 of file Scaling.h.

◆ Scalar

template<typename MatrixType_ >
typedef MatrixType::Scalar Eigen::IterScaling< MatrixType_ >::Scalar

Definition at line 54 of file Scaling.h.

Constructor & Destructor Documentation

◆ IterScaling() [1/2]

template<typename MatrixType_ >
Eigen::IterScaling< MatrixType_ >::IterScaling ( )
inline

Definition at line 58 of file Scaling.h.

58 { init(); }

◆ IterScaling() [2/2]

template<typename MatrixType_ >
Eigen::IterScaling< MatrixType_ >::IterScaling ( const MatrixType matrix)
inline

Definition at line 60 of file Scaling.h.

61  {
62  init();
63  compute(matrix);
64  }
void compute(const MatrixType &mat)
Definition: Scaling.h:75

◆ ~IterScaling()

template<typename MatrixType_ >
Eigen::IterScaling< MatrixType_ >::~IterScaling ( )
inline

Definition at line 66 of file Scaling.h.

66 { }

Member Function Documentation

◆ compute()

template<typename MatrixType_ >
void Eigen::IterScaling< MatrixType_ >::compute ( const MatrixType mat)
inline

Compute the left and right diagonal matrices to scale the input matrix mat

FIXME This algorithm will be modified such that the diagonal elements are permuted on the diagonal.

See also
LeftScaling() RightScaling()

Definition at line 75 of file Scaling.h.

76  {
77  using std::abs;
78  int m = mat.rows();
79  int n = mat.cols();
80  eigen_assert((m>0 && m == n) && "Please give a non - empty matrix");
81  m_left.resize(m);
82  m_right.resize(n);
83  m_left.setOnes();
84  m_right.setOnes();
85  m_matrix = mat;
86  VectorXd Dr, Dc, DrRes, DcRes; // Temporary Left and right scaling vectors
87  Dr.resize(m); Dc.resize(n);
88  DrRes.resize(m); DcRes.resize(n);
89  double EpsRow = 1.0, EpsCol = 1.0;
90  int its = 0;
91  do
92  { // Iterate until the infinite norm of each row and column is approximately 1
93  // Get the maximum value in each row and column
94  Dr.setZero(); Dc.setZero();
95  for (int k=0; k<m_matrix.outerSize(); ++k)
96  {
97  for (typename MatrixType::InnerIterator it(m_matrix, k); it; ++it)
98  {
99  if ( Dr(it.row()) < abs(it.value()) )
100  Dr(it.row()) = abs(it.value());
101 
102  if ( Dc(it.col()) < abs(it.value()) )
103  Dc(it.col()) = abs(it.value());
104  }
105  }
106  for (int i = 0; i < m; ++i)
107  {
108  Dr(i) = std::sqrt(Dr(i));
109  }
110  for (int i = 0; i < n; ++i)
111  {
112  Dc(i) = std::sqrt(Dc(i));
113  }
114  // Save the scaling factors
115  for (int i = 0; i < m; ++i)
116  {
117  m_left(i) /= Dr(i);
118  }
119  for (int i = 0; i < n; ++i)
120  {
121  m_right(i) /= Dc(i);
122  }
123  // Scale the rows and the columns of the matrix
124  DrRes.setZero(); DcRes.setZero();
125  for (int k=0; k<m_matrix.outerSize(); ++k)
126  {
127  for (typename MatrixType::InnerIterator it(m_matrix, k); it; ++it)
128  {
129  it.valueRef() = it.value()/( Dr(it.row()) * Dc(it.col()) );
130  // Accumulate the norms of the row and column vectors
131  if ( DrRes(it.row()) < abs(it.value()) )
132  DrRes(it.row()) = abs(it.value());
133 
134  if ( DcRes(it.col()) < abs(it.value()) )
135  DcRes(it.col()) = abs(it.value());
136  }
137  }
138  DrRes.array() = (1-DrRes.array()).abs();
139  EpsRow = DrRes.maxCoeff();
140  DcRes.array() = (1-DcRes.array()).abs();
141  EpsCol = DcRes.maxCoeff();
142  its++;
143  }while ( (EpsRow >m_tol || EpsCol > m_tol) && (its < m_maxits) );
144  m_isInitialized = true;
145  }
Matrix3f m
int n
int i
#define eigen_assert(x)
MatrixXf mat
VectorXd m_left
Definition: Scaling.h:189
MatrixType m_matrix
Definition: Scaling.h:186
bool m_isInitialized
Definition: Scaling.h:188
VectorXd m_right
Definition: Scaling.h:190
Derived & setOnes(Index rows, Index cols)
constexpr void resize(Index rows, Index cols)
Eigen::AutoDiffScalar< EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(Eigen::internal::remove_all_t< DerType >, typename Eigen::internal::traits< Eigen::internal::remove_all_t< DerType >>::Scalar, product) > sqrt(const Eigen::AutoDiffScalar< DerType > &x)
const Eigen::CwiseUnaryOp< Eigen::internal::scalar_abs_op< typename Derived::Scalar >, const Derived > abs(const Eigen::ArrayBase< Derived > &x)
adouble abs(const adouble &x)
Definition: AdolcForward:74

◆ computeRef()

template<typename MatrixType_ >
void Eigen::IterScaling< MatrixType_ >::computeRef ( MatrixType mat)
inline

Compute the left and right vectors to scale the vectors the input matrix is scaled with the computed vectors at output

See also
compute()

Definition at line 151 of file Scaling.h.

152  {
153  compute (mat);
154  mat = m_matrix;
155  }

◆ init()

template<typename MatrixType_ >
void Eigen::IterScaling< MatrixType_ >::init ( )
inlineprotected

Definition at line 179 of file Scaling.h.

180  {
181  m_tol = 1e-10;
182  m_maxits = 5;
183  m_isInitialized = false;
184  }
Array< double, 1, 3 > e(1./3., 0.5, 2.)

◆ LeftScaling()

template<typename MatrixType_ >
VectorXd& Eigen::IterScaling< MatrixType_ >::LeftScaling ( )
inline

Get the vector to scale the rows of the matrix

Definition at line 158 of file Scaling.h.

159  {
160  return m_left;
161  }

◆ RightScaling()

template<typename MatrixType_ >
VectorXd& Eigen::IterScaling< MatrixType_ >::RightScaling ( )
inline

Get the vector to scale the columns of the matrix

Definition at line 165 of file Scaling.h.

166  {
167  return m_right;
168  }

◆ setTolerance()

template<typename MatrixType_ >
void Eigen::IterScaling< MatrixType_ >::setTolerance ( double  tol)
inline

Set the tolerance for the convergence of the iterative scaling algorithm

Definition at line 172 of file Scaling.h.

173  {
174  m_tol = tol;
175  }

Member Data Documentation

◆ m_info

template<typename MatrixType_ >
ComputationInfo Eigen::IterScaling< MatrixType_ >::m_info
mutableprotected

Definition at line 187 of file Scaling.h.

◆ m_isInitialized

template<typename MatrixType_ >
bool Eigen::IterScaling< MatrixType_ >::m_isInitialized
protected

Definition at line 188 of file Scaling.h.

◆ m_left

template<typename MatrixType_ >
VectorXd Eigen::IterScaling< MatrixType_ >::m_left
protected

Definition at line 189 of file Scaling.h.

◆ m_matrix

template<typename MatrixType_ >
MatrixType Eigen::IterScaling< MatrixType_ >::m_matrix
protected

Definition at line 186 of file Scaling.h.

◆ m_maxits

template<typename MatrixType_ >
int Eigen::IterScaling< MatrixType_ >::m_maxits
protected

Definition at line 192 of file Scaling.h.

◆ m_right

template<typename MatrixType_ >
VectorXd Eigen::IterScaling< MatrixType_ >::m_right
protected

Definition at line 190 of file Scaling.h.

◆ m_tol

template<typename MatrixType_ >
double Eigen::IterScaling< MatrixType_ >::m_tol
protected

Definition at line 191 of file Scaling.h.


The documentation for this class was generated from the following file: