v0.14.0
DGProjection.hpp
Go to the documentation of this file.
1 /*/**
2  * @file DGProjection.hpp
3  * @author your name (you@domain.com)
4  * @brief
5  * @version 0.1
6  * @date 2023-08-23
7  *
8  * @copyright Copyright (c) 2023
9  *
10  */
11 
12 #ifndef DGPROJECTION_HPP
13 #define DGPROJECTION_HPP
14 
15 namespace MoFEM {
16 
17 /**
18  * @brief Execute "this" element in the operator
19  *
20  * \code An example of usage in post-processing
21 
22 using InteriorElem = VolumeElementForcesAndSourcesCore;
23 using PostProcFe = PostProcBrokenMeshInMoab<InteriorElem>;
24 auto fe_post_proc = boost::make_shared<PostProcFe>(mField);
25 auto op_this = new OpLoopThis<InteriorElem>(m_field, "DomainFE", Sev::noisy);
26 fe_post_proc->getOpPtrVector()->push_back(op_this);
27 
28 auto fe_physics = op_this->getThisFEPtr();
29 fe_physics->getRuleHook = [](int, int, int o) {
30  return 2 * o; };
31 
32 int order = 1;
33 
34 // entity data shared between physical and post proc elements
35 auto entity_data_l2 = boost::make_shared<EntitiesFieldData>(MBENTITYSET);
36 // integrated mass matrix of post proc element
37 auto mass_ptr = boost::make_shared<MatrixDouble>();
38 // vector of coeffs shared between physical and post proc elements
39 auto coeffs_ptr = boost::make_shared<MatrixDouble>();
40 // data stored at
41 // integration points of the physical element and evaluated at integration
42 // points of the post proc element
43 auto data_ptr = boost::make_shared<MatrixDouble>();
44 fe_physics->getOpPtrVector()->push_back(new
45  OpDGProjectionMassMatrix(order, mass_ptr, entity_data_l2,
46  AINSWORTH_LEGENDRE_BASE, L2));
47 
48 // You need to call operatpor which will evalaute data_ptr
49 
50 fe_physics->getOpPtrVector()->push_back(new
51  OpCalculateVectorFieldValues("V", data_ptr);
52 fe_physics->getOpPtrVector()->push_back( new
53  OpDGProjectionCoefficients(data_ptr, coeffs_ptr, mass_ptr, entity_data_l2,
54  AINSWORTH_LEGENDRE_BASE, L2));
55 
56 fe_post_proc->getOpPtrVector()->push_back(new
57  OpDGProjectionEvaluation(data_ptr, coeffs_ptr,
58  entity_data_l2, AINSWORTH_LEGENDRE_BASE, L2));
59 
60  * \endcode
61  *
62  * @tparam E template for "this" element type
63  *
64  *
65  *
66  */
67 template <typename E>
69 
71 
72  /**
73  * @brief Construct a new OpThis object
74  *
75  * @param m_field
76  * @param fe_name name of "this" element (domain element)
77  */
78  OpLoopThis(MoFEM::Interface &m_field, const std::string fe_name,
79  const LogManager::SeverityLevel sev = Sev::noisy)
80  : UserDataOperator(NOSPACE, OPSPACE), thisFEPtr(new E(m_field)),
81  feName(fe_name), sevLevel(sev) {}
82 
83  MoFEMErrorCode doWork(int side, EntityType type,
88  };
89 
90  /**
91  * Here user will push operator evaluating field, or data at "this" element
92  * integration points.
93  */
94  boost::ptr_deque<UserDataOperator> &getOpPtrVector() {
95  return thisFEPtr->getOpPtrVector();
96  }
97 
98  boost::shared_ptr<E> &getThisFEPtr() { return thisFEPtr; }
99 
100 protected:
101  const std::string feName;
102  boost::shared_ptr<E> thisFEPtr;
104 };
105 
107  OpDGProjectionMassMatrix(int order, boost::shared_ptr<MatrixDouble> mass_ptr,
108  boost::shared_ptr<EntitiesFieldData> data_l2,
109  const FieldApproximationBase b, const FieldSpace s,
110  int verb = QUIET, Sev sev = Sev::verbose);
111 
112  MoFEMErrorCode doWork(int side, EntityType type,
114 
115 private:
117 };
118 
120 
121  OpDGProjectionCoefficients(boost::shared_ptr<MatrixDouble> data_ptr,
122  boost::shared_ptr<MatrixDouble> coeffs_ptr,
123  boost::shared_ptr<MatrixDouble> mass_ptr,
124  boost::shared_ptr<EntitiesFieldData> data_l2,
125  const FieldApproximationBase b, const FieldSpace s,
126  const LogManager::SeverityLevel sev = Sev::noisy);
127 
128  MoFEMErrorCode doWork(int side, EntityType type,
130 
131 protected:
132  boost::shared_ptr<MatrixDouble> dataPtr;
133  boost::shared_ptr<MatrixDouble> coeffsPtr;
134 };
135 
137 
138  OpDGProjectionEvaluation(boost::shared_ptr<MatrixDouble> data_ptr,
139  boost::shared_ptr<MatrixDouble> coeffs_ptr,
140  boost::shared_ptr<EntitiesFieldData> data_l2,
141  const FieldApproximationBase b, const FieldSpace s,
142  const LogManager::SeverityLevel sev = Sev::noisy);
143 
144  MoFEMErrorCode doWork(int side, EntityType type,
146 };
147 }; // namespace MoFEM
148 
149 #endif // DGPROJECTION_HPP
NOSPACE
@ NOSPACE
Definition: definitions.h:83
MoFEM::OpLoopThis::sevLevel
const LogManager::SeverityLevel sevLevel
Definition: DGProjection.hpp:103
MoFEM::EntitiesFieldData::EntData
Data on single entity (This is passed as argument to DataOperator::doWork)
Definition: EntitiesFieldData.hpp:128
MoFEM::OpDGProjectionMassMatrix::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: DGProjection.cpp:15
MoFEM::OpLoopThis::feName
const std::string feName
Definition: DGProjection.hpp:101
MoFEM::OpLoopThis
Execute "this" element in the operator.
Definition: DGProjection.hpp:68
MoFEM::ForcesAndSourcesCore::UserDataOperator::OPSPACE
@ OPSPACE
operator do Work is execute on space data
Definition: ForcesAndSourcesCore.hpp:570
MoFEM::OpDGProjectionEvaluation::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: DGProjection.cpp:103
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
E
MoFEM::OpDGProjectionCoefficients
Definition: DGProjection.hpp:119
order
constexpr int order
Definition: dg_projection.cpp:18
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
MoFEM::OpLoopThis::getOpPtrVector
boost::ptr_deque< UserDataOperator > & getOpPtrVector()
Definition: DGProjection.hpp:94
FieldSpace
FieldSpace
approximation spaces
Definition: definitions.h:82
VERBOSE
@ VERBOSE
Definition: definitions.h:222
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::OpLoopThis::OpLoopThis
OpLoopThis(MoFEM::Interface &m_field, const std::string fe_name, const LogManager::SeverityLevel sev=Sev::noisy)
Construct a new OpThis object.
Definition: DGProjection.hpp:78
MoFEM::OpDGProjectionCoefficients::coeffsPtr
boost::shared_ptr< MatrixDouble > coeffsPtr
Definition: DGProjection.hpp:133
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::OpDGProjectionCoefficients::dataPtr
boost::shared_ptr< MatrixDouble > dataPtr
Definition: DGProjection.hpp:132
MoFEM::ForcesAndSourcesCore::UserDataOperator
Definition: ForcesAndSourcesCore.hpp:549
convert.type
type
Definition: convert.py:64
MoFEM::OpLoopThis::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: DGProjection.hpp:83
MoFEM::LogManager::SeverityLevel
SeverityLevel
Severity levels.
Definition: LogManager.hpp:33
MoFEM::ForcesAndSourcesCore::UserDataOperator
friend class UserDataOperator
Definition: ForcesAndSourcesCore.hpp:482
MoFEM::OpBaseDerivativesBase
Definition: BaseDerivativesDataOperators.hpp:13
MoFEM::OpLoopThis::getThisFEPtr
boost::shared_ptr< E > & getThisFEPtr()
Definition: DGProjection.hpp:98
MoFEM::OpDGProjectionEvaluation::OpDGProjectionEvaluation
OpDGProjectionEvaluation(boost::shared_ptr< MatrixDouble > data_ptr, boost::shared_ptr< MatrixDouble > coeffs_ptr, boost::shared_ptr< EntitiesFieldData > data_l2, const FieldApproximationBase b, const FieldSpace s, const LogManager::SeverityLevel sev=Sev::noisy)
Definition: DGProjection.cpp:92
MoFEM::OpLoopThis::thisFEPtr
boost::shared_ptr< E > thisFEPtr
Definition: DGProjection.hpp:102
MoFEM::OpDGProjectionEvaluation
Definition: DGProjection.hpp:136
MoFEM::OpDGProjectionMassMatrix::OpDGProjectionMassMatrix
OpDGProjectionMassMatrix(int order, boost::shared_ptr< MatrixDouble > mass_ptr, boost::shared_ptr< EntitiesFieldData > data_l2, const FieldApproximationBase b, const FieldSpace s, int verb=QUIET, Sev sev=Sev::verbose)
Definition: DGProjection.cpp:7
FieldApproximationBase
FieldApproximationBase
approximation base
Definition: definitions.h:58
MoFEM::OpDGProjectionMassMatrix::baseOrder
int baseOrder
Definition: DGProjection.hpp:116
QUIET
@ QUIET
Definition: definitions.h:221
MoFEM::OpDGProjectionCoefficients::OpDGProjectionCoefficients
OpDGProjectionCoefficients(boost::shared_ptr< MatrixDouble > data_ptr, boost::shared_ptr< MatrixDouble > coeffs_ptr, boost::shared_ptr< MatrixDouble > mass_ptr, boost::shared_ptr< EntitiesFieldData > data_l2, const FieldApproximationBase b, const FieldSpace s, const LogManager::SeverityLevel sev=Sev::noisy)
Definition: DGProjection.cpp:34
MoFEM::OpDGProjectionCoefficients::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: DGProjection.cpp:45
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MoFEM::ForcesAndSourcesCore::UserDataOperator::loopThis
MoFEMErrorCode loopThis(const string &fe_name, ForcesAndSourcesCore *this_fe, const int verb=QUIET, const LogManager::SeverityLevel sev=Sev::noisy)
User calls this function to loop over the same element using a different set of integration points....
Definition: ForcesAndSourcesCore.cpp:1792
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
MoFEM::OpDGProjectionMassMatrix
Definition: DGProjection.hpp:106