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 BasisDofHandler > | getCFEBasisDofHandler () 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::Point > | getAtomCoordinatesVec () const |
std::string | getFieldName () const |
std::shared_ptr< const TriangulationBase > | getTriangulation () 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 TriangulationBase > | d_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::Point > | d_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 |
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:
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.
[in] | atomCoordinates | Vector of Coordinates of the atoms |
[in] | atomPartitionTolerance | set the tolerance for partitioning the atomids |
[in] | comm | MPI_Comm object if defined with MPI |
[in] | cfeBasisDataStorageRhs | FEBasisDataStorage object for RHS \(d_{I}\) |
[in] | cfeBasisDataStorageOverlapMatrix | FEBasisDataStorage object for the OverlapMatrix \(M^{cc}\) |
[in] | linAlgOpContext | The linearAlgebraOperator context for solving the linear equation |
[in] | fieldName | The 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);
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.
[in] | atomCoordinates | Vector of Coordinates of the atoms |
[in] | atomPartitionTolerance | set the tolerance for partitioning the atomids |
[in] | comm | MPI_Comm object if defined with MPI |
[in] | fieldName | The fieldname of the enrichment function |
|
default |
Destructor for the class.
std::vector< utils::Point > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getAtomCoordinatesVec |
std::shared_ptr< const AtomIdsPartition< dim > > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getAtomIdsPartition |
std::shared_ptr< const atoms::AtomSphericalDataContainer > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getAtomSphericalDataContainer |
Function to return AtomSphericalDataContainerObject.
std::vector< std::string > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getAtomSymbolVec |
std::shared_ptr< const BasisDofHandler > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getCFEBasisDofHandler |
std::shared_ptr< const BasisManager< ValueTypeBasisData, memorySpace > > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getCFEBasisManager |
const std::unordered_map< global_size_type, std::vector< ValueTypeBasisData > > & dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getClassicalComponentCoeffMap |
const std::unordered_map< global_size_type, utils::OptimizedIndexSet< size_type > > & dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getClassicalComponentLocalIdsMap |
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].
std::shared_ptr< const EnrichmentIdsPartition< dim > > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getEnrichmentIdsPartition |
Function to return EnrichmentIDsObject.
size_type dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getEnrichmentLocalId | ( | global_size_type | enrichmentId | ) | const |
size_type dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getEnrichmentLocalId | ( | size_type | cellId, |
size_type | enrichmentCellLocalId | ||
) | const |
size_type dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getFEOrder |
std::string dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getFieldName |
std::shared_ptr< linearAlgebra::LinAlgOpContext< memorySpace > > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getLinAlgOpContext |
std::shared_ptr< const TriangulationBase > dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::getTriangulation |
bool dftefe::basis::EnrichmentClassicalInterfaceSpherical< ValueTypeBasisData, memorySpace, dim >::isOrthogonalized |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |
|
private |