v0.14.0
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
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface reference 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
dynamic_first_order_con_law.cpp, elasticity.cpp, field_evaluator.cpp, lorentz_force.cpp, seepage.cpp, and thermo_elastic.cpp.

Definition at line 21 of file FieldEvaluator.hpp.

Constructor & Destructor Documentation

◆ FieldEvaluatorInterface()

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

Definition at line 8 of file FieldEvaluator.cpp.

9  : cOre(const_cast<MoFEM::Core &>(core)) {
10 
11  if (!LogManager::checkIfChannelExist("FieldEvaluatorWorld")) {
12  auto core_log = logging::core::get();
13 
15  "FieldEvaluatorWorld"));
17  "FieldEvaluatorSync"));
19  "FieldEvaluatorSelf"));
20 
21  LogManager::setLog("FieldEvaluatorWorld");
22  LogManager::setLog("FieldEvaluatorSync");
23  LogManager::setLog("FieldEvaluatorSelf");
24 
25  MOFEM_LOG_TAG("FieldEvaluatorWorld", "FieldEvaluator");
26  MOFEM_LOG_TAG("FieldEvaluatorSync", "FieldEvaluator");
27  MOFEM_LOG_TAG("FieldEvaluatorSelf", "FieldEvaluator");
28  }
29 
30  MOFEM_LOG("FieldEvaluatorWorld", Sev::noisy) << "Field evaluator intreface";
31 }

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 42 of file FieldEvaluator.cpp.

43  {
44  CoreInterface &m_field = cOre;
46  EntityHandle fe_meshset;
47  fe_meshset = m_field.get_finite_element_meshset(finite_element);
48  Range entities;
49  CHKERR m_field.get_moab().get_entities_by_dimension(fe_meshset, D, entities,
50  true);
51  spd_ptr->treePtr.reset(new AdaptiveKDTree(&m_field.get_moab()));
52  CHKERR spd_ptr->treePtr->build_tree(entities, &spd_ptr->rooTreeSet);
54 }

◆ 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 63 of file FieldEvaluator.cpp.

64  {
65  return buildTree<2>(spd_ptr, finite_element);
66 }

◆ 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 57 of file FieldEvaluator.cpp.

58  {
59  return buildTree<3>(spd_ptr, finite_element);
60 }

◆ 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);
}
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 161 of file FieldEvaluator.cpp.

165  {
166  CoreInterface &m_field = cOre;
168 
169  std::vector<EntityHandle> leafs_out;
170  CHKERR data_ptr->treePtr->distance_search(point, distance, leafs_out);
171  Range tree_ents;
172  for (auto lit : leafs_out)//{
173  CHKERR m_field.get_moab().get_entities_by_dimension(lit, D, tree_ents,
174  true);
175 
176  if (verb >= NOISY)
177  MOFEM_LOG("SELF", Sev::noisy) << "tree entities: " << tree_ents;
178 
179  data_ptr->evalPointEntityHandle.resize(data_ptr->nbEvalPoints);
180  std::fill(data_ptr->evalPointEntityHandle.begin(),
181  data_ptr->evalPointEntityHandle.end(), 0);
182 
183  std::vector<double> local_coords;
184  std::vector<double> shape;
185 
186  auto nb_eval_points = data_ptr->nbEvalPoints;
187 
188  for (auto tet : tree_ents) {
189  const EntityHandle *conn;
190  int num_nodes;
191  CHKERR m_field.get_moab().get_connectivity(tet, conn, num_nodes, true);
192 
193  if constexpr (D == 3) {
194 
195  local_coords.resize(3 * nb_eval_points);
196  shape.resize(4 * nb_eval_points);
197 
198  std::array<double, 12> coords;
199  CHKERR m_field.get_moab().get_coords(conn, num_nodes, coords.data());
200 
202  coords.data(), &data_ptr->evalPoints[0], nb_eval_points,
203  &local_coords[0]);
204  CHKERR Tools::shapeFunMBTET<3>(&shape[0], &local_coords[0],
205  &local_coords[1], &local_coords[2],
206  nb_eval_points);
207 
210  &shape[0], &shape[1], &shape[2], &shape[3]};
212  &data_ptr->shapeFunctions(0, 0), &data_ptr->shapeFunctions(0, 1),
213  &data_ptr->shapeFunctions(0, 2), &data_ptr->shapeFunctions(0, 3)};
214 
217  &local_coords[0], &local_coords[1], &local_coords[2]};
219  &(data_ptr->localCoords(0, 0)), &(data_ptr->localCoords(0, 1)),
220  &(data_ptr->localCoords(0, 2))};
221 
222  for (int n = 0; n != nb_eval_points; ++n) {
223 
224  const double eps = data_ptr->eps;
225  if (t_shape(0) >= 0 - eps && t_shape(0) <= 1 + eps &&
226 
227  t_shape(1) >= 0 - eps && t_shape(1) <= 1 + eps &&
228 
229  t_shape(2) >= 0 - eps && t_shape(2) <= 1 + eps &&
230 
231  t_shape(3) >= 0 - eps && t_shape(3) <= 1 + eps) {
232 
233  data_ptr->evalPointEntityHandle[n] = tet;
234  t_shape_data(i4) = t_shape(i4);
235  t_local_data(j3) = t_local(j3);
236  }
237 
238  ++t_shape;
239  ++t_shape_data;
240  ++t_local;
241  ++t_local_data;
242  }
243  }
244 
245  if constexpr (D == 2) {
246 
247  local_coords.resize(2 * nb_eval_points);
248  shape.resize(3 * nb_eval_points);
249 
250  std::array<double, 9> coords;
251  CHKERR m_field.get_moab().get_coords(conn, num_nodes, coords.data());
252 
254  coords.data(), &data_ptr->evalPoints[0], nb_eval_points,
255  &local_coords[0]);
256  CHKERR Tools::shapeFunMBTRI<2>(&shape[0], &local_coords[0],
257  &local_coords[1], nb_eval_points);
258 
261  &shape[0], &shape[1], &shape[2]};
263  &data_ptr->shapeFunctions(0, 0), &data_ptr->shapeFunctions(0, 1),
264  &data_ptr->shapeFunctions(0, 2)};
265 
268  &local_coords[0], &local_coords[1]};
269  data_ptr->localCoords.resize(nb_eval_points, 2);
271  &(data_ptr->localCoords(0, 0)), &(data_ptr->localCoords(0, 1))};
272 
273  for (int n = 0; n != nb_eval_points; ++n) {
274 
275  const double eps = data_ptr->eps;
276  if (t_shape(0) >= 0 - eps && t_shape(0) <= 1 + eps &&
277 
278  t_shape(1) >= 0 - eps && t_shape(1) <= 1 + eps &&
279 
280  t_shape(2) >= 0 - eps && t_shape(2) <= 1 + eps) {
281 
282  data_ptr->evalPointEntityHandle[n] = tet;
283  t_shape_data(i3) = t_shape(i3);
284  t_local_data(j2) = t_local(j2);
285  }
286 
287  ++t_shape;
288  ++t_shape_data;
289  ++t_local;
290  ++t_local_data;
291  }
292  }
293  }
294 
295  const Problem *prb_ptr;
296  CHKERR m_field.get_problem(problem, &prb_ptr);
297  boost::shared_ptr<NumeredEntFiniteElement_multiIndex> numered_fes(
299 
300  Range in_tets;
301  in_tets.insert_list(data_ptr->evalPointEntityHandle.begin(),
302  data_ptr->evalPointEntityHandle.end());
303  in_tets = in_tets.subset_by_dimension(D);
304 
305  if (verb >= NOISY)
306  MOFEM_LOG("SELF", Sev::noisy) << "in tets: " << in_tets << endl;
307 
308  auto fe_miit =
309  m_field.get_finite_elements()->get<FiniteElement_name_mi_tag>().find(
310  finite_element);
311  if (fe_miit !=
312  m_field.get_finite_elements()->get<FiniteElement_name_mi_tag>().end()) {
313 
314  for (auto peit = in_tets.pair_begin(); peit != in_tets.pair_end(); ++peit) {
315 
316  auto lo =
317  prb_ptr->numeredFiniteElementsPtr->get<Unique_mi_tag>().lower_bound(
319  peit->first, (*fe_miit)->getFEUId()));
320  auto hi =
321  prb_ptr->numeredFiniteElementsPtr->get<Unique_mi_tag>().upper_bound(
323  peit->second, (*fe_miit)->getFEUId()));
324  numered_fes->insert(lo, hi);
325 
326  if (verb >= VERY_NOISY)
327  std::cerr << "numered elements:" << std::endl;
328  for (; lo != hi; ++lo)
329  if (verb >= VERY_NOISY)
330  std::cerr << **lo << endl;
331  }
332  if (verb >= VERY_NOISY)
333  std::cerr << std::endl;
334 
335  if (auto fe_ptr = data_ptr->feMethodPtr.lock()) {
336 
337  MOFEM_LOG_C("FieldEvaluatorSync", Sev::verbose,
338  "Number elements %d to evaluate at proc %d",
339  numered_fes->size(), m_field.get_comm_rank());
340  MOFEM_LOG_SYNCHRONISE(m_field.get_comm());
341 
342  if (!cache_ptr) {
343  cache_ptr = boost::make_shared<CacheTuple>();
344  CHKERR m_field.cache_problem_entities(prb_ptr->getName(), cache_ptr);
345 
346  MOFEM_LOG("FieldEvaluatorSync", Sev::noisy)
347  << "If you call that function many times in the loop consider to "
348  "set "
349  "cache_ptr outside of the loop. Otherwise code can be slow.";
350  }
351 
352  CHKERR m_field.loop_finite_elements(prb_ptr, finite_element, *fe_ptr,
353  lower_rank, upper_rank, numered_fes,
354  bh, cache_ptr, verb);
355 
356  } else
357  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
358  "Pointer to element does not exists");
359  }
360 
362 }

◆ 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 373 of file FieldEvaluator.cpp.

377  {
378  return evalFEAtThePoint<2>(point, distance, problem, finite_element, data_ptr,
379  lower_rank, upper_rank, cache_ptr, bh, verb);
380 }

◆ 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 364 of file FieldEvaluator.cpp.

368  {
369  return evalFEAtThePoint<3>(point, distance, problem, finite_element, data_ptr,
370  lower_rank, upper_rank, cache_ptr, bh, verb);
371 }

◆ 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 
)
inline

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
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 
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  }

◆ query_interface()

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

Implements MoFEM::UnknownInterface.

Definition at line 33 of file FieldEvaluator.cpp.

34  {
36  *iface = const_cast<FieldEvaluatorInterface *>(this);
38 }

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:
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
MoFEM::LogManager::checkIfChannelExist
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
Definition: LogManager.cpp:404
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
FTensor::Tensor1
Definition: Tensor1_value.hpp:8
EntityHandle
NOISY
@ NOISY
Definition: definitions.h:224
VERY_NOISY
@ VERY_NOISY
Definition: definitions.h:225
MoFEM::FieldEvaluatorInterface::FieldEvaluatorInterface
FieldEvaluatorInterface(const MoFEM::Core &core)
Definition: FieldEvaluator.cpp:8
MoFEM::LogManager::createSink
static boost::shared_ptr< SinkType > createSink(boost::shared_ptr< std::ostream > stream_ptr, std::string comm_filter)
Create a sink object.
Definition: LogManager.cpp:298
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
MoFEM::Tools::getLocalCoordinatesOnReferenceThreeNodeTri
static MoFEMErrorCode getLocalCoordinatesOnReferenceThreeNodeTri(const double *elem_coords, const double *glob_coords, const int nb_nodes, double *local_coords)
Get the local coordinates on reference three node tri object.
Definition: Tools.cpp:188
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MOFEM_LOG_C
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:311
MoFEM::EntFiniteElement::getLocalUniqueIdCalculate
UId getLocalUniqueIdCalculate() const
Generate UId for finite element entity.
Definition: FEMultiIndices.hpp:528
MOFEM_LOG_SYNCHRONISE
#define MOFEM_LOG_SYNCHRONISE(comm)
Synchronise "SYNC" channel.
Definition: LogManager.hpp:345
MoFEM::LogManager::getStrmSync
static boost::shared_ptr< std::ostream > getStrmSync()
Get the strm sync object.
Definition: LogManager.cpp:348
MoFEM::FieldEvaluatorInterface::cOre
MoFEM::Core & cOre
Definition: FieldEvaluator.hpp:26
MoFEM::FieldEvaluatorInterface::evalFEAtThePoint3D
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.
Definition: FieldEvaluator.cpp:364
MoFEM::LogManager::getStrmWorld
static boost::shared_ptr< std::ostream > getStrmWorld()
Get the strm world object.
Definition: LogManager.cpp:344
MoFEM::LogManager::getStrmSelf
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Definition: LogManager.cpp:340
MoFEM::Tools::getLocalCoordinatesOnReferenceFourNodeTet
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:88
MOFEM_LOG_TAG
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:339
MoFEM::VolumeElementForcesAndSourcesCore
Volume finite element base.
Definition: VolumeElementForcesAndSourcesCore.hpp:26
NumeredEntFiniteElement_multiIndex
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< Ent_mi_tag >, 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.
Definition: FEMultiIndices.hpp:830
FTensor::Index
Definition: Index.hpp:23
convert.n
n
Definition: convert.py:82
Range
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
eps
static const double eps
Definition: check_base_functions_derivatives_on_tet.cpp:11
MyOp
Operator used to check consistency between local coordinates and global cooridnates for integrated po...
Definition: field_evaluator.cpp:21
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
ReactionDiffusionEquation::D
const double D
diffusivity
Definition: reaction_diffusion.cpp:20
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
MoFEM::LogManager::setLog
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:389
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359