MBDyn-1.7.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups
Functions | Variables
dgeequtest.cc File Reference
#include "mbconfig.h"
#include <iostream>
#include <iomanip>
#include <typeinfo>
#include "dgeequ.h"
#include "spmapmh.h"
#include "ccmh.h"
#include "dirccmh.h"
#include "naivemh.h"
Include dependency graph for dgeequtest.cc:

Go to the source code of this file.

Functions

void ReportMatScale (const char *title, bool fOK, const MatrixScaleBase &matScale, const MatrixHandler &mh, const double cond[2])
 
template<typename T >
void ScaleMatrix (const char *title, MatrixScale< T > &matScale, T &mh)
 
int main (int argc, char *argv[])
 

Variables

static doublereal mat [5][5]
 

Function Documentation

int main ( int  argc,
char *  argv[] 
)

Definition at line 74 of file dgeequtest.cc.

References SolutionManager::ScaleOpt::dTol, SolutionManager::ScaleOpt::iMaxIter, SpMapMatrixHandler::MakeCompressedColumnForm(), mat, NaiveMatrixHandler::Reset(), FullMatrixHandler::Reset(), SpMapMatrixHandler::Reset(), SolutionManager::SCALEF_VERBOSE, SolutionManager::SCALEF_WARN, ScaleMatrix(), grad::sqrt(), and SolutionManager::ScaleOpt::uFlags.

75 {
76  SolutionManager::ScaleOpt scale;
77  scale.uFlags |= SolutionManager::SCALEF_VERBOSE | SolutionManager::SCALEF_WARN;
78  scale.iMaxIter = argc >= 2 ? atoi(argv[1]) : 100;
79  scale.dTol = argc >= 3 ? atof(argv[2]) : sqrt(std::numeric_limits<doublereal>::epsilon());
80 
81  struct {
82  MatrixScale<NaiveMatrixHandler>* pNaive;
83  MatrixScale<NaivePermMatrixHandler>* pNaivePerm;
84  MatrixScale<FullMatrixHandler>* pFull;
85  MatrixScale<CColMatrixHandler<0> >* pCCol0;
86  MatrixScale<CColMatrixHandler<1> >* pCCol1;
87  MatrixScale<DirCColMatrixHandler<0> >* pDirCCol0;
88  MatrixScale<DirCColMatrixHandler<1> >* pDirCCol1;
89  MatrixScale<SpMapMatrixHandler>* pMap;
90  } matScale[] = {
91  { new RowMaxMatrixScale<NaiveMatrixHandler>(scale),
92  new RowMaxMatrixScale<NaivePermMatrixHandler>(scale),
93  new RowMaxMatrixScale<FullMatrixHandler>(scale),
94  new RowMaxMatrixScale<CColMatrixHandler<0> >(scale),
95  new RowMaxMatrixScale<CColMatrixHandler<1> >(scale),
96  new RowMaxMatrixScale<DirCColMatrixHandler<0> >(scale),
97  new RowMaxMatrixScale<DirCColMatrixHandler<1> >(scale),
98  new RowMaxMatrixScale<SpMapMatrixHandler>(scale)},
99  { new RowSumMatrixScale<NaiveMatrixHandler>(scale),
100  new RowSumMatrixScale<NaivePermMatrixHandler>(scale),
101  new RowSumMatrixScale<FullMatrixHandler>(scale),
102  new RowSumMatrixScale<CColMatrixHandler<0> >(scale),
103  new RowSumMatrixScale<CColMatrixHandler<1> >(scale),
104  new RowSumMatrixScale<DirCColMatrixHandler<0> >(scale),
105  new RowSumMatrixScale<DirCColMatrixHandler<1> >(scale),
106  new RowSumMatrixScale<SpMapMatrixHandler>(scale)},
107  { new ColMaxMatrixScale<NaiveMatrixHandler>(scale),
108  new ColMaxMatrixScale<NaivePermMatrixHandler>(scale),
109  new ColMaxMatrixScale<FullMatrixHandler>(scale),
110  new ColMaxMatrixScale<CColMatrixHandler<0> >(scale),
111  new ColMaxMatrixScale<CColMatrixHandler<1> >(scale),
112  new ColMaxMatrixScale<DirCColMatrixHandler<0> >(scale),
113  new ColMaxMatrixScale<DirCColMatrixHandler<1> >(scale),
114  new ColMaxMatrixScale<SpMapMatrixHandler>(scale)},
115  { new ColSumMatrixScale<NaiveMatrixHandler>(scale),
116  new ColSumMatrixScale<NaivePermMatrixHandler>(scale),
117  new ColSumMatrixScale<FullMatrixHandler>(scale),
118  new ColSumMatrixScale<CColMatrixHandler<0> >(scale),
119  new ColSumMatrixScale<CColMatrixHandler<1> >(scale),
120  new ColSumMatrixScale<DirCColMatrixHandler<0> >(scale),
121  new ColSumMatrixScale<DirCColMatrixHandler<1> >(scale),
122  new ColSumMatrixScale<SpMapMatrixHandler>(scale)},
123  { new LapackMatrixScale<NaiveMatrixHandler>(scale),
124  new LapackMatrixScale<NaivePermMatrixHandler>(scale),
125  new LapackMatrixScale<FullMatrixHandler>(scale),
126  new LapackMatrixScale<CColMatrixHandler<0> >(scale),
127  new LapackMatrixScale<CColMatrixHandler<1> >(scale),
128  new LapackMatrixScale<DirCColMatrixHandler<0> >(scale),
129  new LapackMatrixScale<DirCColMatrixHandler<1> >(scale),
130  new LapackMatrixScale<SpMapMatrixHandler>(scale)},
131  { new IterativeMatrixScale<NaiveMatrixHandler>(scale),
132  new IterativeMatrixScale<NaivePermMatrixHandler>(scale),
133  new IterativeMatrixScale<FullMatrixHandler>(scale),
134  new IterativeMatrixScale<CColMatrixHandler<0> >(scale),
135  new IterativeMatrixScale<CColMatrixHandler<1> >(scale),
136  new IterativeMatrixScale<DirCColMatrixHandler<0> >(scale),
137  new IterativeMatrixScale<DirCColMatrixHandler<1> >(scale),
138  new IterativeMatrixScale<SpMapMatrixHandler>(scale)},
139  { new RowMaxColMaxMatrixScale<NaiveMatrixHandler>(scale),
140  new RowMaxColMaxMatrixScale<NaivePermMatrixHandler>(scale),
141  new RowMaxColMaxMatrixScale<FullMatrixHandler>(scale),
142  new RowMaxColMaxMatrixScale<CColMatrixHandler<0> >(scale),
143  new RowMaxColMaxMatrixScale<CColMatrixHandler<1> >(scale),
144  new RowMaxColMaxMatrixScale<DirCColMatrixHandler<0> >(scale),
145  new RowMaxColMaxMatrixScale<DirCColMatrixHandler<1> >(scale),
146  new RowMaxColMaxMatrixScale<SpMapMatrixHandler>(scale)}
147  };
148 
149  const int N = sizeof(matScale)/sizeof(matScale[0]);
150 
151  for (int iMatScale = 0; iMatScale < N; ++iMatScale) {
152  std::vector<integer> perm(5), invperm(5);
153  perm[0] = 4;
154  perm[1] = 3;
155  perm[2] = 2;
156  perm[3] = 1;
157  perm[4] = 0;
158  for (int i = 0; i < 5; i++) {
159  invperm[perm[i]] = i;
160  }
161 
162  NaiveMatrixHandler nm(5);
163  NaivePermMatrixHandler npm(5, perm, invperm);
164  FullMatrixHandler fm(5);
165  SpMapMatrixHandler spm(5, 5);
166 
167  nm.Reset();
168  npm.Reset();
169  fm.Reset();
170  spm.Reset();
171  for (unsigned ir = 0; ir < 5; ir++) {
172  for (unsigned ic = 0; ic < 5; ic++) {
173  if (mat[ir][ic] != 0.) {
174  nm(ir + 1, ic + 1) = mat[ir][ic];
175  npm(ir + 1, ic + 1) = mat[ir][ic];
176  fm(ir + 1, ic + 1) = mat[ir][ic];
177  spm(ir + 1, ic + 1) = mat[ir][ic];
178  }
179  }
180  }
181 
182  std::vector<doublereal> Ax0, Ax1, Axd0, Axd1;
183  std::vector<integer> Ai0, Ai1, Ap0, Ap1, Aid0, Apd0, Aid1, Apd1;
184 
185  spm.MakeCompressedColumnForm(Ax0, Ai0, Ap0, 0);
186  spm.MakeCompressedColumnForm(Ax1, Ai1, Ap1, 1);
187  spm.MakeCompressedColumnForm(Axd0, Aid0, Apd0, 0);
188  spm.MakeCompressedColumnForm(Axd1, Aid1, Apd1, 1);
189 
190  CColMatrixHandler<0> ccm0(Ax0, Ai0, Ap0);
191  CColMatrixHandler<1> ccm1(Ax1, Ai1, Ap1);
192  DirCColMatrixHandler<0> dirccm0(Axd0, Aid0, Apd0);
193  DirCColMatrixHandler<1> dirccm1(Axd1, Aid1, Apd1);
194 
195  ScaleMatrix("Naive", *matScale[iMatScale].pNaive, nm);
196  ScaleMatrix("NaivePerm", *matScale[iMatScale].pNaivePerm, npm);
197  ScaleMatrix("Full", *matScale[iMatScale].pFull, fm);
198  ScaleMatrix("Dir0", *matScale[iMatScale].pDirCCol0, dirccm0);
199  ScaleMatrix("Dir1", *matScale[iMatScale].pDirCCol1, dirccm1);
200  ScaleMatrix("CC0", *matScale[iMatScale].pCCol0, ccm0);
201  ScaleMatrix("CC1", *matScale[iMatScale].pCCol1, ccm1);
202  ScaleMatrix("Map", *matScale[iMatScale].pMap, spm);
203  }
204 
205  for (unsigned i = 0; i < sizeof(matScale)/sizeof(matScale[0]); ++i) {
206  delete matScale[i].pNaive;
207  delete matScale[i].pNaivePerm;
208  delete matScale[i].pFull;
209  delete matScale[i].pCCol0;
210  delete matScale[i].pCCol1;
211  delete matScale[i].pDirCCol0;
212  delete matScale[i].pDirCCol1;
213  delete matScale[i].pMap;
214  }
215 
216  return 0;
217 }
ScaleMatrix
void ScaleMatrix(const char *title, MatrixScale< T > &matScale, T &mh)
Definition: dgeequtest.cc:56
mat
static doublereal mat[5][5]
Definition: dgeequtest.cc:45
grad::sqrt
GradientExpression< UnaryExpr< FuncSqrt, Expr > > sqrt(const GradientExpression< Expr > &u)
Definition: gradient.h:2974
doublereal
double doublereal
Definition: colamd.c:52

Here is the call graph for this function:

void ReportMatScale ( const char *  title,
bool  fOK,
const MatrixScaleBase &  matScale,
const MatrixHandler &  mh,
const double  cond[2] 
)

Definition at line 219 of file dgeequtest.cc.

References c, MatrixScaleBase::GetColScale(), MatrixScaleBase::GetRowScale(), MatrixHandler::iGetNumCols(), MatrixHandler::iGetNumRows(), and MatrixScaleBase::Report().

Referenced by ScaleMatrix().

219  {
220  std::cout << "------------------------------------------------------------" << std::endl;
221  std::cout << title << ": " << (fOK ? "OK" : "NOK") << " : " << typeid(matScale).name() << std::endl;
222 
223  std::cout << "condition number before scaling:" << cond[0] << std::endl;
224  std::cout << "condition number after scaling:" << cond[1] << std::endl;
225 
226  matScale.Report(std::cout);
227  const std::vector<doublereal>& r = matScale.GetRowScale();
228  const std::vector<doublereal>& c = matScale.GetColScale();
229 
230  const int N = std::min(mh.iGetNumRows(), mh.iGetNumCols());
231 
232  for (int i = 0; i < N; ++i) {
233  std::cout
234  << " r[" << i << "]=" << std::setw(12) << (r.empty() ? 1. : r[i])
235  << " c[" << i << "]=" << std::setw(12) << (c.empty() ? 1. : c[i])
236  << std::endl;
237  }
238 
239  std::cout << "matrix after scaling:" << std::endl << mh << std::endl;
240  std::cout << "------------------------------------------------------------" << std::endl;
241 }
MatrixHandler::iGetNumCols
virtual integer iGetNumCols(void) const =0
MatrixScaleBase::GetRowScale
const std::vector< doublereal > & GetRowScale() const
Definition: dgeequ.h:57
MatrixScaleBase::Report
std::ostream & Report(std::ostream &os) const
Definition: dgeequ.cc:52
MatrixScaleBase::GetColScale
const std::vector< doublereal > & GetColScale() const
Definition: dgeequ.h:58
c
static std::stack< cleanup * > c
Definition: cleanup.cc:59
MatrixHandler::iGetNumRows
virtual integer iGetNumRows(void) const =0

Here is the call graph for this function:

template<typename T >
void ScaleMatrix ( const char *  title,
MatrixScale< T > &  matScale,
T &  mh 
)

Definition at line 56 of file dgeequtest.cc.

References MatrixScale< T >::ComputeScaleFactors(), ReportMatScale(), and MatrixScale< T >::ScaleMatrix().

Referenced by main().

57 {
58  doublereal cond[2];
59 
60  //std::cout << "matrix before scaling:" << std::endl << mh << std::endl;
61 
62  cond[0] = mh.ConditionNumber();
63 
64  const bool fOK = matScale.ComputeScaleFactors(mh);
65 
66  matScale.ScaleMatrix(mh);
67 
68  cond[1] = mh.ConditionNumber();
69 
70  ReportMatScale(title, fOK, matScale, mh, cond);
71 }
MatrixScale::ComputeScaleFactors
bool ComputeScaleFactors(const T &mh)
Definition: dgeequ.h:305
ReportMatScale
void ReportMatScale(const char *title, bool fOK, const MatrixScaleBase &matScale, const MatrixHandler &mh, const double cond[2])
Definition: dgeequtest.cc:219
doublereal
double doublereal
Definition: colamd.c:52
MatrixScale::ScaleMatrix
T & ScaleMatrix(T &mh) const
Definition: dgeequ.h:274

Here is the call graph for this function:

Variable Documentation

doublereal mat[5][5]
static
Initial value:
= {
{ 11., 0., 13., 0., 15. },
{ 0., 22., 0., 24., 0. },
{ 31., 0., 33., 0., 35. },
{ 0., 42., 0., 44., 0. },
{ 51., 0., 53., 0., 55. }
}

Definition at line 45 of file dgeequtest.cc.

Referenced by ANN_vector_matrix_ass(), and main().