DofManager.hpp

Go to the documentation of this file.

/*
 * ------------------------------------------------------------------------------------------------------------
 * SPDX-License-Identifier: LGPL-2.1-only
 *
 * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC
 * Copyright (c) 2018-2024 Total, S.A
 * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University
 * Copyright (c) 2023-2024 Chevron
 * Copyright (c) 2019-     GEOS/GEOSX Contributors
 * All rights reserved
 *
 * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details.
 * ------------------------------------------------------------------------------------------------------------
 */
 
#ifndef GEOS_LINEARALGEBRA_DOFMANAGER_HPP_
#define GEOS_LINEARALGEBRA_DOFMANAGER_HPP_
 
#include "common/DataTypes.hpp"
#include "linearAlgebra/utilities/ComponentMask.hpp"
#include "mesh/FieldIdentifiers.hpp"
 
#include <numeric>
 
namespace geos
{
 
class DomainPartition;
class MeshLevel;
class ObjectManagerBase;
class FluxApproximationBase;
 
class DofManager
{
public:
 
  static int constexpr MAX_COMP = 32;
 
  using CompMask = ComponentMask< MAX_COMP >;
 
  struct SubComponent
  {
    string fieldName;  
    CompMask mask;     
  };
 
  struct FieldSupport
  {
public:
    string meshBodyName;
    string meshLevelName;
    std::set< string > regionNames;
 
    bool add( FieldSupport const & input )
    {
      bool added = false;
      if( meshBodyName  == input.meshBodyName && meshLevelName == input.meshLevelName )
      {
        regionNames.insert( input.regionNames.begin(), input.regionNames.end() );
        added = true;
      }
 
      return added;
    }
  };
 
  enum class Connector
  {
    Elem, 
    Face, 
    Edge, 
    Node, 
    None,  
    Stencil 
  };
 
  enum class LocalReorderingType
  {
    None,    
    ReverseCutHillMcKee, 
  };
 
  explicit DofManager( string name );
 
  DofManager( DofManager const & ) = delete;
 
  DofManager( DofManager && ) = default;
 
  DofManager & operator=( DofManager const & ) = delete;
 
  DofManager & operator=( DofManager && ) = default;
 
  void clear();
 
  void setDomain( DomainPartition & domain );
 
  void addField( string const & fieldName,
                 FieldLocation location,
                 integer components,
                 std::vector< FieldSupport > const & regions = {} );
 
  void addField( string const & fieldName,
                 FieldLocation location,
                 integer components,
                 map< std::pair< string, string >, array1d< string > > const & regions );
 
  void setLocalReorderingType( string const & fieldName,
                               LocalReorderingType const reorderingType );
 
  void disableGlobalCouplingForEquation( string const & fieldName,
                                         integer const c );
 
  void disableGlobalCouplingForEquations( string const & fieldName,
                                          arrayView1d< integer const > const components );
 
  void addCoupling( string const & rowFieldName,
                    string const & colFieldName,
                    Connector connectivity,
                    std::vector< FieldSupport > const & regions = {},
                    bool symmetric = true );
  void addCoupling( string const & rowFieldName,
                    string const & colFieldName,
                    Connector connectivity,
                    map< std::pair< string, string >, array1d< string > > const & regions,
                    bool symmetric = true );
 
  void addCoupling( string const & fieldName,
                    FluxApproximationBase const & stencils );
 
  void reorderByRank();
 
  bool fieldExists( string const & name ) const;
 
  string const & getKey( string const & fieldName ) const;
 
  globalIndex numGlobalDofs( string const & fieldName ) const;
 
  globalIndex numGlobalDofs() const;
 
  localIndex numLocalDofs( string const & fieldName ) const;
 
  localIndex numLocalDofs() const;
 
  globalIndex rankOffset( string const & fieldName ) const;
 
  globalIndex rankOffset() const;
 
  integer numComponents( string const & fieldName = "" ) const;
 
  integer numComponents() const;
 
  FieldLocation location( string const & fieldName ) const;
 
  globalIndex globalOffset( string const & fieldName ) const;
 
  array1d< integer > numComponentsPerField() const;
 
  template< typename CONTAINER >
  void getLocalDofComponentLabels( CONTAINER & labels ) const
  {
    labels.resize( numLocalDofs() );
    typename CONTAINER::value_type labelStart = 0;
    auto it = labels.begin();
    for( FieldDescription const & field : m_fields )
    {
      localIndex const numComp = field.numComponents;
      localIndex const numSupp = field.numLocalDof / numComp;
      for( localIndex i = 0; i < numSupp; ++i, it += numComp )
      {
        std::iota( it, it + numComp, labelStart );
      }
      labelStart += numComp;
    }
  }
 
  void setSparsityPattern( SparsityPattern< globalIndex > & pattern ) const;
 
  void copyVectorToField( arrayView1d< real64 const > const & localVector,
                          string const & srcFieldName,
                          string const & dstFieldName,
                          real64 scalingFactor,
                          CompMask mask = CompMask( MAX_COMP, true ) ) const;
 
  template< typename SCALING_FACTOR_TYPE >
  void addVectorToField( arrayView1d< real64 const > const & localVector,
                         string const & srcFieldName,
                         string const & dstFieldName,
                         SCALING_FACTOR_TYPE const & scalingFactor,
                         CompMask mask = CompMask( MAX_COMP, true ) ) const;
 
  void copyFieldToVector( arrayView1d< real64 > const & localVector,
                          string const & srcFieldName,
                          string const & dstFieldName,
                          real64 scalingFactor,
                          CompMask mask = CompMask( MAX_COMP, true ) ) const;
 
  void addFieldToVector( arrayView1d< real64 > const & localVector,
                         string const & srcFieldName,
                         string const & dstFieldName,
                         real64 scalingFactor,
                         CompMask mask = CompMask( MAX_COMP, true ) ) const;
 
  std::vector< SubComponent >
  filterDofs( std::vector< SubComponent > const & excluded ) const;
 
  void setupFrom( DofManager const & source,
                  std::vector< SubComponent > const & selection );
 
  template< typename MATRIX >
  void makeRestrictor( std::vector< SubComponent > const & selection,
                       MPI_Comm const & comm,
                       bool transpose,
                       MATRIX & restrictor ) const;
 
  void printFieldInfo( std::ostream & os = std::cout ) const;
 
private:
 
  struct FieldDescription
  {
    string name;                   
    string key;                    
    string docstring;              
    std::vector< FieldSupport > support;
    FieldLocation location;             
    integer numComponents = 1;     
    CompMask globallyCoupledComponents; 
    localIndex numLocalDof = 0;    
    globalIndex numGlobalDof = 0;  
    globalIndex blockOffset = 0;   
    globalIndex rankOffset = 0;    
    globalIndex globalOffset = 0;  
    LocalReorderingType reorderingType = LocalReorderingType::None; 
  };
 
  struct CouplingDescription
  {
    Connector connector = Connector::None;  
    std::vector< FieldSupport > support; 
    FluxApproximationBase const * stencils = nullptr; 
  };
 
  localIndex getFieldIndex( string const & name ) const;
 
  void computeFieldDimensions( localIndex fieldIndex );
 
  void createIndexArray( FieldDescription const & field,
                         arrayView1d< localIndex const > const permutation );
 
  void removeIndexArray( FieldDescription const & field );
 
  array1d< localIndex > computePermutation( FieldDescription & field );
 
  void computePermutation( FieldDescription const & field,
                           arrayView1d< localIndex > const permutation );
 
 
  void countRowLengthsOneBlock( arrayView1d< localIndex > const & rowLengths,
                                localIndex rowFieldIndex,
                                localIndex colFieldIndex ) const;
 
  void countRowLengthsFromStencil( arrayView1d< localIndex > const & rowLengths,
                                   localIndex fieldIndex ) const;
 
  void setSparsityPatternOneBlock( SparsityPatternView< globalIndex > const & pattern,
                                   localIndex rowFieldIndex,
                                   localIndex colFieldIndex ) const;
 
  void setSparsityPatternFromStencil( SparsityPatternView< globalIndex > const & pattern,
                                      localIndex fieldIndex ) const;
 
  template< int DIMS_PER_DOF >
  void setFiniteElementSparsityPattern( SparsityPattern< globalIndex > & pattern,
                                        localIndex fieldIndex ) const;
 
  template< typename FIELD_OP, typename POLICY, typename SCALING_FACTOR_TYPE >
  void vectorToField( arrayView1d< real64 const > const & localVector,
                      string const & srcFieldName,
                      string const & dstFieldName,
                      SCALING_FACTOR_TYPE const & scalingFactor,
                      CompMask mask ) const;
 
  template< typename FIELD_OP, typename POLICY >
  void fieldToVector( arrayView1d< real64 > const & localVector,
                      string const & srcFieldName,
                      string const & dstFieldName,
                      real64 scalingFactor,
                      CompMask mask ) const;
 
  string m_name;
 
  DomainPartition * m_domain = nullptr;
 
  std::vector< FieldDescription > m_fields;
 
  std::map< std::pair< localIndex, localIndex >, CouplingDescription > m_coupling;
 
  bool m_reordered = false;
};
 
} /* namespace geos */
 
#endif /*GEOS_LINEARALGEBRA_DOFMANAGER_HPP_*/