relaxation_block.h

// ---------------------------------------------------------------------
//
// Copyright (C) 2010 - 2018 by the deal.II authors
//
// This file is part of the deal.II library.
//
// The deal.II library is free software; you can use it, redistribute
// it, and/or modify it under the terms of the GNU Lesser General
// Public License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
// The full text of the license can be found in the file LICENSE at
// the top level of the deal.II distribution.
//
// ---------------------------------------------------------------------

#ifndef dealii_relaxation_block_h
#define dealii_relaxation_block_h

#include <deal.II/base/subscriptor.h>
#include <deal.II/base/smartpointer.h>
#include <deal.II/lac/vector.h>
#include <deal.II/lac/precondition_block_base.h>
#include <deal.II/lac/sparsity_pattern.h>

#include <vector>
#include <set>

DEAL_II_NAMESPACE_OPEN

template <typename MatrixType,
          typename InverseNumberType = typename MatrixType::value_type,
          typename VectorType = Vector<double> >
class RelaxationBlock :
  protected PreconditionBlockBase<InverseNumberType>
{
private:
  typedef typename MatrixType::value_type number;

  typedef InverseNumberType value_type;

public:
  typedef types::global_dof_index size_type;

  class AdditionalData : public Subscriptor
  {
  public:
    AdditionalData (const double relaxation = 1.,
                    const bool invert_diagonal = true,
                    const bool same_diagonal = false,
                    const typename PreconditionBlockBase<InverseNumberType>::Inversion inversion
                    = PreconditionBlockBase<InverseNumberType>::gauss_jordan,
                    const double threshold = 0.,
                    VectorType *temp_ghost_vector = nullptr);

    SparsityPattern block_list;

    double relaxation;

    bool invert_diagonal;

    bool same_diagonal;

    typename PreconditionBlockBase<InverseNumberType>::Inversion inversion;

    double threshold = 0.;

    unsigned int kernel_size = 0;

    std::vector<std::vector<unsigned int> > order;

    mutable VectorType *temp_ghost_vector;

    std::size_t memory_consumption() const;
  };

  void initialize (const MatrixType     &A,
                   const AdditionalData &parameters);

  void clear();

  void invert_diagblocks();

protected:
  void do_step (
    VectorType       &dst,
    const VectorType &prev,
    const VectorType &src,
    const bool backward) const;

  SmartPointer<const MatrixType,RelaxationBlock<MatrixType, InverseNumberType, VectorType> > A;

  SmartPointer<const AdditionalData, RelaxationBlock<MatrixType, InverseNumberType, VectorType> > additional_data;

private:
  void block_kernel(const size_type block_begin, const size_type block_end);
};


template <typename MatrixType,
          typename InverseNumberType = typename MatrixType::value_type,
          typename VectorType = Vector<double> >
class RelaxationBlockJacobi : public virtual Subscriptor,
  protected RelaxationBlock<MatrixType, InverseNumberType, VectorType>
{
public:
//    RelaxationBlockJacobi();

  typedef typename MatrixType::value_type number;

  using typename RelaxationBlock<MatrixType, InverseNumberType, VectorType>::AdditionalData;

  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::initialize;

  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::clear;

  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::size;
  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::inverse;
  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::inverse_householder;
  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::inverse_svd;
  using PreconditionBlockBase<InverseNumberType>::log_statistics;
  void step (VectorType &dst, const VectorType &rhs) const;

  void Tstep (VectorType &dst, const VectorType &rhs) const;

  void vmult (VectorType &dst, const VectorType &rhs) const;

  void Tvmult (VectorType &dst, const VectorType &rhs) const;
};


template <typename MatrixType,
          typename InverseNumberType = typename MatrixType::value_type,
          typename VectorType = Vector<double> >
class RelaxationBlockSOR : public virtual Subscriptor,
  protected RelaxationBlock<MatrixType, InverseNumberType, VectorType>
{
public:
//    RelaxationBlockSOR();

  typedef typename MatrixType::value_type number;

  using typename RelaxationBlock<MatrixType, InverseNumberType, VectorType>::AdditionalData;

  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::initialize;

  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::clear;

  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::size;
  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::inverse;
  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::inverse_householder;
  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::inverse_svd;
  using PreconditionBlockBase<InverseNumberType>::log_statistics;
  void step (VectorType &dst, const VectorType &rhs) const;

  void Tstep (VectorType &dst, const VectorType &rhs) const;

  void vmult (VectorType &dst, const VectorType &rhs) const;

  void Tvmult (VectorType &dst, const VectorType &rhs) const;
};


template <typename MatrixType,
          typename InverseNumberType = typename MatrixType::value_type,
          typename VectorType = Vector<double> >
class RelaxationBlockSSOR : public virtual Subscriptor,
  protected RelaxationBlock<MatrixType, InverseNumberType, VectorType>
{
public:
  typedef typename MatrixType::value_type number;

  using typename RelaxationBlock<MatrixType, InverseNumberType, VectorType>::AdditionalData;

  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::initialize;

  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::clear;

  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::size;
  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::inverse;
  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::inverse_householder;
  using RelaxationBlock<MatrixType, InverseNumberType, VectorType>::inverse_svd;
  using PreconditionBlockBase<InverseNumberType>::log_statistics;
  void step (VectorType &dst, const VectorType &rhs) const;

  void Tstep (VectorType &dst, const VectorType &rhs) const;

  void vmult (VectorType &dst, const VectorType &rhs) const;

  void Tvmult (VectorType &dst, const VectorType &rhs) const;
};


DEAL_II_NAMESPACE_CLOSE

#endif