v0.14.0
FieldEvaluator.cpp
Go to the documentation of this file.
1 /** \file FieldEvaluator.cpp
2  * \brief Field evaluator
3  *
4  */
5 
6 namespace MoFEM {
7 
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 }
32 
34  boost::typeindex::type_index type_index, UnknownInterface **iface) const {
36  *iface = const_cast<FieldEvaluatorInterface *>(this);
38 }
39 
40 template <int D>
42 FieldEvaluatorInterface::buildTree(boost::shared_ptr<SetPtsData> spd_ptr,
43  const std::string finite_element) {
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 }
55 
57 FieldEvaluatorInterface::buildTree3D(boost::shared_ptr<SetPtsData> spd_ptr,
58  const std::string finite_element) {
59  return buildTree<3>(spd_ptr, finite_element);
60 }
61 
63 FieldEvaluatorInterface::buildTree2D(boost::shared_ptr<SetPtsData> spd_ptr,
64  const std::string finite_element) {
65  return buildTree<2>(spd_ptr, finite_element);
66 }
67 
70  int order_row, int order_col,
71  int order_data) {
73 
74  if (auto data_ptr = dataPtr.lock()) {
75 
76  const auto nb_eval_points = data_ptr->nbEvalPoints;
77  const auto &shape_functions = data_ptr->shapeFunctions;
78  const auto &eval_pointentity_handle = data_ptr->evalPointEntityHandle;
79 
80  if (auto fe_ptr = data_ptr->feMethodPtr.lock()) {
81 
82  auto &fe = static_cast<ForcesAndSourcesCore &>(*fe_ptr);
83 
84 #ifndef NDEBUG
85  if (fe_ptr.get() != fe_raw_ptr)
86  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong FE ptr");
87 #endif
88 
89  const auto fe_ent = fe.numeredEntFiniteElementPtr->getEnt();
90  const auto fe_type = type_from_handle(fe_ent);
91  const auto fe_dim = moab::CN::Dimension(fe_type);
92 
93  auto &gauss_pts = fe.gaussPts;
94  int nb_gauss_pts;
95 
96  if (fe_dim == 3) {
97  gauss_pts.resize(4, nb_eval_points, false);
99  &shape_functions(0, 0), &shape_functions(0, 1),
100  &shape_functions(0, 2), &shape_functions(0, 3)};
102  &gauss_pts(0, 0), &gauss_pts(1, 0), &gauss_pts(2, 0)};
103 
104  nb_gauss_pts = 0;
105  for (int nn = 0; nn != nb_eval_points; ++nn) {
106  if (eval_pointentity_handle[nn] == fe_ent) {
107  for (const int i : {0, 1, 2}) {
108  t_gauss_pts(i) = t_shape(i + 1);
109  }
110  gauss_pts(3, nb_gauss_pts) = nn;
111  ++t_gauss_pts;
112  ++nb_gauss_pts;
113  }
114  ++t_shape;
115  }
116  gauss_pts.resize(4, nb_gauss_pts, true);
117  } else if (fe_dim == 2) {
118  gauss_pts.resize(3, nb_eval_points, false);
120  &shape_functions(0, 0), &shape_functions(0, 1),
121  &shape_functions(0, 2)};
123  &gauss_pts(0, 0), &gauss_pts(1, 0)};
124  nb_gauss_pts = 0;
125  for (int nn = 0; nn != nb_eval_points; ++nn) {
126  if (eval_pointentity_handle[nn] == fe_ent) {
127  for (const int i : {0, 1}) {
128  t_gauss_pts(i) = t_shape(i + 1);
129  }
130  gauss_pts(2, nb_gauss_pts) = nn;
131  ++t_gauss_pts;
132  ++nb_gauss_pts;
133  }
134  ++t_shape;
135  }
136  gauss_pts.resize(3, nb_gauss_pts, true);
137  } else {
138  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
139  "Dimension not implemented");
140  }
141 
142 #ifndef NDEBUG
143  MOFEM_LOG("SELF", Sev::noisy)
144  << "nbEvalOPoints / nbGau_sspt_s: " << nb_eval_points << " / "
145  << nb_gauss_pts;
146  MOFEM_LOG("SELF", Sev::noisy) << "gauss pts: " << gauss_pts;
147 #endif
148 
149  } else
150  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
151  "Pointer to element does not exists");
152 
153  } else
154  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
155  "Pointer to data does not exists");
156 
158 }
159 
160 template <int D>
162  const double *const point, const double distance, const std::string problem,
163  const std::string finite_element, boost::shared_ptr<SetPtsData> data_ptr,
164  int lower_rank, int upper_rank, boost::shared_ptr<CacheTuple> cache_ptr,
165  MoFEMTypes bh, VERBOSITY_LEVELS verb) {
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 }
363 
365  const double *const point, const double distance, const std::string problem,
366  const std::string finite_element, boost::shared_ptr<SetPtsData> data_ptr,
367  int lower_rank, int upper_rank, boost::shared_ptr<CacheTuple> cache_ptr,
368  MoFEMTypes bh, VERBOSITY_LEVELS verb) {
369  return evalFEAtThePoint<3>(point, distance, problem, finite_element, data_ptr,
370  lower_rank, upper_rank, cache_ptr, bh, verb);
371 }
372 
374  const double *const point, const double distance, const std::string problem,
375  const std::string finite_element, boost::shared_ptr<SetPtsData> data_ptr,
376  int lower_rank, int upper_rank, boost::shared_ptr<CacheTuple> cache_ptr,
377  MoFEMTypes bh, VERBOSITY_LEVELS verb) {
378  return evalFEAtThePoint<2>(point, distance, problem, finite_element, data_ptr,
379  lower_rank, upper_rank, cache_ptr, bh, verb);
380 }
381 
382 } // namespace MoFEM
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
MoFEM::LogManager::checkIfChannelExist
static bool checkIfChannelExist(const std::string channel)
Check if channel exist.
Definition: LogManager.cpp:404
MoFEM::CoreInterface::loop_finite_elements
virtual MoFEMErrorCode loop_finite_elements(const std::string problem_name, const std::string &fe_name, FEMethod &method, boost::shared_ptr< NumeredEntFiniteElement_multiIndex > fe_ptr=nullptr, MoFEMTypes bh=MF_EXIST, CacheTupleWeakPtr cache_ptr=CacheTupleSharedPtr(), int verb=DEFAULT_VERBOSITY)=0
Make a loop over finite elements.
MoFEM::CoreInterface::get_finite_element_meshset
virtual EntityHandle get_finite_element_meshset(const std::string name) const =0
MoFEM::CoreTmp< 0 >
Core (interface) class.
Definition: Core.hpp:82
FTensor::Tensor1
Definition: Tensor1_value.hpp:8
MoFEM::FieldEvaluatorInterface::query_interface
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: FieldEvaluator.cpp:33
EntityHandle
NOISY
@ NOISY
Definition: definitions.h:224
MoFEM::CoreInterface::get_comm
virtual MPI_Comm & get_comm() const =0
MoFEM::FieldEvaluatorInterface::buildTree
MoFEMErrorCode buildTree(boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
Build spatial tree.
Definition: FieldEvaluator.cpp:42
MoFEM::FieldEvaluatorInterface::SetPts::operator()
MoFEMErrorCode operator()(ForcesAndSourcesCore *fe_raw_ptr, int order_row, int order_col, int order_data)
Definition: FieldEvaluator.cpp:69
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::CoreInterface::get_comm_rank
virtual int get_comm_rank() const =0
VERY_NOISY
@ VERY_NOISY
Definition: definitions.h:225
MoFEM::FieldEvaluatorInterface::FieldEvaluatorInterface
FieldEvaluatorInterface(const MoFEM::Core &core)
Definition: FieldEvaluator.cpp:8
MoFEM::CoreInterface::get_problem
virtual const Problem * get_problem(const std::string problem_name) const =0
Get the problem object.
MoFEM::CoreInterface::get_finite_elements
virtual const FiniteElement_multiIndex * get_finite_elements() const =0
Get the finite elements object.
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::FieldEvaluatorInterface::SetPts::dataPtr
boost::weak_ptr< SetPtsData > dataPtr
Definition: FieldEvaluator.hpp:83
MoFEM::FieldEvaluatorInterface::buildTree3D
MoFEMErrorCode buildTree3D(boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
Build spatial tree.
Definition: FieldEvaluator.cpp:57
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::CoreInterface::get_moab
virtual moab::Interface & get_moab()=0
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
VERBOSITY_LEVELS
VERBOSITY_LEVELS
Verbosity levels.
Definition: definitions.h:219
MoFEM::FieldEvaluatorInterface::buildTree2D
MoFEMErrorCode buildTree2D(boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
Build spatial tree.
Definition: FieldEvaluator.cpp:63
MOFEM_LOG_C
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:311
MoFEM::FieldEvaluatorInterface
Field evaluator interface.
Definition: FieldEvaluator.hpp:21
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::CoreInterface::cache_problem_entities
virtual MoFEMErrorCode cache_problem_entities(const std::string prb_name, CacheTupleWeakPtr cache_ptr)=0
Cache variables.
MoFEM::LogManager::getStrmSelf
static boost::shared_ptr< std::ostream > getStrmSelf()
Get the strm self object.
Definition: LogManager.cpp:340
MoFEM::type_from_handle
auto type_from_handle(const EntityHandle h)
get type from entity handle
Definition: Templates.hpp:1898
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
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
MoFEM::Problem::numeredFiniteElementsPtr
boost::shared_ptr< NumeredEntFiniteElement_multiIndex > numeredFiniteElementsPtr
store finite elements
Definition: ProblemsMultiIndices.hpp:77
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
MoFEM::UnknownInterface
base class for all interface classes
Definition: UnknownInterface.hpp:34
MoFEM::ForcesAndSourcesCore
structure to get information form mofem into EntitiesFieldData
Definition: ForcesAndSourcesCore.hpp:22
Range
MOFEM_LOG
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
MoFEM::FiniteElement_name_mi_tag
Definition: TagMultiIndices.hpp:26
eps
static const double eps
Definition: check_base_functions_derivatives_on_tet.cpp:11
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MoFEM::FieldEvaluatorInterface::evalFEAtThePoint
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.
Definition: FieldEvaluator.cpp:161
MoFEMTypes
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:110
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::Unique_mi_tag
Definition: TagMultiIndices.hpp:18
MoFEM::Problem
keeps basic data about problem
Definition: ProblemsMultiIndices.hpp:54
MoFEM::FieldEvaluatorInterface::evalFEAtThePoint2D
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.
Definition: FieldEvaluator.cpp:373
MoFEM::CoreInterface
Interface.
Definition: Interface.hpp:27
MoFEM::LogManager::setLog
static LoggerType & setLog(const std::string channel)
Set ans resset chanel logger.
Definition: LogManager.cpp:389
MoFEM::Problem::getName
auto getName() const
Definition: ProblemsMultiIndices.hpp:372
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359