Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
AssemblyMap.h
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File AssemblyMap.h
4 //
5 // For more information, please see: http://www.nektar.info
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description: Assembly (e.g. local to global) base mapping routines
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 #ifndef MULTIREGIONS_ASSEMBLY_MAP_H
37 #define MULTIREGIONS_ASSEMBLY_MAP_H
38 
42 #include <vector>
45 
46 
47 namespace Nektar
48 {
49  namespace MultiRegions
50  {
51  // Forward declarations
52  class AssemblyMap;
53  class ExpList;
54  typedef boost::shared_ptr<AssemblyMap> AssemblyMapSharedPtr;
55  static AssemblyMapSharedPtr NullAssemblyMapSharedPtr;
56 
57  /// Base class for constructing local to global mapping of degrees of
58  /// freedom.
60  {
61  public:
62  /// Default constructor.
64  /// Constructor with a communicator
67  const std::string variable = "DefaultVar");
68 
69  /// Constructor for next level in multi-level static condensation.
70  MULTI_REGIONS_EXPORT AssemblyMap(AssemblyMap* oldLevelMap,
71  const BottomUpSubStructuredGraphSharedPtr& multiLevelGraph);
72  /// Destructor.
74 
75  /// Retrieves the communicator
77 
78  /// Retrieves the hash of this map
79  MULTI_REGIONS_EXPORT size_t GetHash() const;
80 
81  MULTI_REGIONS_EXPORT int GetLocalToGlobalMap(const int i) const;
82 
83  MULTI_REGIONS_EXPORT int GetGlobalToUniversalMap(const int i) const;
84 
86 
88 
90 
92 
94 
96 
99  Array<OneD, NekDouble>& global) const;
100 
102  const NekVector<NekDouble>& loc,
103  NekVector< NekDouble>& global) const;
104 
106  const Array<OneD, const NekDouble>& global,
107  Array<OneD, NekDouble>& loc) const;
108 
110  const NekVector<NekDouble>& global,
111  NekVector< NekDouble>& loc) const;
112 
114  const Array<OneD, const NekDouble> &loc,
115  Array<OneD, NekDouble> &global) const;
116 
118  const NekVector<NekDouble>& loc,
119  NekVector< NekDouble>& global) const;
120 
122  Array<OneD, NekDouble>& pGlobal) const;
123 
125  NekVector< NekDouble>& pGlobal) const;
126 
128  Array<OneD, NekDouble>& pGlobal,
129  int offset) const;
130 
131  /// Retrieve the global index of a given local boundary mode.
132  MULTI_REGIONS_EXPORT int GetLocalToGlobalBndMap(const int i) const;
133  /// Retrieve the global indices of the local boundary modes.
135 
137 
139 
140  /// Returns true if using a modal expansion requiring a change of
141  /// sign of some modes.
143 
144  /// Retrieve the sign change of a given local boundary mode.
146  /// Retrieve the sign change for all local boundary modes.
148  /// Retrieves the global index corresponding to a boundary expansion
149  /// mode.
151  /// Retrieves the global indices corresponding to the boundary
152  /// expansion modes.
155  /// Returns the modal sign associated with a given boundary
156  /// expansion mode.
158 
159  /// Returns the global index of the boundary trace giving the
160  /// index on the boundary expansion
164 
165  /// Returns the number of global Dirichlet boundary coefficients.
167  /// Returns the number of local Dirichlet boundary coefficients.
169  /// Returns the total number of global boundary coefficients.
171  /// Returns the total number of local boundary coefficients.
173  /// Returns the total number of local coefficients.
175  /// Returns the total number of global coefficients.
177  /// Retrieves if the system is singular (true) or not (false)
179 
180  ///
182  const NekVector<NekDouble>& global,
184  int offset) const;
185 
187  const NekVector<NekDouble>& global,
188  NekVector<NekDouble>& loc) const;
189 
191  const Array<OneD, const NekDouble>& global,
193  int offset) const;
194 
196  const Array<OneD, const NekDouble>& global,
197  Array<OneD,NekDouble>& loc) const;
198 
200  const NekVector<NekDouble>& loc,
201  NekVector<NekDouble>& global,
202  int offset) const;
203 
205  const NekVector<NekDouble>& loc,
206  NekVector<NekDouble>& global) const;
207 
209  const Array<OneD, const NekDouble>& loc,
210  Array<OneD,NekDouble>& global,
211  int offset) const;
212 
214  const Array<OneD, const NekDouble>& loc,
215  Array<OneD,NekDouble>& global) const;
216 
218  NekVector<NekDouble>& global, int offset) const;
219 
221  NekVector<NekDouble>& global) const;
222 
224  Array<OneD, NekDouble>& global, int offset) const;
225 
227  Array<OneD, NekDouble>& global) const;
228 
230  Array<OneD, NekDouble>& pGlobal) const;
231 
233  NekVector< NekDouble>& pGlobal) const;
234 
236  Array<OneD, NekDouble>& pGlobal,
237  int offset) const;
238 
240 
242 
244 
246 
248 
250 
252 
254 
255  MULTI_REGIONS_EXPORT void PrintStats(std::ostream &out, std::string variable) const;
256 
259 
260  MULTI_REGIONS_EXPORT boost::shared_ptr<AssemblyMap> LinearSpaceMap(const ExpList &locexp, GlobalSysSolnType solnType);
261 
262  /// Returns the bandwidth of the boundary system.
264  /// Returns the level of static condensation for this map.
266  /// Returns the number of patches in this static condensation level.
268  /// Returns the number of local boundary coefficients in each patch.
271  /// Returns the number of local interior coefficients in each patch.
274  /// Returns the local to global mapping for the next level in the
275  /// multi-level static condensation.
276  MULTI_REGIONS_EXPORT const AssemblyMapSharedPtr
278 
279  MULTI_REGIONS_EXPORT void SetNextLevelLocalToGlobalMap( AssemblyMapSharedPtr pNextLevelLocalToGlobalMap);
280 
281  /// Returns the patch map from the previous level
282  /// of the multi-level static condensation.
284  GetPatchMapFromPrevLevel(void) const;
285 
286  /// Returns true if this is the last level in the multi-level
287  /// static condensation.
288  MULTI_REGIONS_EXPORT bool AtLastLevel() const;
289  /// Returns the method of solving global systems.
295 
297  {
299  }
300 
301  protected:
302  /// Session object
304 
305  /// Communicator
307 
308  /// Hash for map
309  size_t m_hash;
310 
311  /// Number of local boundary coefficients
313  /// Total number of global boundary coefficients
315  /// Number of Local Dirichlet Boundary Coefficients
317  /// Number of Global Dirichlet Boundary Coefficients
319  /// Flag indicating if the system is singular or not
321 
322  /// Total number of local coefficients
323  /** This corresponds to the number of total number of coefficients
324  * - For CG this corresponds to the total of bnd + int DOFs
325  * - For DG this corresponds to the number of bnd DOFs.
326  * This means that #m_numLocalCoeffs = #m_numLocalBndCoeffs
327  * This way, we can consider the trace-system solve as a
328  * statically condensed solve without interior DOFs. This allows
329  * us to use the same global system solver for both cases.
330  */
332 
333  /// Total number of global coefficients
334  /** This corresponds to the number of total number of coefficients
335  * - For CG this corresponds to the total of bnd + int DOFs.
336  * - For DG this corresponds to the number of bnd DOFs.
337  * This means that #m_numGlobalCoeffs = #m_numGlobalBndCoeffs
338  * This way, we can consider the trace-system solve as a
339  * statically condensed solve without interior DOFs. This allows
340  * us to use the same global system solver for both cases.
341  */
343 
344  /// Flag indicating if modes require sign reversal.
346 
347  /// Integer map of local boundary coeffs to global space
349  /// Integer sign of local boundary coeffs to global space
351  /// Integer map of bnd cond coeffs to global coefficients
353  /// Integer map of bnd cond coeffs to global coefficients
355  /// Integer map of bnd cond trace number to global trace number
357  /// Integer map of process coeffs to universal space
359  /// Integer map of unique process coeffs to universal space (signed)
361 
362  /// The solution type of the global system
364  /// The bandwith of the global bnd system
366 
367  /// Type type of preconditioner to use in iterative solver.
369 
370  /// Maximum iterations for iterative solver
372 
373  /// Tolerance for iterative solver
375 
376  /// sucessive RHS for iterative solver
378 
381 
382  /// The level of recursion in the case of multi-level static
383  /// condensation.
385  /// The number of patches (~elements) in the current level
387  /// The number of bnd dofs per patch
389  /// The number of int dofs per patch
391  /// Map from the patches of the previous level to the patches of
392  /// the current level
393 
394  /// The local to global mapping of the next level of recursion
395  AssemblyMapSharedPtr m_nextLevelLocalToGlobalMap;
396  /// Lowest static condensation level.
398 
399  /// Calculates the bandwidth of the boundary system.
401 
403  const Array<OneD, const NekDouble>& global,
404  Array<OneD,NekDouble>& loc);
405 
406  private:
407  /// Mapping information for previous level in MultiLevel Solver
409 
410  virtual int v_GetLocalToGlobalMap(const int i) const;
411 
412  virtual int v_GetGlobalToUniversalMap(const int i) const;
413 
414  virtual int v_GetGlobalToUniversalMapUnique(const int i) const;
415 
417 
419 
421 
422  virtual NekDouble v_GetLocalToGlobalSign(const int i) const;
423 
424  virtual const Array<OneD, NekDouble>& v_GetLocalToGlobalSign() const;
425 
426  virtual void v_LocalToGlobal(
427  const Array<OneD, const NekDouble>& loc,
428  Array<OneD, NekDouble>& global) const;
429 
430  virtual void v_LocalToGlobal(
431  const NekVector<NekDouble>& loc,
432  NekVector< NekDouble>& global) const;
433 
434  virtual void v_GlobalToLocal(
435  const Array<OneD, const NekDouble>& global,
436  Array<OneD, NekDouble>& loc) const;
437 
438  virtual void v_GlobalToLocal(
439  const NekVector<NekDouble>& global,
440  NekVector< NekDouble>& loc) const;
441 
442  virtual void v_Assemble(
443  const Array<OneD, const NekDouble> &loc,
444  Array<OneD, NekDouble> &global) const;
445 
446  virtual void v_Assemble(
447  const NekVector<NekDouble>& loc,
448  NekVector< NekDouble>& global) const;
449 
450  virtual void v_UniversalAssemble(
451  Array<OneD, NekDouble>& pGlobal) const;
452 
453  virtual void v_UniversalAssemble(
454  NekVector< NekDouble>& pGlobal) const;
455 
456  virtual void v_UniversalAssemble(
457  Array<OneD, NekDouble>& pGlobal,
458  int offset) const;
459 
460  virtual int v_GetFullSystemBandWidth() const;
461 
462  virtual int v_GetNumNonDirVertexModes() const;
463 
464  virtual int v_GetNumNonDirEdgeModes() const;
465 
466  virtual int v_GetNumNonDirFaceModes() const;
467 
468  virtual int v_GetNumDirEdges() const;
469 
470  virtual int v_GetNumDirFaces() const;
471 
472  virtual int v_GetNumNonDirEdges() const;
473 
474  virtual int v_GetNumNonDirFaces() const;
475 
476  virtual const Array<OneD, const int>&
478 
479  /// Generate a linear space mapping from existing mapping
480  virtual boost::shared_ptr<AssemblyMap> v_LinearSpaceMap(
481  const ExpList &locexp, GlobalSysSolnType solnType);
482  };
483 
484 
485  } // end of namespace
486 } // end of namespace
487 
488 #endif //MULTIREGIONS_ASSEMBLY_MAP_H
489 
490 
const Array< OneD, const int > & GetBndCondTraceToGlobalTraceMap()
PreconditionerType GetPreconType() const
bool m_systemSingular
Flag indicating if the system is singular or not.
Definition: AssemblyMap.h:320
const Array< OneD, const int > & GetGlobalToUniversalMap()
bool m_signChange
Flag indicating if modes require sign reversal.
Definition: AssemblyMap.h:345
void PrintStats(std::ostream &out, std::string variable) const
virtual int v_GetNumNonDirFaceModes() const
const Array< OneD, NekDouble > & GetLocalToGlobalSign() const
PreconditionerType m_preconType
Type type of preconditioner to use in iterative solver.
Definition: AssemblyMap.h:368
int m_numGlobalBndCoeffs
Total number of global boundary coefficients.
Definition: AssemblyMap.h:314
LibUtilities::CommSharedPtr GetComm()
Retrieves the communicator.
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: AssemblyMap.h:306
virtual ~AssemblyMap()
Destructor.
boost::shared_ptr< AssemblyMap > AssemblyMapSharedPtr
Definition: AssemblyMap.h:53
int GetNumLocalCoeffs() const
Returns the total number of local coefficients.
virtual int v_GetNumNonDirVertexModes() const
Array< OneD, int > m_bndCondTraceToGlobalTraceMap
Integer map of bnd cond trace number to global trace number.
Definition: AssemblyMap.h:356
int GetBndSystemBandWidth() const
Returns the bandwidth of the boundary system.
virtual int v_GetNumNonDirFaces() const
size_t GetHash() const
Retrieves the hash of this map.
boost::shared_ptr< BottomUpSubStructuredGraph > BottomUpSubStructuredGraphSharedPtr
#define MULTI_REGIONS_EXPORT
bool AtLastLevel() const
Returns true if this is the last level in the multi-level static condensation.
int GetNumPatches() const
Returns the number of patches in this static condensation level.
Array< OneD, const NekDouble > GetLocalToGlobalBndSign() const
Retrieve the sign change for all local boundary modes.
NekDouble GetIterativeTolerance() const
int m_numLocalCoeffs
Total number of local coefficients.
Definition: AssemblyMap.h:331
const Array< OneD, const int > & GetExtraDirEdges()
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
PatchMapSharedPtr m_patchMapFromPrevLevel
Mapping information for previous level in MultiLevel Solver.
Definition: AssemblyMap.h:408
virtual void v_GlobalToLocal(const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
boost::shared_ptr< PatchMap > PatchMapSharedPtr
const Array< OneD, const int > & GetGlobalToUniversalBndMapUnique()
virtual void v_UniversalAssemble(Array< OneD, NekDouble > &pGlobal) const
const Array< OneD, const int > & GetLocalToGlobalBndMap()
Retrieve the global indices of the local boundary modes.
AssemblyMapSharedPtr m_nextLevelLocalToGlobalMap
Map from the patches of the previous level to the patches of the current level.
Definition: AssemblyMap.h:395
Base class for constructing local to global mapping of degrees of freedom.
Definition: AssemblyMap.h:59
int m_successiveRHS
sucessive RHS for iterative solver
Definition: AssemblyMap.h:377
size_t m_hash
Hash for map.
Definition: AssemblyMap.h:309
void UniversalAssemble(Array< OneD, NekDouble > &pGlobal) const
Base class for all multi-elemental spectral/hp expansions.
Definition: ExpList.h:101
boost::shared_ptr< Comm > CommSharedPtr
Pointer to a Communicator object.
Definition: Comm.h:53
int GetNumGlobalDirBndCoeffs() const
Returns the number of global Dirichlet boundary coefficients.
virtual boost::shared_ptr< AssemblyMap > v_LinearSpaceMap(const ExpList &locexp, GlobalSysSolnType solnType)
Generate a linear space mapping from existing mapping.
GlobalSysSolnType GetGlobalSysSolnType() const
Returns the method of solving global systems.
int m_bndSystemBandWidth
The bandwith of the global bnd system.
Definition: AssemblyMap.h:365
NekDouble m_iterativeTolerance
Tolerance for iterative solver.
Definition: AssemblyMap.h:374
Array< OneD, unsigned int > m_numLocalBndCoeffsPerPatch
The number of bnd dofs per patch.
Definition: AssemblyMap.h:388
boost::shared_ptr< AssemblyMap > LinearSpaceMap(const ExpList &locexp, GlobalSysSolnType solnType)
GlobalSysSolnType m_solnType
The solution type of the global system.
Definition: AssemblyMap.h:363
Array< OneD, NekDouble > m_bndCondCoeffsToGlobalCoeffsSign
Integer map of bnd cond coeffs to global coefficients.
Definition: AssemblyMap.h:354
int m_numGlobalDirBndCoeffs
Number of Global Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:318
NekDouble GetBndCondCoeffsToGlobalCoeffsSign(const int i)
Returns the modal sign associated with a given boundary expansion mode.
int GetStaticCondLevel() const
Returns the level of static condensation for this map.
void GlobalToLocal(const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc) const
virtual int v_GetNumDirFaces() const
virtual int v_GetFullSystemBandWidth() const
void SetNextLevelLocalToGlobalMap(AssemblyMapSharedPtr pNextLevelLocalToGlobalMap)
const PatchMapSharedPtr & GetPatchMapFromPrevLevel(void) const
Returns the patch map from the previous level of the multi-level static condensation.
double NekDouble
void LocalToGlobal(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
virtual void v_LocalToGlobal(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
Array< OneD, unsigned int > m_numLocalIntCoeffsPerPatch
The number of int dofs per patch.
Definition: AssemblyMap.h:390
int m_lowestStaticCondLevel
Lowest static condensation level.
Definition: AssemblyMap.h:397
void CalculateBndSystemBandWidth()
Calculates the bandwidth of the boundary system.
Array< OneD, int > m_localToGlobalBndMap
Integer map of local boundary coeffs to global space.
Definition: AssemblyMap.h:348
const Array< OneD, const int > & GetBndCondCoeffsToGlobalCoeffsMap()
Retrieves the global indices corresponding to the boundary expansion modes.
int m_numLocalDirBndCoeffs
Number of Local Dirichlet Boundary Coefficients.
Definition: AssemblyMap.h:316
bool GetSignChange()
Returns true if using a modal expansion requiring a change of sign of some modes. ...
const Array< OneD, const int > & GetLocalToGlobalMap()
virtual void v_Assemble(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
const Array< OneD, const int > & GetGlobalToUniversalBndMap()
virtual int v_GetNumNonDirEdgeModes() const
Array< OneD, int > m_bndCondCoeffsToGlobalCoeffsMap
Integer map of bnd cond coeffs to global coefficients.
Definition: AssemblyMap.h:352
int m_numLocalBndCoeffs
Number of local boundary coefficients.
Definition: AssemblyMap.h:312
static AssemblyMapSharedPtr NullAssemblyMapSharedPtr
Definition: AssemblyMap.h:55
int m_staticCondLevel
The level of recursion in the case of multi-level static condensation.
Definition: AssemblyMap.h:384
Array< OneD, NekDouble > m_localToGlobalBndSign
Integer sign of local boundary coeffs to global space.
Definition: AssemblyMap.h:350
int GetNumGlobalBndCoeffs() const
Returns the total number of global boundary coefficients.
AssemblyMap()
Default constructor.
Definition: AssemblyMap.cpp:79
int m_maxIterations
Maximum iterations for iterative solver.
Definition: AssemblyMap.h:371
const Array< OneD, const unsigned int > & GetNumLocalBndCoeffsPerPatch()
Returns the number of local boundary coefficients in each patch.
void GlobalToLocalBndWithoutSign(const Array< OneD, const NekDouble > &global, Array< OneD, NekDouble > &loc)
const Array< OneD, const int > & GetGlobalToUniversalMapUnique()
int GetNumLocalDirBndCoeffs() const
Returns the number of local Dirichlet boundary coefficients.
LibUtilities::SessionReaderSharedPtr m_session
Session object.
Definition: AssemblyMap.h:303
void LocalBndToGlobal(const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, int offset) const
void AssembleBnd(const NekVector< NekDouble > &loc, NekVector< NekDouble > &global, int offset) const
virtual const Array< OneD, NekDouble > & v_GetLocalToGlobalSign() const
virtual const Array< OneD, const int > & v_GetExtraDirEdges()
int GetNumGlobalCoeffs() const
Returns the total number of global coefficients.
int GetNumLocalBndCoeffs() const
Returns the total number of local boundary coefficients.
Array< OneD, int > m_globalToUniversalBndMap
Integer map of process coeffs to universal space.
Definition: AssemblyMap.h:358
Array< OneD, int > m_globalToUniversalBndMapUnique
Integer map of unique process coeffs to universal space (signed)
Definition: AssemblyMap.h:360
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMapUnique()
void UniversalAssembleBnd(Array< OneD, NekDouble > &pGlobal) const
int m_numGlobalCoeffs
Total number of global coefficients.
Definition: AssemblyMap.h:342
void GlobalToLocalBnd(const NekVector< NekDouble > &global, NekVector< NekDouble > &loc, int offset) const
const Array< OneD, const unsigned int > & GetNumLocalIntCoeffsPerPatch()
Returns the number of local interior coefficients in each patch.
virtual int v_GetNumNonDirEdges() const
virtual const Array< OneD, const int > & v_GetGlobalToUniversalMap()
const AssemblyMapSharedPtr GetNextLevelLocalToGlobalMap() const
Returns the local to global mapping for the next level in the multi-level static condensation.
void Assemble(const Array< OneD, const NekDouble > &loc, Array< OneD, NekDouble > &global) const
virtual const Array< OneD, const int > & v_GetLocalToGlobalMap()
virtual int v_GetNumDirEdges() const
int m_numPatches
The number of patches (~elements) in the current level.
Definition: AssemblyMap.h:386
bool GetSingularSystem() const
Retrieves if the system is singular (true) or not (false)