scalarproducts.hh

Go to the documentation of this file.

// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_ISTL_SCALARPRODUCTS_HH
#define DUNE_ISTL_SCALARPRODUCTS_HH
 
#include <cmath>
#include <complex>
#include <iostream>
#include <iomanip>
#include <string>
#include <memory>
 
#include <dune/common/exceptions.hh>
#include <dune/common/shared_ptr.hh>
 
#include "bvector.hh"
#include "solvercategory.hh"
 
 
namespace Dune {
  template<class X>
  class ScalarProduct {
  public:
    typedef X domain_type;
    typedef typename X::field_type field_type;
    typedef typename FieldTraits<field_type>::real_type real_type;
 
    virtual field_type dot (const X& x, const X& y) const
    {
      return x.dot(y);
    }
 
    virtual real_type norm (const X& x) const
    {
      return x.two_norm();
    }
 
    virtual SolverCategory::Category category() const
    {
      return SolverCategory::sequential;
    }
 
    virtual ~ScalarProduct () {}
  };
 
  template<class X, class C>
  class ParallelScalarProduct : public ScalarProduct<X>
  {
  public:
    typedef X domain_type;
    typedef typename X::field_type field_type;
    typedef typename FieldTraits<field_type>::real_type real_type;
    typedef C communication_type;
 
    ParallelScalarProduct (std::shared_ptr<const communication_type> com, SolverCategory::Category cat)
      : _communication(com), _category(cat)
    {}
 
    ParallelScalarProduct (const communication_type& com, SolverCategory::Category cat)
      : ParallelScalarProduct(stackobject_to_shared_ptr(com), cat)
    {}
 
 
    virtual field_type dot (const X& x, const X& y) const override
    {
      field_type result(0);
      _communication->dot(x,y,result); // explicitly loop and apply masking
      return result;
    }
 
    virtual real_type norm (const X& x) const override
    {
      return _communication->norm(x);
    }
 
    virtual SolverCategory::Category category() const override
    {
      return _category;
    }
 
  private:
    std::shared_ptr<const communication_type> _communication;
    SolverCategory::Category _category;
  };
 
  template<class X>
  class SeqScalarProduct : public ScalarProduct<X>
  {
    using ScalarProduct<X>::ScalarProduct;
  };
 
  template<class X, class C>
  class NonoverlappingSchwarzScalarProduct : public ParallelScalarProduct<X,C>
  {
  public:
    NonoverlappingSchwarzScalarProduct (std::shared_ptr<const C> comm) :
      ParallelScalarProduct<X,C>(comm,SolverCategory::nonoverlapping) {}
 
    NonoverlappingSchwarzScalarProduct (const C& comm) :
      ParallelScalarProduct<X,C>(comm,SolverCategory::nonoverlapping) {}
  };
 
  template<class X, class C>
  class OverlappingSchwarzScalarProduct : public ParallelScalarProduct<X,C>
  {
  public:
    OverlappingSchwarzScalarProduct (std::shared_ptr<const C> comm) :
      ParallelScalarProduct<X,C>(comm, SolverCategory::overlapping) {}
 
    OverlappingSchwarzScalarProduct (const C& comm) :
      ParallelScalarProduct<X,C>(comm,SolverCategory::overlapping) {}
  };
 
  template<class X, class Comm>
  std::shared_ptr<ScalarProduct<X>> makeScalarProduct(std::shared_ptr<const Comm> comm, SolverCategory::Category category)
  {
    switch(category)
    {
      case SolverCategory::sequential:
        return
          std::make_shared<ScalarProduct<X>>();
      default:
        return
          std::make_shared<ParallelScalarProduct<X,Comm>>(comm,category);
    }
  }
 
  template<class X, class Comm>
  std::shared_ptr<ScalarProduct<X>> createScalarProduct(const Comm& comm, SolverCategory::Category category)
  { return makeScalarProduct<X>(stackobject_to_shared_ptr(comm), category); }
 
} // end namespace Dune
 
#endif