Xpetra_BlockedCrsMatrix_decl.hpp

Go to the documentation of this file.

// @HEADER
// *****************************************************************************
//             Xpetra: A linear algebra interface package
//
// Copyright 2012 NTESS and the Xpetra contributors.
// SPDX-License-Identifier: BSD-3-Clause
// *****************************************************************************
// @HEADER

#ifndef XPETRA_BLOCKEDCRSMATRIX_DECL_HPP
#define XPETRA_BLOCKEDCRSMATRIX_DECL_HPP

#include <Tpetra_KokkosCompat_DefaultNode.hpp>

#include <Teuchos_SerialDenseMatrix.hpp>
#include <Teuchos_Hashtable.hpp>

#include "Xpetra_ConfigDefs.hpp"
#include "Xpetra_Exceptions.hpp"

#include "Xpetra_MapFactory.hpp"
#include "Xpetra_MultiVector.hpp"
#include "Xpetra_BlockedMultiVector.hpp"
#include "Xpetra_MultiVectorFactory.hpp"
#include "Xpetra_BlockedVector.hpp"
#include "Xpetra_CrsGraph.hpp"
#include "Xpetra_CrsMatrix.hpp"
#include "Xpetra_CrsMatrixFactory.hpp"

#include "Xpetra_MapExtractor.hpp"
#include "Xpetra_MapExtractorFactory.hpp"

#include "Xpetra_Matrix.hpp"
#include "Xpetra_MatrixFactory.hpp"
#include "Xpetra_CrsMatrixWrap.hpp"

#ifdef HAVE_XPETRA_THYRA
#include <Thyra_ProductVectorSpaceBase.hpp>
#include <Thyra_VectorSpaceBase.hpp>
#include <Thyra_LinearOpBase.hpp>
#include <Thyra_BlockedLinearOpBase.hpp>
#include <Thyra_PhysicallyBlockedLinearOpBase.hpp>
#include "Xpetra_ThyraUtils.hpp"
#endif

#include "Xpetra_VectorFactory.hpp"

namespace Xpetra {

#ifdef HAVE_XPETRA_THYRA
template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
class ThyraUtils;
#endif

typedef std::string viewLabel_t;

template <class Scalar,
          class LocalOrdinal,
          class GlobalOrdinal,
          class Node = Tpetra::KokkosClassic::DefaultNode::DefaultNodeType>
class BlockedCrsMatrix : public Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> {
 public:
  typedef Scalar scalar_type;
  typedef LocalOrdinal local_ordinal_type;
  typedef GlobalOrdinal global_ordinal_type;
  typedef Node node_type;

 private:
#undef XPETRA_BLOCKEDCRSMATRIX_SHORT
#include "Xpetra_UseShortNames.hpp"

 public:



  BlockedCrsMatrix(const Teuchos::RCP<const BlockedMap>& rangeMaps,
                   const Teuchos::RCP<const BlockedMap>& domainMaps,
                   size_t numEntriesPerRow);


  BlockedCrsMatrix(Teuchos::RCP<const MapExtractor>& rangeMapExtractor,
                   Teuchos::RCP<const MapExtractor>& domainMapExtractor,
                   size_t numEntriesPerRow);

#ifdef HAVE_XPETRA_THYRA

  BlockedCrsMatrix(const Teuchos::RCP<const Thyra::BlockedLinearOpBase<Scalar> >& thyraOp, const Teuchos::RCP<const Teuchos::Comm<int> >& /* comm */);

 private:

  Teuchos::RCP<const Map> mergeMaps(std::vector<Teuchos::RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > >& subMaps);

 public:
#endif

  virtual ~BlockedCrsMatrix();





  void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView<const GlobalOrdinal>& cols, const ArrayView<const Scalar>& vals);


  void insertLocalValues(LocalOrdinal localRow, const ArrayView<const LocalOrdinal>& cols, const ArrayView<const Scalar>& vals);

  void removeEmptyProcessesInPlace(const Teuchos::RCP<const Map>& newMap);


  void replaceGlobalValues(GlobalOrdinal globalRow,
                           const ArrayView<const GlobalOrdinal>& cols,
                           const ArrayView<const Scalar>& vals);


  void replaceLocalValues(LocalOrdinal localRow,
                          const ArrayView<const LocalOrdinal>& cols,
                          const ArrayView<const Scalar>& vals);

  virtual void setAllToScalar(const Scalar& alpha);

  void scale(const Scalar& alpha);




  void resumeFill(const RCP<ParameterList>& params = null);

  void fillComplete(const RCP<const Map>& domainMap, const RCP<const Map>& rangeMap, const RCP<ParameterList>& params = null);

  void fillComplete(const RCP<ParameterList>& params = null);



  global_size_t getGlobalNumRows() const;


  global_size_t getGlobalNumCols() const;

  size_t getLocalNumRows() const;

  global_size_t getGlobalNumEntries() const;

  size_t getLocalNumEntries() const;


  size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const;


  size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const;


  size_t getGlobalMaxNumRowEntries() const;


  size_t getLocalMaxNumRowEntries() const;


  bool isLocallyIndexed() const;


  bool isGloballyIndexed() const;

  bool isFillComplete() const;


  virtual void getLocalRowCopy(LocalOrdinal LocalRow,
                               const ArrayView<LocalOrdinal>& Indices,
                               const ArrayView<Scalar>& Values,
                               size_t& NumEntries) const;


  void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView<const GlobalOrdinal>& indices, ArrayView<const Scalar>& values) const;


  void getLocalRowView(LocalOrdinal LocalRow, ArrayView<const LocalOrdinal>& indices, ArrayView<const Scalar>& values) const;


  void getLocalDiagCopy(Vector& diag) const;

  void leftScale(const Vector& x);

  void rightScale(const Vector& x);

  virtual typename ScalarTraits<Scalar>::magnitudeType getFrobeniusNorm() const;

  virtual bool haveGlobalConstants() const;








  virtual void apply(const MultiVector& X, MultiVector& Y, Teuchos::ETransp mode, Scalar alpha, Scalar beta, bool sumInterfaceValues,
                     const RCP<Xpetra::Import<LocalOrdinal, GlobalOrdinal, Node> >& regionInterfaceImporter,
                     const Teuchos::ArrayRCP<LocalOrdinal>& regionInterfaceLIDs) const;


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

  RCP<const Map> getFullDomainMap() const;

  RCP<const BlockedMap> getBlockedDomainMap() const;

  const RCP<const Map> getDomainMap() const;

  RCP<const Map> getDomainMap(size_t i) const;

  RCP<const Map> getDomainMap(size_t i, bool bThyraMode) const;

  RCP<const Map> getFullRangeMap() const;

  RCP<const BlockedMap> getBlockedRangeMap() const;

  const RCP<const Map> getRangeMap() const;

  RCP<const Map> getRangeMap(size_t i) const;

  RCP<const Map> getRangeMap(size_t i, bool bThyraMode) const;

  RCP<const MapExtractor> getRangeMapExtractor() const;

  RCP<const MapExtractor> getDomainMapExtractor() const;


  //{@

  virtual void bgs_apply(
      const MultiVector& X,                                 
      MultiVector& Y,                                       
      size_t row,                                           
      Teuchos::ETransp mode = Teuchos::NO_TRANS,            
      Scalar alpha          = ScalarTraits<Scalar>::one(),  
      Scalar beta           = ScalarTraits<Scalar>::zero()  
  ) const;


  //{@

  const Teuchos::RCP<const Map> getMap() const;

  void doImport(const Matrix& source, const Import& importer, CombineMode CM);

  void doExport(const Matrix& dest, const Import& importer, CombineMode CM);

  void doImport(const Matrix& source, const Export& exporter, CombineMode CM);

  void doExport(const Matrix& dest, const Export& exporter, CombineMode CM);

  // @}



  std::string description() const;

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


  void setObjectLabel(const std::string& objectLabel);

  bool hasCrsGraph() const;

  RCP<const CrsGraph> getCrsGraph() const;




  virtual bool isDiagonal() const;

  virtual size_t Rows() const;

  virtual size_t Cols() const;

  Teuchos::RCP<Matrix> getCrsMatrix() const;

  Teuchos::RCP<Matrix> getInnermostCrsMatrix();

  Teuchos::RCP<Matrix> getMatrix(size_t r, size_t c) const;

  // void setMatrix(size_t r, size_t c, Teuchos::RCP<CrsMatrix> mat) {
  void setMatrix(size_t r, size_t c, Teuchos::RCP<Matrix> mat);

  // NOTE: This is a rather expensive operation, since all blocks are copied
  // into a new big CrsMatrix
  Teuchos::RCP<Matrix> Merge() const;

  typedef typename CrsMatrix::local_matrix_type local_matrix_type;
  local_matrix_type getLocalMatrixDevice() const;
  typename local_matrix_type::HostMirror getLocalMatrixHost() const;

#ifdef HAVE_XPETRA_THYRA
  Teuchos::RCP<Thyra::BlockedLinearOpBase<Scalar> > getThyraOperator();
#endif
  LocalOrdinal GetStorageBlockSize() const;

  void residual(const MultiVector& X,
                const MultiVector& B,
                MultiVector& R) const;

 private:


  void Add(const Matrix& A, const Scalar scalarA, Matrix& B, const Scalar scalarB) const;


  // Default view is created after fillComplete()
  // Because ColMap might not be available before fillComplete().
  void CreateDefaultView();

 private:
  bool is_diagonal_;                             
  Teuchos::RCP<const MapExtractor> domainmaps_;  
  Teuchos::RCP<const MapExtractor> rangemaps_;   

  std::vector<Teuchos::RCP<Matrix> > blocks_;  
#ifdef HAVE_XPETRA_THYRA
  Teuchos::RCP<const Thyra::BlockedLinearOpBase<Scalar> > thyraOp_;  
#endif
  bool bRangeThyraMode_;   
  bool bDomainThyraMode_;  
};

}  // namespace Xpetra

#define XPETRA_BLOCKEDCRSMATRIX_SHORT
#endif /* XPETRA_BLOCKEDCRSMATRIX_DECL_HPP */