46#ifndef MUELU_GENERALGEOMETRICPFACTORY_DECL_HPP
47#define MUELU_GENERALGEOMETRICPFACTORY_DECL_HPP
49#include <Teuchos_SerialDenseVector.hpp>
51#include <Xpetra_MultiVector.hpp>
52#include <Xpetra_Matrix_fwd.hpp>
55#include "MueLu_PFactory.hpp"
62 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
120#undef MUELU_GENERALGEOMETRICPFACTORY_SHORT
205 RCP<const Map> fineCoordsMap, RCP<GeometricData> myGeometry,
206 RCP<NodesIDs> ghostedCoarseNodes,
207 Array<Array<GO> >& lCoarseNodesGIDs)
const;
210 RCP<const Map> fineCoordsMap, RCP<GeometricData> myGeometry,
211 RCP<NodesIDs> ghostedCoarseNodes,
212 Array<Array<GO> >& lCoarseNodesGIDs)
const;
215 const RCP<Xpetra::MultiVector<
typename Teuchos::ScalarTraits<Scalar>::coordinateType,LO,GO,NO> >& fCoords,
216 const LO nnzP,
const LO dofsPerNode,
217 RCP<const Map>& stridedDomainMapP,
218 RCP<Matrix> & Amat, RCP<Matrix>& P,
219 RCP<Xpetra::MultiVector<
typename Teuchos::ScalarTraits<Scalar>::coordinateType,LO,GO,NO> >& cCoords,
220 RCP<NodesIDs> ghostedCoarseNodes, Array<Array<GO> > coarseNodesGIDs,
221 int interpolationOrder)
const;
223 void ComputeStencil(
const LO numDimension,
const Array<GO> currentNodeIndices,
224 const Array<GO> coarseNodeIndices,
const LO rate[3],
225 const Array<Array<
typename Teuchos::ScalarTraits<Scalar>::coordinateType> > coord,
const int interpolationOrder,
226 std::vector<double>& stencil)
const;
229 const Array<GO> currentNodeIndices,
230 const Array<GO> coarseNodeIndices,
231 const LO rate[3], std::vector<double>& stencil)
const;
234 std::vector<double>& stencil)
const;
236 const Teuchos::SerialDenseVector<LO,double> parameters,
237 double functions[4][8])
const;
240 const typename Teuchos::Array<LocalOrdinal>::iterator& first1,
241 const typename Teuchos::Array<LocalOrdinal>::iterator& last1,
242 const typename Teuchos::Array<LocalOrdinal>::iterator& first2,
243 const typename Teuchos::Array<LocalOrdinal>::iterator& last2)
const;
246 const typename Teuchos::Array<LocalOrdinal>::iterator& first1,
247 const typename Teuchos::Array<LocalOrdinal>::iterator& last1,
248 const typename Teuchos::Array<LocalOrdinal>::iterator& first2,
249 const typename Teuchos::Array<LocalOrdinal>::iterator& last2)
const;
252 const Array<LO> coarseNodeFineIndices,
253 const RCP<GeometricData> myGeo,
const LO myRankIndex,
const LO pi,
254 const LO pj,
const LO pk,
255 const typename std::vector<std::vector<GO> >::iterator blockStart,
256 const typename std::vector<std::vector<GO> >::iterator blockEnd,
257 GO& myGID, LO& myPID, LO& myLID)
const;
263#define MUELU_GENERALGEOMETRICPFACTORY_SHORT
MueLu::DefaultLocalOrdinal LocalOrdinal
MueLu::DefaultScalar Scalar
MueLu::DefaultGlobalOrdinal GlobalOrdinal
Prolongator factory performing geometric coarsening.
void MakeGeneralGeometricP(RCP< GeometricData > myGeo, const RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::coordinateType, LO, GO, NO > > &fCoords, const LO nnzP, const LO dofsPerNode, RCP< const Map > &stridedDomainMapP, RCP< Matrix > &Amat, RCP< Matrix > &P, RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::coordinateType, LO, GO, NO > > &cCoords, RCP< NodesIDs > ghostedCoarseNodes, Array< Array< GO > > coarseNodesGIDs, int interpolationOrder) const
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
virtual ~GeneralGeometricPFactory()
Destructor.
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
void sh_sort_permute(const typename Teuchos::Array< LocalOrdinal >::iterator &first1, const typename Teuchos::Array< LocalOrdinal >::iterator &last1, const typename Teuchos::Array< LocalOrdinal >::iterator &first2, const typename Teuchos::Array< LocalOrdinal >::iterator &last2) const
void ComputeLinearInterpolationStencil(const LO numDimension, const Array< Array< typename Teuchos::ScalarTraits< Scalar >::coordinateType > > coord, std::vector< double > &stencil) const
void BuildP(Level &fineLevel, Level &coarseLevel) const
Abstract Build method.
void ComputeStencil(const LO numDimension, const Array< GO > currentNodeIndices, const Array< GO > coarseNodeIndices, const LO rate[3], const Array< Array< typename Teuchos::ScalarTraits< Scalar >::coordinateType > > coord, const int interpolationOrder, std::vector< double > &stencil) const
void GetCoarsePoints(const int interpolationOrder, const LO blkSize, RCP< const Map > fineCoordsMap, RCP< GeometricData > myGeometry, RCP< NodesIDs > ghostedCoarseNodes, Array< Array< GO > > &lCoarseNodesGIDs) const
void GetInterpolationFunctions(const LO numDimension, const Teuchos::SerialDenseVector< LO, double > parameters, double functions[4][8]) const
void GetGIDLocalLexicographic(const GO i, const GO j, const GO k, const Array< LO > coarseNodeFineIndices, const RCP< GeometricData > myGeo, const LO myRankIndex, const LO pi, const LO pj, const LO pk, const typename std::vector< std::vector< GO > >::iterator blockStart, const typename std::vector< std::vector< GO > >::iterator blockEnd, GO &myGID, LO &myPID, LO &myLID) const
void MeshLayoutInterface(const int interpolationOrder, const LO blkSize, RCP< const Map > fineCoordsMap, RCP< GeometricData > myGeometry, RCP< NodesIDs > ghostedCoarseNodes, Array< Array< GO > > &lCoarseNodesGIDs) const
GeneralGeometricPFactory()
Constructor.
void sh_sort2(const typename Teuchos::Array< LocalOrdinal >::iterator &first1, const typename Teuchos::Array< LocalOrdinal >::iterator &last1, const typename Teuchos::Array< LocalOrdinal >::iterator &first2, const typename Teuchos::Array< LocalOrdinal >::iterator &last2) const
void ComputeConstantInterpolationStencil(const LO numDimension, const Array< GO > currentNodeIndices, const Array< GO > coarseNodeIndices, const LO rate[3], std::vector< double > &stencil) const
Class that holds all level-specific information.
Factory that provides an interface for a concrete implementation of a prolongation operator.
Namespace for MueLu classes and methods.
KokkosClassic::DefaultNode::DefaultNodeType DefaultNode
Tpetra::Details::DefaultTypes::scalar_type DefaultScalar
Array< GO > gCoarseNodesPerDir
Array< GO > startGhostedCoarseNode
Array< LO > lFineNodesPerDir
Array< LO > lCoarseNodesPerDir
std::vector< std::vector< GO > > meshData
Array< LO > ghostedCoarseNodesPerDir
Array< GO > gFineNodesPerDir
std::vector< GO > colInds