v0.15.0
Loading...
Searching...
No Matches
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
15namespace MoFEM {
16
17/**
18 * @brief Execute "this" element in the operator
19 *
20 * \code An example of usage in post-processing
21
22using InteriorElem = VolumeElementForcesAndSourcesCore;
23using PostProcFe = PostProcBrokenMeshInMoab<InteriorElem>;
24auto fe_post_proc = boost::make_shared<PostProcFe>(mField);
25auto op_this = new OpLoopThis<InteriorElem>(m_field, "DomainFE", Sev::noisy);
26fe_post_proc->getOpPtrVector()->push_back(op_this);
27
28auto fe_physics = op_this->getThisFEPtr();
29fe_physics->getRuleHook = [](int, int, int o) {
30 return 2 * o; };
31
32int order = 1;
33
34// entity data shared between physical and post proc elements
35auto entity_data_l2 = boost::make_shared<EntitiesFieldData>(MBENTITYSET);
36// integrated mass matrix of post proc element
37auto mass_ptr = boost::make_shared<MatrixDouble>();
38// vector of coeff shared between physical and post proc elements
39auto 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
43auto data_ptr = boost::make_shared<MatrixDouble>();
44fe_physics->getOpPtrVector()->push_back(new
45 OpDGProjectionMassMatrix(order, mass_ptr, entity_data_l2,
46 AINSWORTH_LEGENDRE_BASE, L2));
47
48// You need to call operator which will evaluate data_ptr
49
50fe_physics->getOpPtrVector()->push_back(new
51 OpCalculateVectorFieldValues("V", data_ptr);
52fe_physics->getOpPtrVector()->push_back( new
53 OpDGProjectionCoefficients(data_ptr, coeffs_ptr, mass_ptr, entity_data_l2,
54 AINSWORTH_LEGENDRE_BASE, L2));
55
56fe_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 */
67template <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
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
100protected:
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
115private:
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
131protected:
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
@ QUIET
@ VERBOSE
FieldApproximationBase
approximation base
Definition definitions.h:58
FieldSpace
approximation spaces
Definition definitions.h:82
@ NOSPACE
Definition definitions.h:83
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
constexpr int order
SeverityLevel
Severity levels.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
@ OPSPACE
operator do Work is execute on space data
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....
boost::shared_ptr< MatrixDouble > coeffsPtr
boost::shared_ptr< MatrixDouble > dataPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
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)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
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)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
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)
Execute "this" element in the operator.
const std::string feName
boost::shared_ptr< E > thisFEPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< E > & getThisFEPtr()
OpLoopThis(MoFEM::Interface &m_field, const std::string fe_name, const LogManager::SeverityLevel sev=Sev::noisy)
Construct a new OpThis object.
boost::ptr_deque< UserDataOperator > & getOpPtrVector()
const LogManager::SeverityLevel sevLevel