Operations.hpp

#pragma once
 
#include <boost/numeric/mtl/operation/two_norm.hpp>
 
#include <amdis/linearalgebra/mtl/MatrixBackend.hpp>
#include <amdis/linearalgebra/mtl/VectorBackend.hpp>
 
namespace AMDiS {
 
  // ||b - A*x||
  template <class T1, class T2, class T3>
  auto residuum(MTLSparseMatrix<T1> const& A, MTLVector<T2> const& x, MTLVector<T3> const& b)
  {
    auto r = b.vector();
    r -= A.matrix() * x.vector();
    return two_norm(r);
  }
 
  // ||b - A*x|| / ||b||
  template <class T1, class T2, class T3>
  auto relResiduum(MTLSparseMatrix<T1> const& A, MTLVector<T2> const& x, MTLVector<T3> const& b)
  {
    return residuum(A,x,b) / two_norm(b.vector());
  }
 
} // end namespace AMDiS