MueLu  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
MueLu_IntrepidPCoarsenFactory_decl.hpp
Go to the documentation of this file.
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // MueLu: A package for multigrid based preconditioning
6 // Copyright 2012 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact
39 // Jonathan Hu (jhu@sandia.gov)
40 // Andrey Prokopenko (aprokop@sandia.gov)
41 // Ray Tuminaro (rstumin@sandia.gov)
42 //
43 // ***********************************************************************
44 //
45 // @HEADER
46 #ifndef MUELU_IPCFACTORY_DECL_HPP
47 #define MUELU_IPCFACTORY_DECL_HPP
48 
49 #include <string>
50 #include <vector>
51 
52 #include "MueLu_ConfigDefs.hpp"
54 
55 #include "MueLu_Level_fwd.hpp"
57 #include "MueLu_PerfUtils_fwd.hpp"
58 #include "MueLu_PFactory.hpp"
61 #include "MueLu_Utilities_fwd.hpp"
62 
63 #include "Intrepid2_Basis.hpp"
64 
65 #include "Kokkos_DynRankView.hpp"
66 
67 #include <Xpetra_Import.hpp>
68 
69 namespace MueLu {
70 
111  template <class Scalar = double, class LocalOrdinal = int, class GlobalOrdinal = LocalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
113 #undef MUELU_INTREPIDPCOARSENFACTORY_SHORT
114 #include "MueLu_UseShortNames.hpp"
115 
116  public:
117  typedef Kokkos::DynRankView<LocalOrdinal,typename Node::device_type> LOFieldContainer;
118  typedef Kokkos::DynRankView<double,typename Node::device_type> SCFieldContainer;
119  typedef Intrepid2::Basis<typename Node::device_type::execution_space,double,double> Basis; // Hardwired on purpose
120 
122 
123 
128 
131 
133 
135 
137 
138 
139  void DeclareInput(Level &fineLevel, Level &coarseLevel) const;
140 
142 
144 
145 
151  void Build(Level& fineLevel, Level &coarseLevel) const;
152 
153  void BuildP(Level &fineLevel, Level &coarseLevel) const; //Build()
154 
156  private:
158 
159  // NOTE: This is hardwired to double on purpose.
161  const std::vector<bool> & hi_nodeIsOwned,
162  const SCFieldContainer & hi_DofCoords,
163  const std::vector<size_t> &lo_node_in_hi,
164  const Basis &lo_Basis,
165  const std::vector<LocalOrdinal> & hi_to_lo_map,
166  const Teuchos::RCP<const Map> & lo_colMap,
167  const Teuchos::RCP<const Map> & lo_domainMap,
168  const Teuchos::RCP<const Map> & hi_map,
169  Teuchos::RCP<Matrix>& P) const;
171 
172  // NOTE: This is hardwired to double on purpose.
174  const std::vector<bool> & hi_nodeIsOwned,
175  const SCFieldContainer & hi_DofCoords,
176  const LOFieldContainer & lo_elemToHiRepresentativeNode,
177  const Basis &lo_basis,
178  const std::vector<LocalOrdinal> & hi_to_lo_map,
179  const Teuchos::RCP<const Map> & lo_colMap,
180  const Teuchos::RCP<const Map> & lo_domainMap,
181  const Teuchos::RCP<const Map> & hi_map,
182  Teuchos::RCP<Matrix>& P) const;
183 
184 
185  }; //class IntrepidPCoarsenFactory
186 
187 
188  /* Utility functions for use with Intrepid */
189  namespace MueLuIntrepid {
190 
191  template<class Scalar,class KokkosExecutionSpace>
193 
194  template<class Scalar,class KokkosDeviceType>
195  void IntrepidGetLoNodeInHi(const Teuchos::RCP<Intrepid2::Basis<typename KokkosDeviceType::execution_space,Scalar,Scalar> > &hi_basis,
196  const Teuchos::RCP<Intrepid2::Basis<typename KokkosDeviceType::execution_space,Scalar,Scalar> > &lo_basis,
197  std::vector<size_t> & lo_node_in_hi,
198  Kokkos::DynRankView<Scalar,KokkosDeviceType> & hi_DofCoords);
199 
200  template<class LocalOrdinal, class GlobalOrdinal, class Node, class LOFieldContainer>
201  void GenerateLoNodeInHiViaGIDs(const std::vector<std::vector<size_t> > & candidates,const LOFieldContainer & hi_elemToNode,
202  RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > & hi_columnMap,
203  LOFieldContainer & lo_elemToHiRepresentativeNode);
204 
205  template <class LocalOrdinal, class LOFieldContainer>
206  void BuildLoElemToNode(const LOFieldContainer & hi_elemToNode,
207  const std::vector<bool> & hi_nodeIsOwned,
208  const std::vector<size_t> & lo_node_in_hi,
209  const Teuchos::ArrayRCP<const int> & hi_isDirichlet,
210  LOFieldContainer & lo_elemToNode,
211  std::vector<bool> & lo_nodeIsOwned,
212  std::vector<LocalOrdinal> & hi_to_lo_map,
213  int & lo_numOwnedNodes);
214 
215  template <class LocalOrdinal, class LOFieldContainer>
216  void BuildLoElemToNodeViaRepresentatives(const LOFieldContainer & hi_elemToNode,
217  const std::vector<bool> & hi_nodeIsOwned,
218  const LOFieldContainer & lo_elemToHiRepresentativeNode,
219  LOFieldContainer & lo_elemToNode,
220  std::vector<bool> & lo_nodeIsOwned,
221  std::vector<LocalOrdinal> & hi_to_lo_map,
222  int & lo_numOwnedNodes);
223 
224 
225  template <class LocalOrdinal, class GlobalOrdinal, class Node>
226  void GenerateColMapFromImport(const Xpetra::Import<LocalOrdinal,GlobalOrdinal,Node> & hi_importer,const std::vector<LocalOrdinal> &hi_to_lo_map,const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> & lo_domainMap, const size_t & lo_columnMapLength, RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > & lo_columnMap);
227 
228 
229  template<class Basis, class SCFieldContainer>
230  void GenerateRepresentativeBasisNodes(const Basis & basis, const SCFieldContainer & ReferenceNodeLocations, const double threshold, std::vector<std::vector<size_t> > & representative_node_candidates);
231 
232  // ! Given an element to (global) node map and a basis, determine one global ordinal per geometric entity (vertex, edge, face,
233  // ! interior). On exit, seeds container is of dimension (spaceDim+1), and contains a sorted vector of local ordinals
234  // ! belonging to entities of that dimension. Only locally-owned degrees of freedom (as determined by rowMap and columnMap)
235  // ! will be stored in seeds.
236  template<class Basis, class LOFieldContainer, class LocalOrdinal, class GlobalOrdinal, class Node>
237  void FindGeometricSeedOrdinals(Teuchos::RCP<Basis> basis, const LOFieldContainer &elementToNodeMap,
238  std::vector<std::vector<LocalOrdinal> > &seeds,
241 
242  }//namespace MueLuIntrepid
243 } //namespace MueLu
244 
245 #define MUELU_INTREPIDPCOARSENFACTORY_SHORT
246 #endif // MUELU_INTREPIDPCOARSENFACTORY_DECL_HPP
void GenerateLoNodeInHiViaGIDs(const std::vector< std::vector< size_t > > &candidates, const LOFieldContainer &hi_elemToNode, RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > &hi_columnMap, LOFieldContainer &lo_elemToHiRepresentativeNode)
Teuchos::RCP< Intrepid2::Basis< KokkosExecutionSpace, Scalar, Scalar > > BasisFactory(const std::string &name, int &degree)
void FindGeometricSeedOrdinals(Teuchos::RCP< Basis > basis, const LOFieldContainer &elementToNodeMap, std::vector< std::vector< LocalOrdinal > > &seeds, const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &rowMap, const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &columnMap)
void IntrepidGetLoNodeInHi(const Teuchos::RCP< Intrepid2::Basis< typename KokkosDeviceType::execution_space, Scalar, Scalar > > &hi_basis, const Teuchos::RCP< Intrepid2::Basis< typename KokkosDeviceType::execution_space, Scalar, Scalar > > &lo_basis, std::vector< size_t > &lo_node_in_hi, Kokkos::DynRankView< Scalar, KokkosDeviceType > &hi_DofCoords)
Factory for building transfer operators based on coarsening in polynomial degree, following the Intre...
void BuildP(Level &fineLevel, Level &coarseLevel) const
Abstract Build method.
IntrepidPCoarsenFactory()
Constructor. User can supply a factory for generating the tentative prolongator.
Class that holds all level-specific information.
Definition: MueLu_Level.hpp:99
Kokkos::DynRankView< double, typename Node::device_type > SCFieldContainer
void BuildLoElemToNode(const LOFieldContainer &hi_elemToNode, const std::vector< bool > &hi_nodeIsOwned, const std::vector< size_t > &lo_node_in_hi, const Teuchos::ArrayRCP< const int > &hi_isDirichlet, LOFieldContainer &lo_elemToNode, std::vector< bool > &lo_nodeIsOwned, std::vector< LocalOrdinal > &hi_to_lo_map, int &lo_numOwnedNodes)
void BuildLoElemToNodeViaRepresentatives(const LOFieldContainer &hi_elemToNode, const std::vector< bool > &hi_nodeIsOwned, const LOFieldContainer &lo_elemToHiRepresentativeNode, LOFieldContainer &lo_elemToNode, std::vector< bool > &lo_nodeIsOwned, std::vector< LocalOrdinal > &hi_to_lo_map, int &lo_numOwnedNodes)
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
void GenerateLinearCoarsening_pn_kirby_to_p1(const LOFieldContainer &hi_elemToNode, const std::vector< bool > &hi_nodeIsOwned, const SCFieldContainer &hi_DofCoords, const std::vector< size_t > &lo_node_in_hi, const Basis &lo_Basis, const std::vector< LocalOrdinal > &hi_to_lo_map, const Teuchos::RCP< const Map > &lo_colMap, const Teuchos::RCP< const Map > &lo_domainMap, const Teuchos::RCP< const Map > &hi_map, Teuchos::RCP< Matrix > &P) const
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.
Kokkos::DynRankView< LocalOrdinal, typename Node::device_type > LOFieldContainer
Intrepid2::Basis< typename Node::device_type::execution_space, double, double > Basis
void GenerateColMapFromImport(const Xpetra::Import< LocalOrdinal, GlobalOrdinal, Node > &hi_importer, const std::vector< LocalOrdinal > &hi_to_lo_map, const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &lo_domainMap, const size_t &lo_columnMapLength, RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > &lo_columnMap)
void Build(Level &fineLevel, Level &coarseLevel) const
Build method.
void GenerateRepresentativeBasisNodes(const Basis &basis, const SCFieldContainer &ReferenceNodeLocations, const double threshold, std::vector< std::vector< size_t > > &representative_node_candidates)
void GenerateLinearCoarsening_pn_kirby_to_pm(const LOFieldContainer &hi_elemToNode, const std::vector< bool > &hi_nodeIsOwned, const SCFieldContainer &hi_DofCoords, const LOFieldContainer &lo_elemToHiRepresentativeNode, const Basis &lo_basis, const std::vector< LocalOrdinal > &hi_to_lo_map, const Teuchos::RCP< const Map > &lo_colMap, const Teuchos::RCP< const Map > &lo_domainMap, const Teuchos::RCP< const Map > &hi_map, Teuchos::RCP< Matrix > &P) const
Factory that provides an interface for a concrete implementation of a prolongation operator...