20 #ifndef GEOS_LINEARALGEBRA_DOFMANAGER_HPP_
21 #define GEOS_LINEARALGEBRA_DOFMANAGER_HPP_
33 class DomainPartition;
35 class ObjectManagerBase;
36 class FluxApproximationBase;
238 string const & colFieldName,
241 bool symmetric = true );
246 string const & colFieldName,
249 bool symmetric =
true );
291 string const &
getKey(
string const & fieldName )
const;
370 template<
typename CONTAINER >
374 typename CONTAINER::value_type labelStart = 0;
375 auto it = labels.begin();
376 for( FieldDescription
const & field : m_fields )
378 integer const numComp = field.numComponents;
379 localIndex const numSupp = field.numLocalDof / numComp;
380 for(
localIndex i = 0; i < numSupp; ++i, it += numComp )
382 std::iota( it, it + numComp, labelStart );
384 labelStart += numComp;
404 string const & srcFieldName,
405 string const & dstFieldName,
418 template<
typename SCALING_FACTOR_TYPE >
420 string const & srcFieldName,
421 string const & dstFieldName,
422 SCALING_FACTOR_TYPE
const & scalingFactor,
435 string const & srcFieldName,
436 string const & dstFieldName,
450 string const & srcFieldName,
451 string const & dstFieldName,
488 template<
typename MATRIX >
490 MPI_Comm
const & comm,
492 MATRIX & restrictor )
const;
506 struct FieldDescription
527 struct CouplingDescription
537 integer getFieldIndex(
string const & name )
const;
543 void computeFieldDimensions(
integer const fieldIndex );
550 void createIndexArray( FieldDescription
const & field,
551 arrayView1d< localIndex const >
const permutation );
557 void removeIndexArray( FieldDescription
const & field );
564 array1d< localIndex > computePermutation( FieldDescription & field );
572 void computePermutation( FieldDescription
const & field,
573 arrayView1d< localIndex >
const permutation );
582 void countRowLengthsOneBlock( arrayView1d< localIndex >
const & rowLengths,
586 void countRowLengthsFromStencil( arrayView1d< localIndex >
const & rowLengths,
597 void setSparsityPatternOneBlock( SparsityPatternView< globalIndex >
const & pattern,
601 void setSparsityPatternFromStencil( SparsityPatternView< globalIndex >
const & pattern,
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,
631 template<
typename FIELD_OP,
typename POLICY >
632 void fieldToVector( arrayView1d< real64 >
const & localVector,
633 string const & srcFieldName,
634 string const & dstFieldName,
642 DomainPartition * m_domain =
nullptr;
645 stdVector< FieldDescription > m_fields;
648 std::map< std::pair< integer, integer >, CouplingDescription > m_coupling;
651 bool m_reordered =
false;
The DoFManager is responsible for allocating global dofs, constructing sparsity patterns,...
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 ®ions={})
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.
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 ®ions={}, 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.
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.
@ 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.
Connector
Enumeration of geometric objects for connectivity type. Note that this enum is nearly identical to Fi...
@ 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 ®ions, 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 ®ions)
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.
Base template for ordered and unordered maps.
ArrayView< T, 1 > arrayView1d
Alias for 1D array view.
stdVector< string > string_array
A 1-dimensional array of geos::string types.
GEOS_GLOBALINDEX_TYPE globalIndex
Global index type (for indexing objects across MPI partitions).
LvArray::SparsityPattern< COL_INDEX, INDEX_TYPE, LvArray::ChaiBuffer > SparsityPattern
Alias for Sparsity pattern class.
double real64
64-bit floating point type.
GEOS_LOCALINDEX_TYPE localIndex
Local index type (for indexing objects within an MPI partition).
FieldLocation
Enum defining the possible location of a field on the mesh.
int integer
Signed integer type.
Array< T, 1 > array1d
Alias for 1D array.
internal::StdVectorWrapper< T, Allocator, USE_STD_CONTAINER_BOUNDS_CHECKING > stdVector
Describes field support on a single mesh body/level.
std::set< string > regionNames
list of the region names
bool add(FieldSupport const &input)
add the regionNames contained in input if the meshBodyName and the meshLevelName of input match the o...
string meshLevelName
name of the mesh level
string meshBodyName
name of the mesh body
Describes a selection of components from a DoF field.
CompMask mask
Mask that defines component selection.
string fieldName
Name of the DOF field in DofManager.