|
| v0.14.0
|
Go to the documentation of this file.
10 #ifndef __FIELD_EVALUATOR_HPP__
11 #define __FIELD_EVALUATOR_HPP__
42 SetPtsData(boost::shared_ptr<MoFEM::ForcesAndSourcesCore> fe_method_ptr,
43 const double *eval_points,
const int nb_eval_points,
51 inline void setEvalPoints(
const double *ptr,
const int nb_eval_points) {
69 boost::scoped_ptr<AdaptiveKDTree>
treePtr;
80 int order_col,
int order_data);
105 template <
typename VE,
typename SPD = SetPtsData,
typename SP = SetPts>
106 boost::shared_ptr<SPD>
107 getData(
const double *ptr =
nullptr,
const int nb_eval_points = 0,
110 boost::scoped_ptr<VE> elePtr;
111 boost::scoped_ptr<SPD> setPtsDataPtr;
112 boost::scoped_ptr<SP> setPtsPtr;
114 boost::shared_ptr<PackData> pack_data(
new PackData());
116 pack_data->elePtr.reset(
new VE(m_field));
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());
135 const std::string finite_element);
141 const std::string finite_element);
187 const std::string problem,
188 const std::string finite_element,
189 boost::shared_ptr<SetPtsData> data_ptr,
int lower_rank,
190 int upper_rank, boost::shared_ptr<CacheTuple> cache_ptr,
198 const std::string problem,
199 const std::string finite_element,
200 boost::shared_ptr<SetPtsData> data_ptr,
int lower_rank,
201 int upper_rank, boost::shared_ptr<CacheTuple> cache_ptr,
210 const std::string finite_element);
218 const std::string problem,
219 const std::string finite_element,
220 boost::shared_ptr<SetPtsData> data_ptr,
int lower_rank,
221 int upper_rank, boost::shared_ptr<CacheTuple> cache_ptr,
228 #endif // __FIELD_EVALUATOR_HPP__
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
MoFEMErrorCode buildTree(boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
Build spatial tree.
MoFEMErrorCode operator()(ForcesAndSourcesCore *fe_raw_ptr, int order_row, int order_col, int order_data)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
FieldEvaluatorInterface(const MoFEM::Core &core)
Deprecated interface functions.
MatrixDouble shapeFunctions
boost::weak_ptr< SetPtsData > dataPtr
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.
MoFEMErrorCode buildTree3D(boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
Build spatial tree.
Default evaluator for setting integration points.
implementation of Data Operators for Forces and Sources
VERBOSITY_LEVELS
Verbosity levels.
MoFEMErrorCode buildTree2D(boost::shared_ptr< SetPtsData > spd_ptr, const std::string finite_element)
Build spatial tree.
Field evaluator interface.
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.
SetPts(boost::shared_ptr< SetPtsData > data_ptr)
std::vector< EntityHandle > evalPointEntityHandle
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.
base class for all interface classes
structure to get information form mofem into EntitiesFieldData
const double * evalPoints
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.
MoFEMTypes
Those types control how functions respond on arguments, f.e. error handling.
void setEvalPoints(const double *ptr, const int nb_eval_points)
boost::weak_ptr< MoFEM::ForcesAndSourcesCore > feMethodPtr
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.
boost::scoped_ptr< AdaptiveKDTree > treePtr