17 #ifndef __deal2__solver_qmrs_h
18 #define __deal2__solver_qmrs_h
20 #include <deal.II/base/config.h>
21 #include <deal.II/lac/solver.h>
22 #include <deal.II/lac/solver_control.h>
23 #include <deal.II/base/logstream.h>
25 #include <deal.II/base/subscriptor.h>
29 DEAL_II_NAMESPACE_OPEN
66 template <
class VECTOR = Vector<
double> >
135 template<
class MATRIX,
class PRECONDITIONER>
140 const PRECONDITIONER &precondition);
154 const VECTOR &d)
const;
201 template<
class MATRIX,
class PRECONDITIONER>
203 iterate(
const MATRIX &A,
const PRECONDITIONER &precondition);
215 template<
class VECTOR>
218 const AdditionalData &data)
221 additional_data(data)
226 template<
class VECTOR>
228 const AdditionalData &data)
231 additional_data(data)
236 template<
class VECTOR>
240 return std::sqrt(res2);
245 template<
class VECTOR>
250 const VECTOR &)
const
255 template<
class VECTOR>
256 template<
class MATRIX,
class PRECONDITIONER>
261 const PRECONDITIONER &precondition)
263 deallog.
push(
"QMRS");
266 Vv = this->memory.alloc();
267 Vp = this->memory.alloc();
268 Vq = this->memory.alloc();
269 Vt = this->memory.alloc();
270 Vd = this->memory.alloc();
289 deallog <<
"Restart step " << step << std::endl;
290 state = iterate(A, precondition);
295 this->memory.free(Vv);
296 this->memory.free(Vp);
297 this->memory.free(Vq);
298 this->memory.free(Vt);
299 this->memory.free(Vd);
307 this->control().last_value()));
313 template<
class VECTOR>
314 template<
class MATRIX,
class PRECONDITIONER>
317 const PRECONDITIONER &precondition)
334 const VECTOR &b = *Vb;
338 double tau, rho, theta=0, sigma, alpha, psi, theta_old, rho_old, beta;
344 precondition.vmult(q,x);
355 precondition.vmult(q,p);
370 if (std::fabs(sigma/rho) < additional_data.breakdown)
382 d.sadd(psi*theta_old, psi*alpha, p);
385 print_vectors(step,x,v,d);
387 if (additional_data.exact_residual)
394 res = std::sqrt((it+1)*tau);
395 state = this->control().check(step,res);
400 if (std::fabs(rho) < additional_data.breakdown)
404 precondition.vmult(q,v);
409 precondition.vmult(q,p);
416 DEAL_II_NAMESPACE_CLOSE
Stop iteration, goal not reached.
void vmult(VECTOR &u, const VECTOR &v) const
virtual void print_vectors(const unsigned int step, const VECTOR &x, const VECTOR &r, const VECTOR &d) const
AdditionalData additional_data
#define AssertThrow(cond, exc)
virtual double criterion()
Stop iteration, goal reached.
bool iterate(const MATRIX &A, const PRECONDITIONER &precondition)
void push(const std::string &text)
void solve(const MATRIX &A, VECTOR &x, const VECTOR &b, const PRECONDITIONER &precondition)
AdditionalData(bool exact_residual=false, double breakdown=1.e-16)
SolverQMRS(SolverControl &cn, VectorMemory< VECTOR > &mem, const AdditionalData &data=AdditionalData())