GEOS
DofManager.hpp
Go to the documentation of this file.
1 /*
2  * ------------------------------------------------------------------------------------------------------------
3  * SPDX-License-Identifier: LGPL-2.1-only
4  *
5  * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC
6  * Copyright (c) 2018-2024 TotalEnergies
7  * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University
8  * Copyright (c) 2023-2024 Chevron
9  * Copyright (c) 2019- GEOS/GEOSX Contributors
10  * All rights reserved
11  *
12  * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details.
13  * ------------------------------------------------------------------------------------------------------------
14  */
15 
20 #ifndef GEOS_LINEARALGEBRA_DOFMANAGER_HPP_
21 #define GEOS_LINEARALGEBRA_DOFMANAGER_HPP_
22 
23 #include "common/DataTypes.hpp"
24 #include "common/Span.hpp"
27 
28 #include <numeric>
29 
30 namespace geos
31 {
32 
33 class DomainPartition;
34 class MeshLevel;
35 class ObjectManagerBase;
36 class FluxApproximationBase;
37 
45 {
46 public:
47 
49  static int constexpr maxNumComp = 32;
50 
53 
59  struct SubComponent
60  {
61  string fieldName;
63  };
64 
68  struct FieldSupport
69  {
70 public:
72  string meshBodyName;
74  string meshLevelName;
76  std::set< string > regionNames;
77 
86  bool add( FieldSupport const & input )
87  {
88  bool added = false;
89  if( meshBodyName == input.meshBodyName && meshLevelName == input.meshLevelName )
90  {
91  regionNames.insert( input.regionNames.begin(), input.regionNames.end() );
92  added = true;
93  }
94 
95  return added;
96  }
97  };
98 
104  enum class Connector
105  {
106  Elem,
107  Face,
108  Edge,
109  Node,
110  None,
111  Stencil
112  };
113 
118  {
119  None,
121  };
122 
128  explicit DofManager( string name );
129 
133  DofManager( DofManager const & ) = delete;
134 
138  DofManager( DofManager && ) = default;
139 
144  DofManager & operator=( DofManager const & ) = delete;
145 
150  DofManager & operator=( DofManager && ) = default;
151 
155  void clear();
156 
165  void setDomain( DomainPartition & domain );
166 
175  void addField( string const & fieldName,
177  integer components,
178  stdVector< FieldSupport > const & regions = {} );
179 
185  void addField( string const & fieldName,
187  integer components,
188  map< std::pair< string, string >, string_array > const & regions );
189 
195  void setLocalReorderingType( string const & fieldName,
196  LocalReorderingType const reorderingType );
197 
203  void disableGlobalCouplingForEquation( string const & fieldName,
204  integer const c );
205 
211  void disableGlobalCouplingForEquations( string const & fieldName,
212  arrayView1d< integer const > const components );
213 
237  void addCoupling( string const & rowFieldName,
238  string const & colFieldName,
239  Connector connectivity,
240  stdVector< FieldSupport > const & regions = {},
241  bool symmetric = true );
245  void addCoupling( string const & rowFieldName,
246  string const & colFieldName,
247  Connector connectivity,
248  map< std::pair< string, string >, string_array > const & regions,
249  bool symmetric = true );
250 
259  void addCoupling( string const & fieldName,
260  FluxApproximationBase const & stencils );
261 
278 
284  bool fieldExists( string const & name ) const;
285 
291  string const & getKey( string const & fieldName ) const;
292 
297  globalIndex numGlobalDofs( string const & fieldName ) const;
298 
303 
308  localIndex numLocalDofs( string const & fieldName ) const;
309 
314 
319  globalIndex rankOffset( string const & fieldName ) const;
320 
325 
330  integer numComponents( string const & fieldName ) const;
331 
336 
342  FieldLocation location( string const & fieldName ) const;
343 
348  globalIndex globalOffset( string const & fieldName ) const;
349 
354  Span< FieldSupport const > support( string const & fieldName ) const;
355 
361 
370  template< typename CONTAINER >
371  void getLocalDofComponentLabels( CONTAINER & labels ) const
372  {
373  labels.resize( numLocalDofs() );
374  typename CONTAINER::value_type labelStart = 0;
375  auto it = labels.begin();
376  for( FieldDescription const & field : m_fields )
377  {
378  integer const numComp = field.numComponents;
379  localIndex const numSupp = field.numLocalDof / numComp;
380  for( localIndex i = 0; i < numSupp; ++i, it += numComp )
381  {
382  std::iota( it, it + numComp, labelStart );
383  }
384  labelStart += numComp;
385  }
386  }
387 
393 
404  string const & srcFieldName,
405  string const & dstFieldName,
406  real64 scalingFactor,
407  CompMask mask = CompMask( maxNumComp, true ) ) const;
408 
418  template< typename SCALING_FACTOR_TYPE >
420  string const & srcFieldName,
421  string const & dstFieldName,
422  SCALING_FACTOR_TYPE const & scalingFactor,
423  CompMask mask = CompMask( maxNumComp, true ) ) const;
424 
434  void copyFieldToVector( arrayView1d< real64 > const & localVector,
435  string const & srcFieldName,
436  string const & dstFieldName,
437  real64 scalingFactor,
438  CompMask mask = CompMask( maxNumComp, true ) ) const;
439 
449  void addFieldToVector( arrayView1d< real64 > const & localVector,
450  string const & srcFieldName,
451  string const & dstFieldName,
452  real64 scalingFactor,
453  CompMask mask = CompMask( maxNumComp, true ) ) const;
454 
465  filterDofs( stdVector< SubComponent > const & excluded ) const;
466 
473  void setupFrom( DofManager const & source,
474  stdVector< SubComponent > const & selection );
475 
488  template< typename MATRIX >
489  void makeRestrictor( stdVector< SubComponent > const & selection,
490  MPI_Comm const & comm,
491  bool transpose,
492  MATRIX & restrictor ) const;
493 
499  void printFieldInfo( std::ostream & os = std::cout ) const;
500 
501 private:
502 
506  struct FieldDescription
507  {
508  string name;
509  string key;
510  string docstring;
513  integer numComponents = 1;
514  CompMask globallyCoupledComponents;
516  localIndex numLocalDof = 0;
517  globalIndex numGlobalDof = 0;
518  globalIndex blockOffset = 0;
519  globalIndex rankOffset = 0;
522  };
523 
527  struct CouplingDescription
528  {
529  Connector connector = Connector::None;
531  FluxApproximationBase const * stencils = nullptr;
532  };
533 
537  integer getFieldIndex( string const & name ) const;
538 
543  void computeFieldDimensions( integer const fieldIndex );
544 
550  void createIndexArray( FieldDescription const & field,
551  arrayView1d< localIndex const > const permutation );
552 
557  void removeIndexArray( FieldDescription const & field );
558 
564  array1d< localIndex > computePermutation( FieldDescription & field );
565 
572  void computePermutation( FieldDescription const & field,
573  arrayView1d< localIndex > const permutation );
574 
575 
582  void countRowLengthsOneBlock( arrayView1d< localIndex > const & rowLengths,
583  integer rowFieldIndex,
584  integer colFieldIndex ) const;
585 
586  void countRowLengthsFromStencil( arrayView1d< localIndex > const & rowLengths,
587  integer fieldIndex ) const;
588 
597  void setSparsityPatternOneBlock( SparsityPatternView< globalIndex > const & pattern,
598  integer rowFieldIndex,
599  integer colFieldIndex ) const;
600 
601  void setSparsityPatternFromStencil( SparsityPatternView< globalIndex > const & pattern,
602  integer fieldIndex ) const;
603 
614  template< typename FIELD_OP, typename POLICY, typename SCALING_FACTOR_TYPE >
615  void vectorToField( arrayView1d< real64 const > const & localVector,
616  string const & srcFieldName,
617  string const & dstFieldName,
618  SCALING_FACTOR_TYPE const & scalingFactor,
619  CompMask mask ) const;
620 
631  template< typename FIELD_OP, typename POLICY >
632  void fieldToVector( arrayView1d< real64 > const & localVector,
633  string const & srcFieldName,
634  string const & dstFieldName,
635  real64 scalingFactor,
636  CompMask mask ) const;
637 
639  string m_name;
640 
642  DomainPartition * m_domain = nullptr;
643 
645  stdVector< FieldDescription > m_fields;
646 
648  std::map< std::pair< integer, integer >, CouplingDescription > m_coupling;
649 
651  bool m_reordered = false;
652 };
653 
654 } /* namespace geos */
655 
656 #endif /*GEOS_LINEARALGEBRA_DOFMANAGER_HPP_*/
The DoFManager is responsible for allocating global dofs, constructing sparsity patterns,...
Definition: DofManager.hpp:45
void setSparsityPattern(SparsityPattern< globalIndex > &pattern) const
Populate sparsity pattern of the entire system matrix.
void setupFrom(DofManager const &source, stdVector< SubComponent > const &selection)
Populate this manager from another using a sub-selection of fields/components.
void addField(string const &fieldName, FieldLocation location, integer components, stdVector< FieldSupport > const &regions={})
Add a new field and enumerate its degrees-of-freedom.
void copyFieldToVector(arrayView1d< real64 > const &localVector, string const &srcFieldName, string const &dstFieldName, real64 scalingFactor, CompMask mask=CompMask(maxNumComp, true)) const
Copy values from simulation data arrays to vectors.
static constexpr int maxNumComp
Maximum number of components in a field.
Definition: DofManager.hpp:49
array1d< integer > numComponentsPerField() const
Return an array of number of components per field, sorted by field registration order.
void addVectorToField(arrayView1d< real64 const > const &localVector, string const &srcFieldName, string const &dstFieldName, SCALING_FACTOR_TYPE const &scalingFactor, CompMask mask=CompMask(maxNumComp, true)) const
Add values from LA vectors to simulation data arrays.
void copyVectorToField(arrayView1d< real64 const > const &localVector, string const &srcFieldName, string const &dstFieldName, real64 scalingFactor, CompMask mask=CompMask(maxNumComp, true)) const
Copy values from LA vectors to simulation data arrays.
DofManager(DofManager const &)=delete
Deleted copy constructor.
globalIndex rankOffset(string const &fieldName) const
void disableGlobalCouplingForEquations(string const &fieldName, arrayView1d< integer const > const components)
Disable the global coupling for a set of equations.
integer numComponents() const
localIndex numLocalDofs() const
void setLocalReorderingType(string const &fieldName, LocalReorderingType const reorderingType)
Set the local reodering of the dof numbers.
void addCoupling(string const &rowFieldName, string const &colFieldName, Connector connectivity, stdVector< FieldSupport > const &regions={}, bool symmetric=true)
Add coupling between two fields.
globalIndex numGlobalDofs(string const &fieldName) const
void addCoupling(string const &fieldName, FluxApproximationBase const &stencils)
Special interface for self-connectivity through a stencil.
DofManager(string name)
Constructor.
void getLocalDofComponentLabels(CONTAINER &labels) const
Fill a container with unique dof labels for each local dof.
Definition: DofManager.hpp:371
void clear()
Remove all fields and couplings and re-enable addition of new fields.
void reorderByRank()
Finish populating fields and apply appropriate dof renumbering.
bool fieldExists(string const &name) const
Check if string key is already being used.
void makeRestrictor(stdVector< SubComponent > const &selection, MPI_Comm const &comm, bool transpose, MATRIX &restrictor) const
Create a matrix that restricts vectors and matrices to a subset of DOFs.
globalIndex rankOffset() const
globalIndex globalOffset(string const &fieldName) const
FieldLocation location(string const &fieldName) const
Get the support location type of the field.
DofManager & operator=(DofManager &&)=default
Defaulted move assignment.
LocalReorderingType
Indicates the type of (local to a rank) reordering applied to a given field.
Definition: DofManager.hpp:118
@ ReverseCutHillMcKee
Use reverve CutHill-McKee reordering algorithm.
@ None
Do not reorder the variables.
globalIndex numGlobalDofs() const
localIndex numLocalDofs(string const &fieldName) const
void disableGlobalCouplingForEquation(string const &fieldName, integer const c)
Disable the global coupling for a given equation.
ComponentMask< maxNumComp > CompMask
Type of component mask used by DofManager.
Definition: DofManager.hpp:52
Connector
Enumeration of geometric objects for connectivity type. Note that this enum is nearly identical to Fi...
Definition: DofManager.hpp:105
@ None
there is no connectivity (self connected field, like a lumped mass matrix)
@ Node
connectivity is node (like in finite volumes MPFA)
@ Face
connectivity is face (like in finite volumes TPFA)
@ Elem
connectivity is element (like in finite elements)
@ Stencil
connectivity is through a (set of) user-provided stencil(s)
@ Edge
connectivity is edge (like fracture element connectors)
DofManager & operator=(DofManager const &)=delete
Deleted copy assignment.
Span< FieldSupport const > support(string const &fieldName) const
void setDomain(DomainPartition &domain)
Assign a domain.
stdVector< SubComponent > filterDofs(stdVector< SubComponent > const &excluded) const
Create a dof selection by filtering out excluded components.
void addCoupling(string const &rowFieldName, string const &colFieldName, Connector connectivity, map< std::pair< string, string >, string_array > const &regions, bool symmetric=true)
Add coupling between two fields.
void printFieldInfo(std::ostream &os=std::cout) const
Print the summary of declared fields and coupling.
void addFieldToVector(arrayView1d< real64 > const &localVector, string const &srcFieldName, string const &dstFieldName, real64 scalingFactor, CompMask mask=CompMask(maxNumComp, true)) const
Add values from a simulation data array to a DOF vector.
DofManager(DofManager &&)=default
Move constructor.
string const & getKey(string const &fieldName) const
Return the key used to record the field in the DofManager.
integer numComponents(string const &fieldName) const
void addField(string const &fieldName, FieldLocation location, integer components, map< std::pair< string, string >, string_array > const &regions)
Add a new field and enumerate its degrees-of-freedom.
Partition of the decomposed physical domain. It also manages the connexion information to its neighbo...
Lightweight non-owning wrapper over a contiguous range of elements.
Definition: Span.hpp:42
Base template for ordered and unordered maps.
ArrayView< T, 1 > arrayView1d
Alias for 1D array view.
Definition: DataTypes.hpp:179
stdVector< string > string_array
A 1-dimensional array of geos::string types.
Definition: DataTypes.hpp:361
GEOS_GLOBALINDEX_TYPE globalIndex
Global index type (for indexing objects across MPI partitions).
Definition: DataTypes.hpp:87
LvArray::SparsityPattern< COL_INDEX, INDEX_TYPE, LvArray::ChaiBuffer > SparsityPattern
Alias for Sparsity pattern class.
Definition: DataTypes.hpp:297
double real64
64-bit floating point type.
Definition: DataTypes.hpp:98
GEOS_LOCALINDEX_TYPE localIndex
Local index type (for indexing objects within an MPI partition).
Definition: DataTypes.hpp:84
FieldLocation
Enum defining the possible location of a field on the mesh.
int integer
Signed integer type.
Definition: DataTypes.hpp:81
Array< T, 1 > array1d
Alias for 1D array.
Definition: DataTypes.hpp:175
internal::StdVectorWrapper< T, Allocator, USE_STD_CONTAINER_BOUNDS_CHECKING > stdVector
Describes field support on a single mesh body/level.
Definition: DofManager.hpp:69
std::set< string > regionNames
list of the region names
Definition: DofManager.hpp:76
bool add(FieldSupport const &input)
add the regionNames contained in input if the meshBodyName and the meshLevelName of input match the o...
Definition: DofManager.hpp:86
string meshLevelName
name of the mesh level
Definition: DofManager.hpp:74
string meshBodyName
name of the mesh body
Definition: DofManager.hpp:72
Describes a selection of components from a DoF field.
Definition: DofManager.hpp:60
CompMask mask
Mask that defines component selection.
Definition: DofManager.hpp:62
string fieldName
Name of the DOF field in DofManager.
Definition: DofManager.hpp:61