1 #ifndef DUNE_PARDISO_HH
2 #define DUNE_PARDISO_HH
10 #define F77_FUNC(func) func
12 #define F77_FUNC(func) func ## _
19 (
void *,
int *,
int *);
22 (
void *,
int *,
int *,
int *,
int *,
int *,
23 double *,
int *,
int *,
int *,
int *,
int *,
24 int *,
double *,
double *,
int *);
34 template<
class M,
class X,
class Y>
78 if (A_.rowdim(i.index()) != 1)
79 DUNE_THROW(NotImplemented,
"SeqPardiso: row blocksize != 1.");
81 for (
ColIterator j = (*i).begin(); j != endj; ++j) {
82 if (A_.coldim(j.index()) != 1)
83 DUNE_THROW(NotImplemented,
"SeqPardiso: column blocksize != 1.");
88 std::cout <<
"dimension = " << n_ <<
", number of nonzeros = " << nnz << std::endl;
97 ia_[i.index()] = count+1;
99 for (
ColIterator j = (*i).begin(); j != endj; ++j) {
101 ja_[
count] = j.index()+1;
108 F77_FUNC(pardisoinit) (pt_, &mtype_, iparm_);
113 iparm_[2] = num_procs_;
115 F77_FUNC(pardiso) (pt_, &maxfct_, &mnum_, &mtype_, &phase,
116 &n_, a_, ia_, ja_, &idum, &nrhs_,
117 iparm_, &msglvl_, &ddum, &ddum, &error_);
120 DUNE_THROW(MathError,
"Constructor SeqPardiso: Factorization failed. Error code " << error_);
122 std::cout <<
"Constructor SeqPardiso: Factorization completed." << std::endl;
125 DUNE_THROW(NotImplemented,
"no Pardiso library available, reconfigure with correct --with-pardiso options");
134 virtual void pre (X& x, Y& b) {}
141 virtual void apply (X& v,
const Y& d)
151 for (
int i = 0; i < n_; i++) {
156 F77_FUNC(pardiso) (pt_, &maxfct_, &mnum_, &mtype_, &phase,
157 &n_, a_, ia_, ja_, &idum, &nrhs_,
158 iparm_, &msglvl_, b, x, &error_);
161 DUNE_THROW(MathError,
"SeqPardiso.apply: Backsolve failed. Error code " << error_);
163 for (
int i = 0; i < n_; i++)
166 std::cout <<
"SeqPardiso: Backsolve completed." << std::endl;
184 F77_FUNC(pardiso) (pt_, &maxfct_, &mnum_, &mtype_, &phase,
185 &n_, &ddum, ia_, ja_, &idum, &nrhs_,
186 iparm_, &msglvl_, &ddum, &ddum, &error_);
210 template<
class M,
class X,
class Y>