ir.h

Go to the documentation of this file.

// $Id: ir.h,v 1.7.2.3 2019/05/28 11:13:02 garloff Exp $

#ifndef TBCI_SOLVER_IR_H
#define TBCI_SOLVER_IR_H

#include "tbci/basics.h"

NAMESPACE_TBCI

INST3(template <typename T, Matrix<T>, Vector<T> > class NN friend \
int IR (const Matrix<T>&, Vector<T>&, const Vector<T>&,\
        const Preconditioner_Sig<T, Matrix<T> >&, unsigned int&, double&);)
INST3(template <typename T, BdMatrix<T>, Vector<T> > class NN friend \
int IR (const BdMatrix<T>&, Vector<T>&, const Vector<T>&,\
        const Preconditioner_Sig<T, BdMatrix<T> >&, unsigned int&, double&);)


//Aenderung Marco:
template < typename T, typename SysMatrix, typename SysVector >
int IR (const SysMatrix &A, SysVector &x, const SysVector &b,
        const Preconditioner_Sig<T, SysMatrix> &M, unsigned int &max_iter, double &tol)
{
  double resid;
  unsigned int dim = A.rows();
  SysVector z(dim);

  double normb = b.fabs();
  if (normb == 0.0)  normb = 1.0;
  SysVector r(b - A*x);

  
  if ((resid = r.fabs() / normb) <= tol) 
  {
    tol = resid;
    max_iter = 0;
    return 0;
  }
  
  for (unsigned int i = 1; i <= max_iter; i++) 
  {
    z = M.solve(r);
    x += z;
    r = b - A * x;
    if ((resid = r.fabs() / normb) <= tol) 
    {
      tol = resid;
      max_iter = i;
      return 0;
    }
  }

  tol = resid;
  return 1;
}

NAMESPACE_END

#endif