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