v0.15.0
Loading...
Searching...
No Matches
FieldEvaluator.hpp
Go to the documentation of this file.
1/** \file FieldEvaluator.hpp
2 * \brief Field Evaluator
3 *
4 * Evaluate field at given coordinate
5 *
6 *
7 * \ingroup field_evaluator
8 */
9
10#ifndef __FIELD_EVALUATOR_HPP__
11 #define __FIELD_EVALUATOR_HPP__
12
13 #include "UnknownInterface.hpp"
14
15namespace MoFEM {
16
17/** \brief Field evaluator interface
18
19 * \ingroup field_evaluator
20 */
22
23 MoFEMErrorCode query_interface(boost::typeindex::type_index type_index,
24 UnknownInterface **iface) const;
25
28
29 struct SetPtsData {
30
31 SetPtsData() = delete;
32
33 /**
34 * @brief Set the Gauss Pts data
35 *
36 * @param fe_method_ptr pointer to finite element instance
37 * @param eval_points pointer to array with evaluation points
38 * @param nb_eval_points number of evaluated points
39 * @param eps tolerance used to find if point is in the element
40 * @param verb
41 */
42 SetPtsData(boost::shared_ptr<MoFEM::ForcesAndSourcesCore> fe_method_ptr,
43 const double *eval_points, const int nb_eval_points,
44 const double eps, VERBOSITY_LEVELS verb = QUIET)
45 : feMethodPtr(fe_method_ptr), evalPoints(eval_points),
46 nbEvalPoints(nb_eval_points), eps(eps), verb(verb) {
47 localCoords.resize(nbEvalPoints, 3);
49 }
50
51 inline void setEvalPoints(const double *ptr, const int nb_eval_points) {
52 evalPoints = ptr;
53 nbEvalPoints = nb_eval_points;
54 localCoords.resize(nbEvalPoints, 3, false);
55 shapeFunctions.resize(nbEvalPoints, 4, false);
56 }
57
58 boost::weak_ptr<MoFEM::ForcesAndSourcesCore> feMethodPtr;
59 const double *evalPoints;
61 double eps;
63
66 std::vector<EntityHandle> evalPointEntityHandle;
67
69 boost::scoped_ptr<AdaptiveKDTree> treePtr;
70 };
71
72 /**
73 * @brief Default evaluator for setting integration points
74 *
75 */
76 struct SetPts {
77 SetPts() = delete;
78 SetPts(boost::shared_ptr<SetPtsData> data_ptr) : dataPtr(data_ptr) {}
79 MoFEMErrorCode operator()(ForcesAndSourcesCore *fe_raw_ptr, int order_row,
80 int order_col, int order_data);
81
82 private:
83 boost::weak_ptr<SetPtsData> dataPtr;
84 };
85
86 /**
87 * @brief Get the Data object
88 *
89 * Pack pointers with data structures for field evaluator and finite
90 * element. Function return shared pointer if returned shared pointer
91 * is reset; all data are destroyed. The idea is to pack all data in
92 * one structure, create shared pointer to it and return aliased shared
93 * pointer to one of the elements of the data structure. It is a bit
94 * complicated, but has great flexibility.
95 *
96 * @tparam VE
97 * @tparam SetPtsData
98 * @tparam SetPts
99 * @param ptr
100 * @param nb_eval_points
101 * @param eps
102 * @param verb
103 * @return boost::shared_ptr<SPD>
104 */
105 template <typename VE, typename SPD = SetPtsData, typename SP = SetPts>
106 boost::shared_ptr<SPD>
107 getData(const double *ptr = nullptr, const int nb_eval_points = 0,
108 const double eps = 1e-12, VERBOSITY_LEVELS verb = QUIET) {
109 struct PackData {
110 boost::scoped_ptr<VE> elePtr;
111 boost::scoped_ptr<SPD> setPtsDataPtr;
112 boost::scoped_ptr<SP> setPtsPtr;
113 };
114 boost::shared_ptr<PackData> pack_data(new PackData());
115 MoFEM::Interface &m_field = cOre;
116 pack_data->elePtr.reset(new VE(m_field));
117
118 pack_data->setPtsDataPtr.reset(
119 new SPD(boost::shared_ptr<VE>(pack_data, pack_data->elePtr.get()), ptr,
120 nb_eval_points, eps, verb));
121 pack_data->setPtsPtr.reset(new SP(
122 boost::shared_ptr<SPD>(pack_data, pack_data->setPtsDataPtr.get())));
123 pack_data->elePtr->setRuleHook = boost::ref(*pack_data->setPtsPtr);
124 boost::shared_ptr<SPD> data(pack_data, pack_data->setPtsDataPtr.get());
125 return data;
126 }
127
128 /** @deprecated use buildTree<SPACE_DIM> */
130 MoFEMErrorCode buildTree3D(boost::shared_ptr<SetPtsData> spd_ptr,
131 const std::string finite_element);
132
133 /** @deprecated use buildTree<SPACE_DIM> */
135 MoFEMErrorCode buildTree2D(boost::shared_ptr<SetPtsData> spd_ptr,
136 const std::string finite_element);
137
138 /** @deprecated use evalFEAtThePoint<SPACE_DIM> */
140 const double *const point, const double distance,
141 const std::string problem, const std::string finite_element,
142 boost::shared_ptr<SetPtsData> data_ptr, int lower_rank, int upper_rank,
143 boost::shared_ptr<CacheTuple> cache_ptr, MoFEMTypes bh = MF_EXIST,
144 VERBOSITY_LEVELS verb = QUIET);
145
146 /** @deprecated use evalFEAtThePoint<SPACE_DIM> */
148 const double *const point, const double distance,
149 const std::string problem, const std::string finite_element,
150 boost::shared_ptr<SetPtsData> data_ptr, int lower_rank, int upper_rank,
151 boost::shared_ptr<CacheTuple> cache_ptr, MoFEMTypes bh = MF_EXIST,
152 VERBOSITY_LEVELS verb = QUIET);
153
154 /**
155 * @brief Build spatial tree
156 *
157 * @param finite_element finite element name
158 * @return MoFEMErrorCode
159 */
160 template <int D>
161 MoFEMErrorCode buildTree(boost::shared_ptr<SetPtsData> spd_ptr,
162 const std::string finite_element);
163
164 /**
165 * @brief Evaluate field at artbitray position
166 *
167 * \code
168
169 std::array<double, 3> point = {0, 0, 0};
170 const double dist = 0.3;
171 std::array<double, 6> eval_points = {-1., -1., -1., 1., 1., 1. };
172
173 using VolEle = VolumeElementForcesAndSourcesCore;
174 auto data_ptr =
175 m_field.getInterface<FieldEvaluatorInterface>()->
176 getData<VolEle>(point.data(), point.size/3);
177
178 if(auto vol_ele = data_ptr->feMethod.lock()) {
179 // push operators to finite element instance, e.g.
180 vol_ele->getOpPtrVector().push_back(new MyOp());
181 // iterate over elemnts with evaluated points
182 auto cache_ptr = boost::make_shared<CacheTuple>();
183 CHKERR m_field.cache_problem_entities(prb_ptr->getName(), cache_ptr);
184 // SPACE_DIM is 3 in this example
185 CHKERR m_field.getInterface<FieldEvaluatorInterface>()
186 ->evalFEAtThePoint<3>(point.data(), dist, prb_ptr->getName(),
187 "FINITE_ELEMENT_NAME",
188 data_ptr, m_field.get_comm_rank(),
189 m_field.get_comm_rank(), cache_ptr);
190 }
191
192 * \endcode
193 *
194 * @param point point used to find tetrahedrons
195 * @param distance distance from the point where tetrahedrons are searched
196 * @param problem problem name
197 * @param finite_element finite element name
198 * @param data_ptr pointer to data abut gauss points
199 * @param lower_rank lower processor
200 * @param upper_rank upper process
201 * @param cache_ptr cache or problem entities
202 * @param bh control looping over entities, e.g. throwing error if element not
203 found
204 * @param verb verbosity level
205 * @return MoFEMErrorCode
206 */
207 template <int D>
209 evalFEAtThePoint(const double *const point, const double distance,
210 const std::string problem, const std::string finite_element,
211 boost::shared_ptr<SetPtsData> data_ptr, int lower_rank,
212 int upper_rank, boost::shared_ptr<CacheTuple> cache_ptr,
214
215private:
216 template <int D>
217 MoFEMErrorCode buildTreeImpl(boost::shared_ptr<SetPtsData> spd_ptr,
218 const std::string finite_element);
219
220 template <int D>
222 evalFEAtThePointImpl(const double *const point, const double distance,
223 const std::string problem,
224 const std::string finite_element,
225 boost::shared_ptr<SetPtsData> data_ptr, int lower_rank,
226 int upper_rank, boost::shared_ptr<CacheTuple> cache_ptr,
228};
229
230template <>
232FieldEvaluatorInterface::buildTree<2>(boost::shared_ptr<SetPtsData> spd_ptr,
233 const std::string finite_element);
234
235template <>
237FieldEvaluatorInterface::buildTree<3>(boost::shared_ptr<SetPtsData> spd_ptr,
238 const std::string finite_element);
239
240template <>
242 const double *const point, const double distance, const std::string problem,
243 const std::string finite_element, boost::shared_ptr<SetPtsData> data_ptr,
244 int lower_rank, int upper_rank, boost::shared_ptr<CacheTuple> cache_ptr,
246
247template <>
249 const double *const point, const double distance, const std::string problem,
250 const std::string finite_element, boost::shared_ptr<SetPtsData> data_ptr,
251 int lower_rank, int upper_rank, boost::shared_ptr<CacheTuple> cache_ptr,
253
254} // namespace MoFEM
255
256#endif // __FIELD_EVALUATOR_HPP__
257
258/**
259 * \defgroup field_evaluator Field Evaluator
260 * \brief Evaluate field at the point
261 *
262 * \ingroup mofem
263 */
MoFEM interface.
static const double eps
VERBOSITY_LEVELS
Verbosity levels.
@ QUIET
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
@ MF_EXIST
#define DEPRECATED
Definition definitions.h:17
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
Core (interface) class.
Definition Core.hpp:82
Deprecated interface functions.
boost::weak_ptr< MoFEM::ForcesAndSourcesCore > feMethodPtr
void setEvalPoints(const double *ptr, const int nb_eval_points)
boost::scoped_ptr< AdaptiveKDTree > treePtr
SetPtsData(boost::shared_ptr< MoFEM::ForcesAndSourcesCore > fe_method_ptr, const double *eval_points, const int nb_eval_points, const double eps, VERBOSITY_LEVELS verb=QUIET)
Set the Gauss Pts data.
std::vector< EntityHandle > evalPointEntityHandle
Default evaluator for setting integration points.
MoFEMErrorCode operator()(ForcesAndSourcesCore *fe_raw_ptr, int order_row, int order_col, int order_data)
SetPts(boost::shared_ptr< SetPtsData > data_ptr)
boost::weak_ptr< SetPtsData > dataPtr
Field evaluator interface.
DEPRECATED MoFEMErrorCode buildTree3D(boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
DEPRECATED MoFEMErrorCode buildTree2D(boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
MoFEMErrorCode evalFEAtThePointImpl(const double *const point, const double distance, const std::string problem, const std::string finite_element, boost::shared_ptr< SetPtsData > data_ptr, int lower_rank, int upper_rank, boost::shared_ptr< CacheTuple > cache_ptr, MoFEMTypes bh=MF_EXIST, VERBOSITY_LEVELS verb=QUIET)
DEPRECATED MoFEMErrorCode evalFEAtThePoint2D(const double *const point, const double distance, const std::string problem, const std::string finite_element, boost::shared_ptr< SetPtsData > data_ptr, int lower_rank, int upper_rank, boost::shared_ptr< CacheTuple > cache_ptr, MoFEMTypes bh=MF_EXIST, VERBOSITY_LEVELS verb=QUIET)
boost::shared_ptr< SPD > getData(const double *ptr=nullptr, const int nb_eval_points=0, const double eps=1e-12, VERBOSITY_LEVELS verb=QUIET)
Get the Data object.
FieldEvaluatorInterface(const MoFEM::Core &core)
MoFEMErrorCode buildTreeImpl(boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
DEPRECATED MoFEMErrorCode evalFEAtThePoint3D(const double *const point, const double distance, const std::string problem, const std::string finite_element, boost::shared_ptr< SetPtsData > data_ptr, int lower_rank, int upper_rank, boost::shared_ptr< CacheTuple > cache_ptr, MoFEMTypes bh=MF_EXIST, VERBOSITY_LEVELS verb=QUIET)
MoFEMErrorCode buildTree(boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
Build spatial tree.
MoFEMErrorCode evalFEAtThePoint(const double *const point, const double distance, const std::string problem, const std::string finite_element, boost::shared_ptr< SetPtsData > data_ptr, int lower_rank, int upper_rank, boost::shared_ptr< CacheTuple > cache_ptr, MoFEMTypes bh=MF_EXIST, VERBOSITY_LEVELS verb=QUIET)
Evaluate field at artbitray position.
structure to get information form mofem into EntitiesFieldData
base class for all interface classes