3 #ifndef DUNE_ISTL_PARDISO_HH
4 #define DUNE_ISTL_PARDISO_HH
12 #define F77_FUNC(func) func
14 #define F77_FUNC(func) func ## _
21 (
void *,
int *,
int *);
24 (
void *,
int *,
int *,
int *,
int *,
int *,
25 double *,
int *,
int *,
int *,
int *,
int *,
26 int *,
double *,
double *,
int *);
36 template<
class M,
class X,
class Y>
77 RowIterator endi = A_.end();
78 for (RowIterator i = A_.begin(); i != endi; ++i)
80 if (A_.rowdim(i.index()) != 1)
81 DUNE_THROW(NotImplemented,
"SeqPardiso: row blocksize != 1.");
82 ColIterator endj = (*i).end();
83 for (ColIterator j = (*i).begin(); j != endj; ++j) {
84 if (A_.coldim(j.index()) != 1)
85 DUNE_THROW(NotImplemented,
"SeqPardiso: column blocksize != 1.");
90 std::cout <<
"dimension = " << n_ <<
", number of nonzeros = " << nnz << std::endl;
97 for (RowIterator i = A_.begin(); i != endi; ++i)
99 ia_[i.index()] = count+1;
100 ColIterator endj = (*i).end();
101 for (ColIterator j = (*i).begin(); j != endj; ++j) {
103 ja_[
count] = j.index()+1;
110 F77_FUNC(pardisoinit) (pt_, &mtype_, iparm_);
115 iparm_[2] = num_procs_;
117 F77_FUNC(pardiso) (pt_, &maxfct_, &mnum_, &mtype_, &phase,
118 &n_, a_, ia_, ja_, &idum, &nrhs_,
119 iparm_, &msglvl_, &ddum, &ddum, &error_);
122 DUNE_THROW(MathError,
"Constructor SeqPardiso: Factorization failed. Error code " << error_);
124 std::cout <<
"Constructor SeqPardiso: Factorization completed." << std::endl;
127 DUNE_THROW(NotImplemented,
"no Pardiso library available, reconfigure with correct --with-pardiso options");
136 virtual void pre (X& x, Y& b) {}
143 virtual void apply (X& v,
const Y& d)
153 for (
int i = 0; i < n_; i++) {
158 F77_FUNC(pardiso) (pt_, &maxfct_, &mnum_, &mtype_, &phase,
159 &n_, a_, ia_, ja_, &idum, &nrhs_,
160 iparm_, &msglvl_, b, x, &error_);
163 DUNE_THROW(MathError,
"SeqPardiso.apply: Backsolve failed. Error code " << error_);
165 for (
int i = 0; i < n_; i++)
168 std::cout <<
"SeqPardiso: Backsolve completed." << std::endl;
186 F77_FUNC(pardiso) (pt_, &maxfct_, &mnum_, &mtype_, &phase,
187 &n_, &ddum, ia_, ja_, &idum, &nrhs_,
188 iparm_, &msglvl_, &ddum, &ddum, &error_);
212 template<
class M,
class X,
class Y>
X domain_type
The domain type of the preconditioner.
Definition: pardiso.hh:42
~SeqPardiso()
Definition: pardiso.hh:179
M::RowIterator RowIterator
Definition: pardiso.hh:48
Base class for matrix free definition of preconditioners.
Definition: preconditioner.hh:26
Y range_type
The range type of the preconditioner.
Definition: pardiso.hh:44
Matrix & A
Definition: matrixmatrix.hh:216
M::ColIterator ColIterator
Definition: pardiso.hh:49
Templates characterizing the type of a solver.
Category for sequential solvers.
Definition: solvercategory.hh:21
X::field_type field_type
The field type of the preconditioner.
Definition: pardiso.hh:46
Whether this is a direct solver.
Definition: solvertype.hh:22
std::size_t count
Definition: matrixmatrix.hh:215
#define F77_FUNC(func)
Definition: pardiso.hh:14
SeqPardiso(const M &A)
Constructor.
Definition: pardiso.hh:62
Definition: solvertype.hh:13
Define general preconditioner interface.
virtual void pre(X &x, Y &b)
Prepare the preconditioner.
Definition: pardiso.hh:136
virtual void post(X &x)
Clean up.
Definition: pardiso.hh:177
virtual void apply(X &v, const Y &d)
Apply the preconditioner.
Definition: pardiso.hh:143
M matrix_type
The matrix type the preconditioner is for.
Definition: pardiso.hh:40
The sequential Pardiso preconditioner.
Definition: pardiso.hh:37
Definition: basearray.hh:19
The category the preconditioner is part of.
Definition: pardiso.hh:54