v0.13.1
Classes | Public Member Functions | Public Attributes | Private Member Functions | List of all members
MoFEM::FieldEvaluatorInterface Struct Reference

Field evaluator interface. More...

#include <src/interfaces/FieldEvaluator.hpp>

Inheritance diagram for MoFEM::FieldEvaluatorInterface:
[legend]
Collaboration diagram for MoFEM::FieldEvaluatorInterface:
[legend]

Classes

struct  SetPts
 Default evaluator for setting integration points. More...
 
struct  SetPtsData
 

Public Member Functions

MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 FieldEvaluatorInterface (const MoFEM::Core &core)
 
template<typename VE , typename SPD = SetPtsData, typename SP = SetPts>
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. More...
 
MoFEMErrorCode buildTree3D (boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
 Build spatial tree. More...
 
MoFEMErrorCode buildTree2D (boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
 Build spatial tree. More...
 
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)
 Evaluate field at artbitray position. More...
 
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)
 Evaluate field at artbitray position. More...
 
- Public Member Functions inherited from MoFEM::UnknownInterface
virtual MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const =0
 
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()=default
 

Public Attributes

MoFEM::CorecOre
 

Private Member Functions

template<int D>
MoFEMErrorCode buildTree (boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
 Build spatial tree. More...
 
template<int D>
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. More...
 

Additional Inherited Members

- Static Public Member Functions inherited from MoFEM::UnknownInterface
static MoFEMErrorCode getLibVersion (Version &version)
 Get library version. More...
 
static MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version)
 Get database major version. More...
 
static MoFEMErrorCode setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
 Get database major version. More...
 
static MoFEMErrorCode getInterfaceVersion (Version &version)
 Get database major version. More...
 

Detailed Description

Field evaluator interface.

Examples
EshelbianPlasticity.cpp, field_evaluator.cpp, and lorentz_force.cpp.

Definition at line 21 of file FieldEvaluator.hpp.

Constructor & Destructor Documentation

◆ FieldEvaluatorInterface()

MoFEM::FieldEvaluatorInterface::FieldEvaluatorInterface ( const MoFEM::Core core)

Definition at line 10 of file FieldEvaluator.cpp.

11 : cOre(const_cast<MoFEM::Core &>(core)) {
12
13 if (!LogManager::checkIfChannelExist("FieldEvaluatorWorld")) {
14 auto core_log = logging::core::get();
15
17 "FieldEvaluatorWorld"));
19 "FieldEvaluatorSync"));
21 "FieldEvaluatorSelf"));
22
23 LogManager::setLog("FieldEvaluatorWorld");
24 LogManager::setLog("FieldEvaluatorSync");
25 LogManager::setLog("FieldEvaluatorSelf");
26
27 MOFEM_LOG_TAG("FieldEvaluatorWorld", "FieldEvaluator");
28 MOFEM_LOG_TAG("FieldEvaluatorSync", "FieldEvaluator");
29 MOFEM_LOG_TAG("FieldEvaluatorSelf", "FieldEvaluator");
30 }
31
32 MOFEM_LOG("FieldEvaluatorWorld", Sev::noisy) << "Field evaluator intreface";
33}
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:364
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:301
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:332
Core (interface) class.
Definition: Core.hpp:82
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:279
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:323
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.
Definition: LogManager.cpp:327
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
Definition: LogManager.cpp:374
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Definition: LogManager.cpp:319

Member Function Documentation

◆ buildTree()

template<int D>
MoFEMErrorCode MoFEM::FieldEvaluatorInterface::buildTree ( boost::shared_ptr< SetPtsData spd_ptr,
const std::string  finite_element 
)
private

Build spatial tree.

Parameters
finite_elementfinite element name
Returns
MoFEMErrorCode

Definition at line 44 of file FieldEvaluator.cpp.

45 {
46 CoreInterface &m_field = cOre;
48 EntityHandle fe_meshset;
49 fe_meshset = m_field.get_finite_element_meshset(finite_element);
50 Range entities;
51 CHKERR m_field.get_moab().get_entities_by_dimension(fe_meshset, D, entities,
52 true);
53 spd_ptr->treePtr.reset(new AdaptiveKDTree(&m_field.get_moab()));
54 CHKERR spd_ptr->treePtr->build_tree(entities, &spd_ptr->rooTreeSet);
56}
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
const double D
diffusivity

◆ buildTree2D()

MoFEMErrorCode MoFEM::FieldEvaluatorInterface::buildTree2D ( boost::shared_ptr< SetPtsData spd_ptr,
const std::string  finite_element 
)

Build spatial tree.

Parameters
finite_elementfinite element name
Returns
MoFEMErrorCode
Examples
field_evaluator.cpp.

Definition at line 65 of file FieldEvaluator.cpp.

66 {
67 return buildTree<2>(spd_ptr, finite_element);
68}

◆ buildTree3D()

MoFEMErrorCode MoFEM::FieldEvaluatorInterface::buildTree3D ( boost::shared_ptr< SetPtsData spd_ptr,
const std::string  finite_element 
)

Build spatial tree.

Parameters
finite_elementfinite element name
Returns
MoFEMErrorCode
Examples
field_evaluator.cpp, and lorentz_force.cpp.

Definition at line 59 of file FieldEvaluator.cpp.

60 {
61 return buildTree<3>(spd_ptr, finite_element);
62}

◆ evalFEAtThePoint()

template<int D>
MoFEMErrorCode MoFEM::FieldEvaluatorInterface::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 
)
private

Evaluate field at artbitray position.

std::array<double, 3> point = {0, 0, 0};
const double dist = 0.3;
std::array<double, 6> eval_points = {-1., -1., -1., 1., 1., 1. };
using VolEle = VolumeElementForcesAndSourcesCore;
auto data_ptr =
getData<VolEle>(point.data(), point.size/3);
if(auto vol_ele = data_ptr->feMethod.lock()) {
// push operators to finite element instance, e.g.
vol_ele->getOpPtrVector().push_back(new MyOp());
// iterate over elemnts with evaluated points
auto cache_ptr = boost::make_shared<CacheTuple>();
CHKERR m_field.cache_problem_entities(prb_ptr->getName(), cache_ptr);
CHKERR m_field.getInterface<FieldEvaluatorInterface>()
->evalFEAtThePoint3D(point.data(), dist, prb_ptr->getName(),
"FINITE_ELEMENT_NAME",
data_ptr, m_field.get_comm_rank(),
m_field.get_comm_rank(), cache_ptr);
}
FieldEvaluatorInterface(const MoFEM::Core &core)
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)
Evaluate field at artbitray position.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
Operator used to check consistency between local coordinates and global cooridnates for integrated po...
Parameters
pointpoint used to find tetrahedrons
distancedistance from the point where tetrahedrons are searched
problemproblem name
finite_elementfinite element name
data_ptrpointer to data abut gauss points
lower_ranklower processor
upper_rankupper process
cache_ptrcache or problem entities
bhcontrol looping over entities, e.g. throwing error if element not found
verbverbosity level
Returns
MoFEMErrorCode

Definition at line 163 of file FieldEvaluator.cpp.

167 {
168 CoreInterface &m_field = cOre;
170
171 std::vector<EntityHandle> leafs_out;
172 CHKERR data_ptr->treePtr->distance_search(point, distance, leafs_out);
173 Range tree_ents;
174 for (auto lit : leafs_out)
175 CHKERR m_field.get_moab().get_entities_by_dimension(lit, D, tree_ents,
176 true);
177 if (verb >= NOISY)
178 MOFEM_LOG("SELF", Sev::noisy) << "tree entities: " << tree_ents;
179
180 data_ptr->evalPointEntityHandle.resize(data_ptr->nbEvalPoints);
181 std::fill(data_ptr->evalPointEntityHandle.begin(),
182 data_ptr->evalPointEntityHandle.end(), 0);
183
184 std::vector<double> local_coords;
185 std::vector<double> shape;
186
187 auto nb_eval_points = data_ptr->nbEvalPoints;
188
189 for (auto tet : tree_ents) {
190 const EntityHandle *conn;
191 int num_nodes;
192 CHKERR m_field.get_moab().get_connectivity(tet, conn, num_nodes, true);
193
194 if constexpr (D == 3) {
195
196 local_coords.resize(3 * nb_eval_points);
197 shape.resize(4 * nb_eval_points);
198
199 std::array<double, 12> coords;
200 CHKERR m_field.get_moab().get_coords(conn, num_nodes, coords.data());
201
203 coords.data(), &data_ptr->evalPoints[0], nb_eval_points,
204 &local_coords[0]);
205 CHKERR Tools::shapeFunMBTET<3>(&shape[0], &local_coords[0],
206 &local_coords[1], &local_coords[2],
207 nb_eval_points);
208
211 &shape[0], &shape[1], &shape[2], &shape[3]};
213 &data_ptr->shapeFunctions(0, 0), &data_ptr->shapeFunctions(0, 1),
214 &data_ptr->shapeFunctions(0, 2), &data_ptr->shapeFunctions(0, 3)};
215
218 &local_coords[0], &local_coords[1], &local_coords[2]};
220 &(data_ptr->localCoords(0, 0)), &(data_ptr->localCoords(0, 1)),
221 &(data_ptr->localCoords(0, 2))};
222
223
224 for (int n = 0; n != nb_eval_points; ++n) {
225
226 const double eps = data_ptr->eps;
227 if (t_shape(0) >= 0 - eps && t_shape(0) <= 1 + eps &&
228
229 t_shape(1) >= 0 - eps && t_shape(1) <= 1 + eps &&
230
231 t_shape(2) >= 0 - eps && t_shape(2) <= 1 + eps &&
232
233 t_shape(3) >= 0 - eps && t_shape(3) <= 1 + eps) {
234
235 data_ptr->evalPointEntityHandle[n] = tet;
236 t_shape_data(i4) = t_shape(i4);
237 t_local_data(j3) = t_local(j3);
238 }
239
240 ++t_shape;
241 ++t_shape_data;
242 ++t_local;
243 ++t_local_data;
244
245 }
246
247 }
248
249 if constexpr (D == 2) {
250
251 local_coords.resize(2 * nb_eval_points);
252 shape.resize(3 * nb_eval_points);
253
254 std::array<double, 9> coords;
255 CHKERR m_field.get_moab().get_coords(conn, num_nodes, coords.data());
256
258 coords.data(), &data_ptr->evalPoints[0], nb_eval_points,
259 &local_coords[0]);
260 CHKERR Tools::shapeFunMBTRI<2>(&shape[0], &local_coords[0],
261 &local_coords[1], nb_eval_points);
262
265 &shape[0], &shape[1], &shape[2]};
267 &data_ptr->shapeFunctions(0, 0), &data_ptr->shapeFunctions(0, 1),
268 &data_ptr->shapeFunctions(0, 2)};
269
272 &local_coords[0], &local_coords[1]};
273 data_ptr->localCoords.resize(nb_eval_points, 2);
275 &(data_ptr->localCoords(0, 0)), &(data_ptr->localCoords(0, 1))};
276
277 for (int n = 0; n != nb_eval_points; ++n) {
278
279 const double eps = data_ptr->eps;
280 if (t_shape(0) >= 0 - eps && t_shape(0) <= 1 + eps &&
281
282 t_shape(1) >= 0 - eps && t_shape(1) <= 1 + eps &&
283
284 t_shape(2) >= 0 - eps && t_shape(2) <= 1 + eps) {
285
286 data_ptr->evalPointEntityHandle[n] = tet;
287 t_shape_data(i3) = t_shape(i3);
288 t_local_data(j2) = t_local(j2);
289 }
290
291 ++t_shape;
292 ++t_shape_data;
293 ++t_local;
294 ++t_local_data;
295 }
296 }
297
298
299 }
300
301 const Problem *prb_ptr;
302 CHKERR m_field.get_problem(problem, &prb_ptr);
303 boost::shared_ptr<NumeredEntFiniteElement_multiIndex> numered_fes(
305
306 Range in_tets;
307 in_tets.insert_list(data_ptr->evalPointEntityHandle.begin(),
308 data_ptr->evalPointEntityHandle.end());
309 in_tets = in_tets.subset_by_dimension(D);
310
311 if (verb >= NOISY)
312 MOFEM_LOG("SELF", Sev::noisy) << "in tets: " << in_tets << endl;
313
314 for (auto peit = in_tets.pair_begin(); peit != in_tets.pair_end(); ++peit) {
315 auto lo =
316 prb_ptr->numeredFiniteElementsPtr->get<Composite_Name_And_Ent_mi_tag>()
317 .lower_bound(boost::make_tuple(finite_element, peit->first));
318 auto hi =
319 prb_ptr->numeredFiniteElementsPtr->get<Composite_Name_And_Ent_mi_tag>()
320 .upper_bound(boost::make_tuple(finite_element, peit->second));
321 numered_fes->insert(lo, hi);
322
323 if (verb >= VERY_NOISY)
324 std::cerr << "numered elements:" << std::endl;
325 for (; lo != hi; ++lo)
326 if (verb >= VERY_NOISY)
327 std::cerr << **lo << endl;
328 }
329 if (verb >= VERY_NOISY)
330 std::cerr << std::endl;
331
332 if (auto fe_ptr = data_ptr->feMethodPtr.lock()) {
333
334 MOFEM_LOG_C("FieldEvaluatorSync", Sev::verbose,
335 "Number elements %d to evaluate at proc %d",
336 numered_fes->size(), m_field.get_comm_rank());
337 MOFEM_LOG_SYNCHRONISE(m_field.get_comm());
338
339 if (!cache_ptr) {
340 cache_ptr = boost::make_shared<CacheTuple>();
341 CHKERR m_field.cache_problem_entities(prb_ptr->getName(), cache_ptr);
342
343 MOFEM_LOG("FieldEvaluatorSync", Sev::noisy)
344 << "If you call that function many times in the loop consider to "
345 "set "
346 "cache_ptr outside of the loop. Otherwise code can be slow.";
347 }
348
349 CHKERR m_field.loop_finite_elements(prb_ptr, finite_element, *fe_ptr,
350 lower_rank, upper_rank, numered_fes, bh,
351 cache_ptr, verb);
352
353 } else
354 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
355 "Pointer to element does not exists");
356
358}
#define MOFEM_LOG_SYNCHRONISE(comm)
Synchronise "SYNC" channel.
Definition: LogManager.hpp:338
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:304
static const double eps
@ NOISY
Definition: definitions.h:211
@ VERY_NOISY
Definition: definitions.h:212
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
FTensor::Index< 'n', SPACE_DIM > n
multi_index_container< boost::shared_ptr< NumeredEntFiniteElement >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< NumeredEntFiniteElement::interface_type_EntFiniteElement, UId, &NumeredEntFiniteElement::getLocalUniqueId > >, ordered_non_unique< tag< Part_mi_tag >, member< NumeredEntFiniteElement, unsigned int, &NumeredEntFiniteElement::part > >, ordered_non_unique< tag< FiniteElement_name_mi_tag >, const_mem_fun< NumeredEntFiniteElement::interface_type_FiniteElement, boost::string_ref, &NumeredEntFiniteElement::getNameRef > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< NumeredEntFiniteElement::interface_type_RefEntity, EntityHandle, &NumeredEntFiniteElement::getEnt > >, ordered_non_unique< tag< Composite_Name_And_Ent_mi_tag >, composite_key< NumeredEntFiniteElement, const_mem_fun< NumeredEntFiniteElement::interface_type_FiniteElement, boost::string_ref, &NumeredEntFiniteElement::getNameRef >, const_mem_fun< NumeredEntFiniteElement::interface_type_RefEntity, EntityHandle, &NumeredEntFiniteElement::getEnt > > >, ordered_non_unique< tag< Composite_Name_And_Part_mi_tag >, composite_key< NumeredEntFiniteElement, const_mem_fun< NumeredEntFiniteElement::interface_type_FiniteElement, boost::string_ref, &NumeredEntFiniteElement::getNameRef >, member< NumeredEntFiniteElement, unsigned int, &NumeredEntFiniteElement::part > > > > > NumeredEntFiniteElement_multiIndex
MultiIndex for entities for NumeredEntFiniteElement.
static MoFEMErrorCode getLocalCoordinatesOnReferenceFourNodeTet(const double *elem_coords, const double *glob_coords, const int nb_nodes, double *local_coords)
Get the Local Coordinates On Reference Four Node Tet object.
Definition: Tools.cpp:90
static MoFEMErrorCode getLocalCoordinatesOnReferenceTriNodeTri(const double *elem_coords, const double *glob_coords, const int nb_nodes, double *local_coords)
Get the Local Coordinates On Reference Four Node Tet object.
Definition: Tools.cpp:142

◆ evalFEAtThePoint2D()

MoFEMErrorCode MoFEM::FieldEvaluatorInterface::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 
)

Evaluate field at artbitray position.

std::array<double, 3> point = {0, 0, 0};
const double dist = 0.3;
std::array<double, 6> eval_points = {-1., -1., -1., 1., 1., 1. };
using VolEle = VolumeElementForcesAndSourcesCore;
auto data_ptr =
getData<VolEle>(point.data(), point.size/3);
if(auto vol_ele = data_ptr->feMethod.lock()) {
// push operators to finite element instance, e.g.
vol_ele->getOpPtrVector().push_back(new MyOp());
// iterate over elemnts with evaluated points
auto cache_ptr = boost::make_shared<CacheTuple>();
CHKERR m_field.cache_problem_entities(prb_ptr->getName(), cache_ptr);
CHKERR m_field.getInterface<FieldEvaluatorInterface>()
->evalFEAtThePoint3D(point.data(), dist, prb_ptr->getName(),
"FINITE_ELEMENT_NAME",
data_ptr, m_field.get_comm_rank(),
m_field.get_comm_rank(), cache_ptr);
}
Parameters
pointpoint used to find tetrahedrons
distancedistance from the point where tetrahedrons are searched
problemproblem name
finite_elementfinite element name
data_ptrpointer to data abut gauss points
lower_ranklower processor
upper_rankupper process
cache_ptrcache or problem entities
bhcontrol looping over entities, e.g. throwing error if element not found
verbverbosity level
Returns
MoFEMErrorCode
Examples
field_evaluator.cpp.

Definition at line 369 of file FieldEvaluator.cpp.

373 {
374 return evalFEAtThePoint<2>(point, distance, problem, finite_element, data_ptr,
375 lower_rank, upper_rank, cache_ptr, bh, verb);
376}

◆ evalFEAtThePoint3D()

MoFEMErrorCode MoFEM::FieldEvaluatorInterface::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 
)

Evaluate field at artbitray position.

std::array<double, 3> point = {0, 0, 0};
const double dist = 0.3;
std::array<double, 6> eval_points = {-1., -1., -1., 1., 1., 1. };
using VolEle = VolumeElementForcesAndSourcesCore;
auto data_ptr =
getData<VolEle>(point.data(), point.size/3);
if(auto vol_ele = data_ptr->feMethod.lock()) {
// push operators to finite element instance, e.g.
vol_ele->getOpPtrVector().push_back(new MyOp());
// iterate over elemnts with evaluated points
auto cache_ptr = boost::make_shared<CacheTuple>();
CHKERR m_field.cache_problem_entities(prb_ptr->getName(), cache_ptr);
CHKERR m_field.getInterface<FieldEvaluatorInterface>()
->evalFEAtThePoint3D(point.data(), dist, prb_ptr->getName(),
"FINITE_ELEMENT_NAME",
data_ptr, m_field.get_comm_rank(),
m_field.get_comm_rank(), cache_ptr);
}
Parameters
pointpoint used to find tetrahedrons
distancedistance from the point where tetrahedrons are searched
problemproblem name
finite_elementfinite element name
data_ptrpointer to data abut gauss points
lower_ranklower processor
upper_rankupper process
cache_ptrcache or problem entities
bhcontrol looping over entities, e.g. throwing error if element not found
verbverbosity level
Returns
MoFEMErrorCode
Examples
field_evaluator.cpp, and lorentz_force.cpp.

Definition at line 360 of file FieldEvaluator.cpp.

364 {
365 return evalFEAtThePoint<3>(point, distance, problem, finite_element, data_ptr,
366 lower_rank, upper_rank, cache_ptr, bh, verb);
367}

◆ getData()

template<typename VE , typename SPD = SetPtsData, typename SP = SetPts>
boost::shared_ptr< SPD > MoFEM::FieldEvaluatorInterface::getData ( const double ptr = nullptr,
const int  nb_eval_points = 0,
const double  eps = 1e-12,
VERBOSITY_LEVELS  verb = QUIET 
)

Get the Data object.

Pack pointers with data structures for field evaluator and finite element. Function return shared pointer if returned shared pointer is reset; all data are destroyed. The idea is to pack all data in one structure, create shared pointer to it and return aliased shared pointer to one of the elements of the data structure. It is a bit complicated, but has great flexibility.

Template Parameters
VE
SetPtsData
SetPts
Parameters
ptr
nb_eval_points
eps
verb
Returns
boost::shared_ptr<SPD>
Examples
EshelbianPlasticity.cpp, field_evaluator.cpp, and lorentz_force.cpp.

Definition at line 107 of file FieldEvaluator.hpp.

108 {
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 pack_data->setPtsDataPtr.reset(
118 new SPD(boost::shared_ptr<VE>(pack_data, pack_data->elePtr.get()), ptr,
119 nb_eval_points, eps, verb));
120 pack_data->setPtsPtr.reset(new SP(
121 boost::shared_ptr<SPD>(pack_data, pack_data->setPtsDataPtr.get())));
122 pack_data->elePtr->setRuleHook = boost::ref(*pack_data->setPtsPtr);
123 boost::shared_ptr<SPD> data(pack_data, pack_data->setPtsDataPtr.get());
124 return data;
125 }
Deprecated interface functions.

◆ query_interface()

MoFEMErrorCode MoFEM::FieldEvaluatorInterface::query_interface ( boost::typeindex::type_index  type_index,
UnknownInterface **  iface 
) const
virtual

Implements MoFEM::UnknownInterface.

Definition at line 35 of file FieldEvaluator.cpp.

36 {
38 *iface = const_cast<FieldEvaluatorInterface *>(this);
40}
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440

Member Data Documentation

◆ cOre

MoFEM::Core& MoFEM::FieldEvaluatorInterface::cOre

Definition at line 26 of file FieldEvaluator.hpp.


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