DFT-EFE
 
Loading...
Searching...
No Matches
dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim > Class Template Reference

Class to get the interface between Classical and Enrichment basis. It takes as the classical basis as input. The main functionalities of the class are: More...

#include <EnrichmentClassicalInterfaceSpherical.h>

Public Member Functions

 EnrichmentClassicalInterfaceSpherical (std::shared_ptr< const FEBasisDataStorage< ValueTypeBasisData, memorySpace > > cfeBasisDataStorageOverlapMatrix, std::shared_ptr< const FEBasisDataStorage< ValueTypeBasisData, memorySpace > > cfeBasisDataStorageRhs, std::shared_ptr< const atoms::AtomSphericalDataContainer > atomSphericalDataContainer, const double atomPartitionTolerance, const std::vector< std::string > &atomSymbolVec, const std::vector< utils::Point > &atomCoordinatesVec, const std::string fieldName, std::shared_ptr< linearAlgebra::LinAlgOpContext< memorySpace > > linAlgOpContext, const utils::mpi::MPIComm &comm)
 This Constructor for augmenting the orthogonalized EFE basis with classical FE basis. and then 1. Partition the enrichment dofs and 2. Form the orthogonalized basis by L2 projecting the pristine enrichment functions on the classicla FE basis. One can write any field \(u^h_{\alpha}(x) = \sum_{i=1}^{n_h}N_i^Cu^C_{\alpha,i} +\sum_{I=1}^{N_e}N_I^{O,u}u_{\alpha,I}\) where \(\alpha\) represents the number of discretized fields, C is the classical part and O is the orthogonalized enriched basis part, \(N_I\) represents the number of enrichment dofs and \(n_h\) represents the number of classical finite element dofs. The orthogonalized basis functions can be written as \(N^{O,u}_I(x) = N^{A,u}_I(x) − N^{B,u}_I(x)\), where \(N^{A,u}_I(x)\) represents the pristine enrichment functions and \(N^{B,u}_I(x)\f is the component of enrichment function along to the classical basis and can be written in terms of classical basis functions as \)@_fakenlN^{B,u}_I(x) = \sum^{n_h}_{l=1}c_{I,l}N^{C}_{l}(x)\f. Since the \(N^{O,u}_I(x)\)'s are orthogonal to classical basis, one can solve the equation \(M^{cc}.c = d\) where M is the discretized overlap matrix in classical FE basis and \(M^{cc}_{j,l} = \integral_{\omega}N^{C}_{l}(x)N^{C}_{j}(x)dx\) and d is the RHS written as \(d_{I,j,k} = \integral_{\omega}N^{A,u}_{i,j}(x)N^{C}_{k}(x)dx\) and \(c\) is the basisInterfaceCoefficient written as a multivector. More...
 
 EnrichmentClassicalInterfaceSpherical (std::shared_ptr< const TriangulationBase > triangulation, size_type feOrder, std::shared_ptr< const atoms::AtomSphericalDataContainer > atomSphericalDataContainer, const double atomPartitionTolerance, const std::vector< std::string > &atomSymbolVec, const std::vector< utils::Point > &atomCoordinatesVec, const std::string fieldName, const utils::mpi::MPIComm &comm)
 This Constructor for augmenting the EFE basis with classical FE basis. More...
 
 ~EnrichmentClassicalInterfaceSpherical ()=default
 Destructor for the class. More...
 
std::shared_ptr< const atoms::AtomSphericalDataContainer > getAtomSphericalDataContainer () const
 Function to return AtomSphericalDataContainerObject. More...
 
std::shared_ptr< const EnrichmentIdsPartition< dim > > getEnrichmentIdsPartition () const
 Function to return EnrichmentIDsObject. More...
 
std::shared_ptr< const AtomIdsPartition< dim > > getAtomIdsPartition () const
 
std::shared_ptr< const BasisManager< ValueTypeBasisData, memorySpace > > getCFEBasisManager () const
 
std::shared_ptr< const BasisDofHandlergetCFEBasisDofHandler () const
 
const std::unordered_map< global_size_type, utils::OptimizedIndexSet< size_type > > & getClassicalComponentLocalIdsMap () const
 
const std::unordered_map< global_size_type, std::vector< ValueTypeBasisData > > & getClassicalComponentCoeffMap () const
 
std::shared_ptr< linearAlgebra::LinAlgOpContext< memorySpace > > getLinAlgOpContext () const
 
bool isOrthogonalized () const
 
std::vector< std::string > getAtomSymbolVec () const
 
std::vector< utils::PointgetAtomCoordinatesVec () const
 
std::string getFieldName () const
 
std::shared_ptr< const TriangulationBasegetTriangulation () const
 
size_type getFEOrder () const
 
global_size_type getEnrichmentId (size_type cellId, size_type enrichmentCellLocalId) const
 The localid is determined by the storage pattern of the components of basisInterfaceCoeff multivector. The multivector has a storage pattern of [LocallyownedEnrichmentIds, GhostEnrichemntIds]. More...
 
size_type getEnrichmentLocalId (global_size_type enrichmentId) const
 
size_type getEnrichmentLocalId (size_type cellId, size_type enrichmentCellLocalId) const
 

Private Attributes

std::shared_ptr< const EnrichmentIdsPartition< dim > > d_enrichmentIdsPartition
 
std::shared_ptr< const AtomIdsPartition< dim > > d_atomIdsPartition
 
std::shared_ptr< const atoms::AtomSphericalDataContainer > d_atomSphericalDataContainer
 
std::shared_ptr< const TriangulationBased_triangulation
 
bool d_isOrthogonalized
 
std::shared_ptr< const FEBasisDofHandler< ValueTypeBasisData, memorySpace, dim > > d_cfeBasisDofHandler
 
std::shared_ptr< const FEBasisManager< ValueTypeBasisData, ValueTypeBasisData, memorySpace, dim > > d_cfeBasisManager
 
const std::vector< std::string > d_atomSymbolVec
 
const std::vector< utils::Pointd_atomCoordinatesVec
 
const std::string d_fieldName
 
std::vector< std::vector< global_size_type > > d_overlappingEnrichmentIdsInCells
 
std::unordered_map< global_size_type, utils::OptimizedIndexSet< size_type > > d_enrichmentIdToClassicalLocalIdMap
 
std::unordered_map< global_size_type, std::vector< ValueTypeBasisData > > d_enrichmentIdToInterfaceCoeffMap
 
std::shared_ptr< linearAlgebra::LinAlgOpContext< memorySpace > > d_linAlgOpContext
 
size_type d_feOrder
 

Detailed Description

template<typename ValueTypeBasisData, utils::MemorySpace memorySpace, size_type dim>
class dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >

Class to get the interface between Classical and Enrichment basis. It takes as the classical basis as input. The main functionalities of the class are:

  1. Partition the enrichment degrees of freedom based on a.) periodic/non-periodic BCs b.) Orthogonalized or pristine enrichment dofs needed.
  2. Acts as input to the EFEBasisDofHandler the basic class for all EFE operations.
  3. For Orthogonalized EFE it carries an extra set of c's from Mc = d obtained from classical dofs.
  4. For periodic BC, there will be contributions to enrichment set from periodic images.

Constructor & Destructor Documentation

◆ EnrichmentClassicalInterfaceSpherical() [1/2]

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::EnrichmentClassicalInterfaceSpherical ( std::shared_ptr< const FEBasisDataStorage< ValueTypeBasisData, memorySpace > >  cfeBasisDataStorageOverlapMatrix,
std::shared_ptr< const FEBasisDataStorage< ValueTypeBasisData, memorySpace > >  cfeBasisDataStorageRhs,
std::shared_ptr< const atoms::AtomSphericalDataContainer >  atomSphericalDataContainer,
const double  atomPartitionTolerance,
const std::vector< std::string > &  atomSymbolVec,
const std::vector< utils::Point > &  atomCoordinatesVec,
const std::string  fieldName,
std::shared_ptr< linearAlgebra::LinAlgOpContext< memorySpace > >  linAlgOpContext,
const utils::mpi::MPIComm comm 
)

This Constructor for augmenting the orthogonalized EFE basis with classical FE basis. and then 1. Partition the enrichment dofs and 2. Form the orthogonalized basis by L2 projecting the pristine enrichment functions on the classicla FE basis. One can write any field \(u^h_{\alpha}(x) = \sum_{i=1}^{n_h}N_i^Cu^C_{\alpha,i} +\sum_{I=1}^{N_e}N_I^{O,u}u_{\alpha,I}\) where \(\alpha\) represents the number of discretized fields, C is the classical part and O is the orthogonalized enriched basis part, \(N_I\) represents the number of enrichment dofs and \(n_h\) represents the number of classical finite element dofs. The orthogonalized basis functions can be written as \(N^{O,u}_I(x) = N^{A,u}_I(x) − N^{B,u}_I(x)\), where \(N^{A,u}_I(x)\) represents the pristine enrichment functions and \(N^{B,u}_I(x)\f is the component of enrichment function along to the classical basis and can be written in terms of classical basis functions as \)@_fakenlN^{B,u}_I(x) = \sum^{n_h}_{l=1}c_{I,l}N^{C}_{l}(x)\f. Since the \(N^{O,u}_I(x)\)'s are orthogonal to classical basis, one can solve the equation \(M^{cc}.c = d\) where M is the discretized overlap matrix in classical FE basis and \(M^{cc}_{j,l} = \integral_{\omega}N^{C}_{l}(x)N^{C}_{j}(x)dx\) and d is the RHS written as \(d_{I,j,k} = \integral_{\omega}N^{A,u}_{i,j}(x)N^{C}_{k}(x)dx\) and \(c\) is the basisInterfaceCoefficient written as a multivector.

Parameters
[in]atomCoordinatesVector of Coordinates of the atoms
[in]atomPartitionToleranceset the tolerance for partitioning the atomids
[in]commMPI_Comm object if defined with MPI
[in]cfeBasisDataStorageRhsFEBasisDataStorage object for RHS \(d_{I}\)
[in]cfeBasisDataStorageOverlapMatrixFEBasisDataStorage object for the OverlapMatrix \(M^{cc}\)
[in]linAlgOpContextThe linearAlgebraOperator context for solving the linear equation
[in]fieldNameThe fieldname of the enrichment function

Can also do via the M^(-1) route withot solving CG.

  linearAlgebra::MultiVector<ValueTypeBasisData, memorySpace> d(
        d_cfeBasisManager->getMPIPatternP2P(),
        linAlgOpContext,
        nTotalEnrichmentIds);
  d.setValue(0.0);

  FEBasisOperations<ValueTypeBasisData, ValueTypeBasisData, memorySpace,
  dim> cfeBasisOperations(cfeBasisDataStorageRhs,
  cellTimesNumVec);
  rootCout << "Begin creating integrateWithBasisValues\n";

Integrate this with different quarature rule. (i.e. adaptive for the enrichment functions) , inp will be in adaptive grid cfeBasisOperations.integrateWithBasisValues(quadValuesEnrichmentFunction, d_cfeBasisManager, d); rootCout << "End creating integrateWithBasisValues\n"; utils::throwException( cfeBasisDataStorageOverlapMatrix->getQuadratureRuleContainer() ->getQuadratureRuleAttributes() .getQuadratureFamily() == quadrature::QuadratureFamily::GLL, "The quadrature rule for integration of Classical FE dofs has to be GLL if Mc = d" "is not solved via a poisson solve. Contact developers if extra options are needed.");

std::shared_ptr<dftefe::linearAlgebra::OperatorContext<ValueTypeBasisData, ValueTypeBasisData, memorySpace>> MInvContext = std::make_shared<dftefe::basis::CFEOverlapInverseOpContextGLL<ValueTypeBasisData, ValueTypeBasisData, memorySpace, dim>> (*d_cfeBasisManager, cfeBasisDataStorageOverlapMatrix, linAlgOpContext);

MInvContext->apply(d,*basisInterfaceCoeff);

Here is the call graph for this function:

◆ EnrichmentClassicalInterfaceSpherical() [2/2]

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::EnrichmentClassicalInterfaceSpherical ( std::shared_ptr< const TriangulationBase triangulation,
size_type  feOrder,
std::shared_ptr< const atoms::AtomSphericalDataContainer >  atomSphericalDataContainer,
const double  atomPartitionTolerance,
const std::vector< std::string > &  atomSymbolVec,
const std::vector< utils::Point > &  atomCoordinatesVec,
const std::string  fieldName,
const utils::mpi::MPIComm comm 
)

This Constructor for augmenting the EFE basis with classical FE basis.

Parameters
[in]atomCoordinatesVector of Coordinates of the atoms
[in]atomPartitionToleranceset the tolerance for partitioning the atomids
[in]commMPI_Comm object if defined with MPI
[in]fieldNameThe fieldname of the enrichment function
Here is the call graph for this function:

◆ ~EnrichmentClassicalInterfaceSpherical()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::~EnrichmentClassicalInterfaceSpherical ( )
default

Destructor for the class.

Member Function Documentation

◆ getAtomCoordinatesVec()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::vector< utils::Point > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getAtomCoordinatesVec

◆ getAtomIdsPartition()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr< const AtomIdsPartition< dim > > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getAtomIdsPartition

◆ getAtomSphericalDataContainer()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr< const atoms::AtomSphericalDataContainer > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getAtomSphericalDataContainer

Function to return AtomSphericalDataContainerObject.

◆ getAtomSymbolVec()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::vector< std::string > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getAtomSymbolVec

◆ getCFEBasisDofHandler()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr< const BasisDofHandler > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getCFEBasisDofHandler
Here is the call graph for this function:

◆ getCFEBasisManager()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr< const BasisManager< ValueTypeBasisData, memorySpace > > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getCFEBasisManager
Here is the call graph for this function:

◆ getClassicalComponentCoeffMap()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
const std::unordered_map< global_size_type, std::vector< ValueTypeBasisData > > & dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getClassicalComponentCoeffMap
Here is the call graph for this function:

◆ getClassicalComponentLocalIdsMap()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
const std::unordered_map< global_size_type, utils::OptimizedIndexSet< size_type > > & dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getClassicalComponentLocalIdsMap
Here is the call graph for this function:

◆ getEnrichmentId()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
global_size_type dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getEnrichmentId ( size_type  cellId,
size_type  enrichmentCellLocalId 
) const

The localid is determined by the storage pattern of the components of basisInterfaceCoeff multivector. The multivector has a storage pattern of [LocallyownedEnrichmentIds, GhostEnrichemntIds].

Here is the call graph for this function:

◆ getEnrichmentIdsPartition()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr< const EnrichmentIdsPartition< dim > > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getEnrichmentIdsPartition

Function to return EnrichmentIDsObject.

Here is the caller graph for this function:

◆ getEnrichmentLocalId() [1/2]

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
size_type dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getEnrichmentLocalId ( global_size_type  enrichmentId) const
Here is the call graph for this function:

◆ getEnrichmentLocalId() [2/2]

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
size_type dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getEnrichmentLocalId ( size_type  cellId,
size_type  enrichmentCellLocalId 
) const
Here is the call graph for this function:

◆ getFEOrder()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
size_type dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getFEOrder

◆ getFieldName()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::string dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getFieldName

◆ getLinAlgOpContext()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr< linearAlgebra::LinAlgOpContext< memorySpace > > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getLinAlgOpContext
Here is the call graph for this function:

◆ getTriangulation()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr< const TriangulationBase > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getTriangulation

◆ isOrthogonalized()

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
bool dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::isOrthogonalized

Member Data Documentation

◆ d_atomCoordinatesVec

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
const std::vector<utils::Point> dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_atomCoordinatesVec
private

◆ d_atomIdsPartition

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr<const AtomIdsPartition<dim> > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_atomIdsPartition
private

◆ d_atomSphericalDataContainer

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr<const atoms::AtomSphericalDataContainer> dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_atomSphericalDataContainer
private

◆ d_atomSymbolVec

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
const std::vector<std::string> dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_atomSymbolVec
private

◆ d_cfeBasisDofHandler

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr< const FEBasisDofHandler<ValueTypeBasisData, memorySpace, dim> > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_cfeBasisDofHandler
private

◆ d_cfeBasisManager

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr<const FEBasisManager<ValueTypeBasisData, ValueTypeBasisData, memorySpace, dim> > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_cfeBasisManager
private

◆ d_enrichmentIdsPartition

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr<const EnrichmentIdsPartition<dim> > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_enrichmentIdsPartition
private

◆ d_enrichmentIdToClassicalLocalIdMap

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::unordered_map<global_size_type, utils::OptimizedIndexSet<size_type> > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_enrichmentIdToClassicalLocalIdMap
private

◆ d_enrichmentIdToInterfaceCoeffMap

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::unordered_map<global_size_type, std::vector<ValueTypeBasisData> > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_enrichmentIdToInterfaceCoeffMap
private

◆ d_feOrder

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
size_type dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_feOrder
private

◆ d_fieldName

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
const std::string dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_fieldName
private

◆ d_isOrthogonalized

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
bool dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_isOrthogonalized
private

◆ d_linAlgOpContext

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr<linearAlgebra::LinAlgOpContext<memorySpace> > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_linAlgOpContext
private

◆ d_overlappingEnrichmentIdsInCells

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::vector<std::vector<global_size_type> > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_overlappingEnrichmentIdsInCells
private

◆ d_triangulation

template<typename ValueTypeBasisData , utils::MemorySpace memorySpace, size_type dim>
std::shared_ptr<const TriangulationBase> dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::d_triangulation
private

The documentation for this class was generated from the following files: