Ifpack2_BlockRelaxation_decl.hpp

Go to the documentation of this file.

/*@HEADER
// ***********************************************************************
//
//       Ifpack2: Templated Object-Oriented Algebraic Preconditioner Package
//                 Copyright (2009) Sandia Corporation
//
// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
// license for use of this work by or on behalf of the U.S. Government.
//
// Redistribution and use in source and binary forms, with or without
// modification, are permitted provided that the following conditions are
// met:
//
// 1. Redistributions of source code must retain the above copyright
// notice, this list of conditions and the following disclaimer.
//
// 2. Redistributions in binary form must reproduce the above copyright
// notice, this list of conditions and the following disclaimer in the
// documentation and/or other materials provided with the distribution.
//
// 3. Neither the name of the Corporation nor the names of the
// contributors may be used to endorse or promote products derived from
// this software without specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
//
// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
//
// ***********************************************************************
//@HEADER
*/

#ifndef IFPACK2_BLOCKRELAXATION_DECL_HPP
#define IFPACK2_BLOCKRELAXATION_DECL_HPP


#include "Ifpack2_Preconditioner.hpp"
#include "Ifpack2_Partitioner.hpp"
#include "Ifpack2_Details_CanChangeMatrix.hpp"
#include "Ifpack2_ContainerFactory.hpp"
#include "Tpetra_BlockCrsMatrix.hpp"
#include <type_traits>

namespace Ifpack2 {

template<class MatrixType, class ContainerType = Container<MatrixType> >
class BlockRelaxation :
    virtual public Ifpack2::Preconditioner<typename MatrixType::scalar_type,
                                           typename MatrixType::local_ordinal_type,
                                           typename MatrixType::global_ordinal_type,
                                           typename MatrixType::node_type>,
    virtual public Ifpack2::Details::CanChangeMatrix<Tpetra::RowMatrix<typename MatrixType::scalar_type,
                                                                       typename MatrixType::local_ordinal_type,
                                                                       typename MatrixType::global_ordinal_type,
                                                                       typename MatrixType::node_type> >
{
public:


  typedef typename MatrixType::scalar_type scalar_type;

  typedef typename MatrixType::local_ordinal_type local_ordinal_type;

  typedef typename MatrixType::global_ordinal_type global_ordinal_type;
  typedef typename MatrixType::node_type node_type;

  typedef typename Teuchos::ScalarTraits<scalar_type>::magnitudeType magnitude_type;

  typedef Tpetra::RowMatrix<scalar_type, local_ordinal_type, global_ordinal_type, node_type> row_matrix_type;

  static_assert (std::is_same<MatrixType, row_matrix_type>::value,
                 "Ifpack2::BlockRelaxation: Please use MatrixType = Tpetra::RowMatrix.");

  static_assert (std::is_same<ContainerType, Container<row_matrix_type> >::value,
                 "Ifpack2::BlockRelaxation: Do NOT specify the (second) "
                 "ContainerType template parameter explicitly.  The default "
                 "value is fine.  Please instead specify the container type to "
                 "use by setting the \"relaxation: container\" parameter.");

  typedef Tpetra::Import<local_ordinal_type, global_ordinal_type, node_type> import_type;

private:

  void setParametersImpl(Teuchos::ParameterList& params);

  void computeImporter() const;


  typedef Tpetra::MultiVector<scalar_type, local_ordinal_type,
                              global_ordinal_type, node_type> MV;
  typedef Tpetra::Vector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> vector_type;
  typedef Teuchos::ScalarTraits<scalar_type> STS;
  typedef Teuchos::ScalarTraits<magnitude_type> STM;
  typedef Tpetra::CrsMatrix<scalar_type, local_ordinal_type,
                                 global_ordinal_type, node_type> crs_matrix_type;
  typedef Tpetra::BlockCrsMatrix<scalar_type, local_ordinal_type,
                                 global_ordinal_type, node_type> block_crs_matrix_type;
  typedef Tpetra::Map<local_ordinal_type, global_ordinal_type, node_type> map_type;
public:

  // \name Constructors and Destructors

  explicit BlockRelaxation (const Teuchos::RCP<const row_matrix_type>& Matrix);

  virtual ~BlockRelaxation ();




  void setParameters(const Teuchos::ParameterList& params);

  bool supportsZeroStartingSolution() { return true; }

  void setZeroStartingSolution (bool zeroStartingSolution) { ZeroStartingSolution_ = zeroStartingSolution; };

  Teuchos::RCP<const Teuchos::ParameterList>
  getValidParameters () const;

  void initialize();

  inline bool isInitialized() const {
    return(IsInitialized_);
  }

  void compute();

  inline bool isComputed() const {
    return(IsComputed_);
  }



  virtual void
  setMatrix (const Teuchos::RCP<const row_matrix_type>& A);




  void apply(const MV& X,
             MV& Y,
             Teuchos::ETransp mode = Teuchos::NO_TRANS,
             scalar_type alpha = Teuchos::ScalarTraits<scalar_type>::one(),
             scalar_type beta = Teuchos::ScalarTraits<scalar_type>::zero()) const;

  Teuchos::RCP<const map_type> getDomainMap() const;

  Teuchos::RCP<const map_type> getRangeMap() const;

  bool hasTransposeApply() const;


  void applyMat(const MV& X,
                MV& Y,
                Teuchos::ETransp mode = Teuchos::NO_TRANS) const;



  Teuchos::RCP<const Teuchos::Comm<int> > getComm() const;

  Teuchos::RCP<const row_matrix_type> getMatrix() const;

  int getNumInitialize() const;

  int getNumCompute() const;

  int getNumApply() const;

  double getInitializeTime() const;

  double getComputeTime() const;

  double getApplyTime() const;

  size_t getNodeSmootherComplexity() const;  



  std::string description() const;

  void
  describe (Teuchos::FancyOStream& out,
            const Teuchos::EVerbosityLevel verbLevel =
            Teuchos::Describable::verbLevel_default) const;


  Teuchos::RCP<Ifpack2::Partitioner<Tpetra::RowGraph<local_ordinal_type,global_ordinal_type,node_type> > > getPartitioner(){return Partitioner_;}

private:

  BlockRelaxation (const BlockRelaxation<MatrixType, ContainerType> & RHS);

  BlockRelaxation<MatrixType,ContainerType>&
  operator= (const BlockRelaxation<MatrixType, ContainerType>& RHS);

  virtual void ApplyInverseJacobi (const MV& X, MV& Y) const;

  virtual void ApplyInverseGS (const MV& X, MV& Y) const;

  virtual void ApplyInverseSGS (const MV& X, MV& Y) const;

  void ExtractSubmatricesStructure();



  Teuchos::RCP<const row_matrix_type> A_;

  mutable Teuchos::RCP<Container<row_matrix_type> > Container_;

  // FIXME (mfh 06 Oct 2014) This doesn't comply with the naming
  // convention for instance members of a class.  Furthermore, the
  // class should keep the Vector, not the ArrayRCP to the data _in_
  // the Vector.
  // FIXED! (amk 10 Nov 2015)
  mutable Teuchos::RCP<vector_type> DiagRCP_;

  Teuchos::RCP<Ifpack2::Partitioner<Tpetra::RowGraph<local_ordinal_type,global_ordinal_type,node_type> > > Partitioner_;

  std::string PartitionerType_;

  Teuchos::ParameterList List_;

  int NumSweeps_;

  local_ordinal_type NumLocalBlocks_;

  std::string containerType_;

  Details::RelaxationType PrecType_;

  bool IsParallel_;

  bool ZeroStartingSolution_;

  bool hasBlockCrsMatrix_;

  bool DoBackwardGS_;

  int OverlapLevel_;

  //only valid with block Jacobi relaxation and overlapping blocks (e.g., defined by "partitioner: type" "user" and "parts: " or "global ID parts:"). Average solutions in overlapped regions (i.e., after summing different solutions divide by number of blocks contain this dof). When false (the default) symmetric averaging performed (i.e., average residuals and solutions).
  bool nonsymCombine_;

  std::string schwarzCombineMode_;

  scalar_type DampingFactor_;

  bool decouple_;

  bool IsInitialized_;

  bool IsComputed_;

  int NumInitialize_;

  int NumCompute_;

  bool TimerForApply_;

  mutable int NumApply_;

  double InitializeTime_;

  double ComputeTime_;

  mutable double ApplyTime_;

  local_ordinal_type NumLocalRows_;

  global_ordinal_type NumGlobalRows_;

  global_ordinal_type NumGlobalNonzeros_;

  Teuchos::RCP<vector_type> W_;

  mutable Teuchos::RCP<const Tpetra::Import<local_ordinal_type, global_ordinal_type, node_type>> Importer_;

}; //class BlockRelaxation

}//namespace Ifpack2

#endif // IFPACK2_BLOCKRELAXATION_DECL_HPP