v0.11.1
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 /*
11  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
12  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
14  * License for more details.
15  *
16  * You should have received a copy of the GNU Lesser General Public
17  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
18  */
19 
20 #ifndef __FIELD_EVALUATOR_HPP__
21 #define __FIELD_EVALUATOR_HPP__
22 
23 #include "UnknownInterface.hpp"
24 
25 namespace MoFEM {
26 
29 
30 /** \brief Field evaluator interface
31 
32  * \ingroup field_evaluator
33  */
35 
37  UnknownInterface **iface) const;
38 
41 
42  struct SetPtsData {
43 
44  /**
45  * @brief Set the Gauss Pts data
46  *
47  * @param fe_method_ptr pointer to finite element instance
48  * @param eval_points pointer to array with evaluation points
49  * @param nb_eval_points number of evaluated points
50  * @param eps tolerance used to find if point is in the element
51  * @param verb
52  */
53  SetPtsData(boost::shared_ptr<MoFEM::ForcesAndSourcesCore> fe_method_ptr,
54  const double *eval_points, const int nb_eval_points,
55  const double eps, VERBOSITY_LEVELS verb = QUIET)
56  : feMethodPtr(fe_method_ptr), evalPoints(eval_points),
57  nbEvalPoints(nb_eval_points), eps(eps), verb(verb) {
58  localCoords.resize(nbEvalPoints, 3);
59  shapeFunctions.resize(nbEvalPoints, 4);
60  }
61 
62  inline void setEvalPoints(const double *ptr, const int nb_eval_points) {
63  evalPoints = ptr;
64  nbEvalPoints = nb_eval_points;
65  localCoords.resize(nbEvalPoints, 3, false);
66  shapeFunctions.resize(nbEvalPoints, 4, false);
67  }
68 
69  boost::weak_ptr<MoFEM::ForcesAndSourcesCore> feMethodPtr;
70  const double *evalPoints;
72  double eps;
74 
77  std::vector<EntityHandle> evalPointEntityHandle;
78 
79  EntityHandle rooTreeSet;
80  boost::scoped_ptr<AdaptiveKDTree> treePtr;
81  };
82 
83  /**
84  * @brief Default evaluator for setting integration points
85  *
86  */
87  struct SetPts {
88  SetPts() = delete;
89  SetPts(boost::shared_ptr<SetPtsData> data_ptr) : dataPtr(data_ptr) {}
90  MoFEMErrorCode operator()(int order_row, int order_col, int order_data);
91 
92  private:
93  boost::weak_ptr<SetPtsData> dataPtr;
94  };
95 
96  /**
97  * @brief Get the Data object
98  *
99  * Pack pointers with data structures for field evaluator and finite
100  * element. Function return shared pointer if returned shared pointer
101  * is reset; all data are destroyed. The idea is to pack all data in
102  * one structure, create shared pointer to it and return aliased shared
103  * pointer to one of the elements of the data structure. It is a bit
104  * complicated, but has great flexibility.
105  *
106  * @tparam VE
107  * @tparam SetPtsData
108  * @tparam SetPts
109  * @param ptr
110  * @param nb_eval_points
111  * @param eps
112  * @param verb
113  * @return boost::shared_ptr<SPD>
114  */
115  template <typename VE, typename SPD = SetPtsData, typename SP = SetPts>
116  boost::shared_ptr<SPD>
117  getData(const double *ptr = nullptr, const int nb_eval_points = 0,
118  const double eps = 1e-12, VERBOSITY_LEVELS verb = QUIET) {
119  struct PackData {
120  boost::scoped_ptr<VE> volElePtr;
121  boost::scoped_ptr<SPD> setPtsDataPtr;
122  boost::scoped_ptr<SP> setPtsPtr;
123  };
124  boost::shared_ptr<PackData> pack_data(new PackData());
125  MoFEM::Interface &m_field = cOre;
126  pack_data->volElePtr.reset(new VE(m_field));
127  pack_data->setPtsDataPtr.reset(
128  new SPD(boost::shared_ptr<VE>(pack_data, pack_data->volElePtr.get()),
129  ptr, nb_eval_points, eps, verb));
130  pack_data->setPtsPtr.reset(new SP(
131  boost::shared_ptr<SPD>(pack_data, pack_data->setPtsDataPtr.get())));
132  pack_data->volElePtr->setRuleHook = boost::ref(*pack_data->setPtsPtr);
133  boost::shared_ptr<SPD> data(pack_data, pack_data->setPtsDataPtr.get());
134  return data;
135  }
136 
137  /**
138  * @brief Build spatial tree
139  *
140  * @param finite_element finite element name
141  * @return MoFEMErrorCode
142  */
143  MoFEMErrorCode buildTree3D(boost::shared_ptr<SetPtsData> spd_ptr,
144  const std::string finite_element);
145 
146  /**
147  * @brief Evaluate field at artbitray position
148  *
149  * \code
150 
151  std::array<double, 3> point = {0, 0, 0};
152  const double dist = 0.3;
153  std::array<double, 6> eval_points = {-1., -1., -1., 1., 1., 1. };
154 
155  using VolEle = VolumeElementForcesAndSourcesCore;
156  auto data_ptr =
157  m_field.getInterface<FieldEvaluatorInterface>()->
158  getData<VolEle>(point.data(), point.size/3);
159 
160  if(auto vol_ele = data_ptr->feMethod.lock()) {
161  // push operators to finite element instance, e.g.
162  vol_ele->getOpPtrVector().push_back(new MyOp());
163  // iterate over elemnts with evaluated points
164  auto cache_ptr = boost::make_shared<CacheTuple>();
165  CHKERR m_field.cache_problem_entities(prb_ptr->getName(), cache_ptr);
166  CHKERR m_field.getInterface<FieldEvaluatorInterface>()
167  ->evalFEAtThePoint3D(point.data(), dist, prb_ptr->getName(),
168  "FINITE_ELEMENT_NAME",
169  data_ptr, m_field.get_comm_rank(),
170  m_field.get_comm_rank(), cache_ptr);
171  }
172 
173  * \endcode
174  *
175  * @param point point used to find tetrahedrons
176  * @param distance distance from the point where tetrahedrons are searched
177  * @param problem problem name
178  * @param finite_element finite element name
179  * @param data_ptr pointer to data abut gauss points
180  * @param lower_rank lower processor
181  * @param upper_rank upper process
182  * @param cache_ptr cache or problem entities
183  * @param bh control looping over entities, e.g. throwing error if element not
184  found
185  * @param verb verbosity level
186  * @return MoFEMErrorCode
187  */
189  evalFEAtThePoint3D(const double *const point, const double distance,
190  const std::string problem,
191  const std::string finite_element,
192  boost::shared_ptr<SetPtsData> data_ptr, int lower_rank,
193  int upper_rank, boost::shared_ptr<CacheTuple> cache_ptr,
195 };
196 
197 } // namespace MoFEM
198 
199 #endif // __FIELD_EVALUATOR_HPP__
200 
201 /**
202  * \defgroup field_evaluator Field Evaluator
203  * \brief Evaluate field at the point
204  *
205  * \ingroup mofem
206  */
MoFEM interface.
VERBOSITY_LEVELS
Verbosity levels.
Definition: definitions.h:275
@ QUIET
Definition: definitions.h:277
@ FIELDEVALUATOR_INTERFACE
Definition: definitions.h:71
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
Definition: definitions.h:189
@ MF_EXIST
Definition: definitions.h:189
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
static const MOFEMuuid IDD_MOFEMFieldEvaluator
Core (interface) class.
Definition: Core.hpp:77
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()(int order_row, int order_col, int order_data)
SetPts(boost::shared_ptr< SetPtsData > data_ptr)
boost::weak_ptr< SetPtsData > dataPtr
Field evaluator interface.
MoFEMErrorCode buildTree3D(boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
Build spatial tree.
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
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.
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.
MoFEM interface unique ID.
base class for all interface classes