21 #include "libmesh/mesh_function.h" 22 #include "libmesh/dense_vector.h" 23 #include "libmesh/equation_systems.h" 24 #include "libmesh/numeric_vector.h" 25 #include "libmesh/dof_map.h" 26 #include "libmesh/point_locator_tree.h" 27 #include "libmesh/fe_base.h" 28 #include "libmesh/fe_interface.h" 29 #include "libmesh/fe_compute_data.h" 30 #include "libmesh/mesh_base.h" 31 #include "libmesh/point.h" 32 #include "libmesh/elem.h" 33 #include "libmesh/int_range.h" 34 #include "libmesh/fe_map.h" 45 std::vector<unsigned int> vars,
49 _eqn_systems (eqn_systems),
52 _system_vars (
std::move(vars)),
53 _out_of_mesh_mode (false)
62 const unsigned int var,
66 _eqn_systems (eqn_systems),
70 _out_of_mesh_mode (false)
77 _eqn_systems (mf._eqn_systems),
79 _dof_map (mf._dof_map),
80 _system_vars (mf._system_vars),
81 _out_of_mesh_mode (mf._out_of_mesh_mode)
96 std::make_unique<std::set<subdomain_id_type>>
127 #ifdef LIBMESH_ENABLE_DEPRECATED 130 libmesh_deprecated();
137 #endif // LIBMESH_ENABLE_DEPRECATED 155 return std::make_unique<MeshFunction>(*this);
175 std::map<const Elem *, DenseVector<Number>> buffer;
177 std::map<const Elem *, Number> return_value;
178 for (
const auto & [elem, vec] : buffer)
179 return_value[elem] = vec(0);
190 std::vector<Gradient> buf (1);
200 std::map<const Elem *, std::vector<Gradient>> buffer;
202 std::map<const Elem *, Gradient> return_value;
203 for (
const auto & [elem, vec] : buffer)
204 return_value[elem] = vec[0];
210 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 216 std::vector<Tensor> buf (1);
232 const std::set<subdomain_id_type> * subdomain_ids)
249 output.
resize (cast_int<unsigned int>
254 const unsigned int dim = element->
dim();
287 std::vector<dof_id_type> dof_indices;
297 output(index) =
value;
321 const std::set<subdomain_id_type> * subdomain_ids)
329 std::set<const Elem *> candidate_element = this->
find_elements(p,subdomain_ids);
333 for (
const auto & element : candidate_element)
335 const unsigned int dim = element->dim();
369 std::vector<dof_id_type> dof_indices;
379 temp_output(index) =
value;
388 output[element] = temp_output;
396 std::vector<Gradient> & output)
405 std::vector<Gradient> & output,
406 const std::set<subdomain_id_type> * subdomain_ids)
421 std::map<
const Elem *, std::vector<Gradient>> & output)
430 std::map<
const Elem *, std::vector<Gradient>> & output,
431 const std::set<subdomain_id_type> * subdomain_ids)
439 std::set<const Elem *> candidate_element = this->
find_elements(p,subdomain_ids);
443 for (
const auto & element : candidate_element)
446 std::vector<Gradient> temp_output (cast_int<unsigned int>(this->
_system_vars.size()));
451 output.emplace(element, std::move(temp_output));
458 const Elem * element,
459 std::vector<Gradient> & output)
467 const unsigned int dim = element->
dim();
475 std::vector<Point> point_list (1, mapped_point);
494 std::vector<dof_id_type> dof_indices;
506 const std::vector<std::vector<RealGradient>> & dphi = point_fe->get_dphi();
507 point_fe->reinit(element, &point_list);
513 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 529 for (std::size_t v=0; v<
dim; v++)
530 for (std::size_t xyz=0; xyz<LIBMESH_DIM; xyz++)
535 * this->
_vector(dof_indices[i]);
540 output[index] = grad;
545 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 548 std::vector<Tensor> & output)
557 std::vector<Tensor> & output,
558 const std::set<subdomain_id_type> * subdomain_ids)
571 libmesh_warning(
"Warning: Requested the Hessian of an Infinite element." 572 <<
"Second derivatives for Infinite elements" 573 <<
" are not yet implemented!" 582 const unsigned int dim = element->
dim();
591 std::vector<Point> point_list (1, mapped_point);
609 const std::vector<std::vector<RealTensor>> & d2phi =
610 point_fe->get_d2phi();
611 point_fe->reinit(element, &point_list);
614 std::vector<dof_id_type> dof_indices;
623 output[index] = hess;
631 const std::set<subdomain_id_type> * subdomain_ids)
const 642 cast_ptr<const MeshFunction *>(this->
_master);
644 "ERROR: If you use out-of-mesh-mode in connection with master mesh " 645 "functions, you must enable out-of-mesh mode for both the master and the slave mesh function.");
650 const Elem * element = (*_point_locator)(p, subdomain_ids);
657 const std::set<subdomain_id_type> * subdomain_ids)
const 668 cast_ptr<const MeshFunction *>(this->
_master);
670 "ERROR: If you use out-of-mesh-mode in connection with master mesh " 671 "functions, you must enable out-of-mesh mode for both the master and the slave mesh function.");
676 std::set<const Elem *> candidate_elements;
677 std::set<const Elem *> final_candidate_elements;
678 (*_point_locator)(p, candidate_elements, subdomain_ids);
682 for (
const auto & element : candidate_elements)
685 return final_candidate_elements;
716 std::set<const Elem *> point_neighbors;
719 for (
const auto & elem : point_neighbors)
777 _subdomain_ids = std::make_unique<std::set<subdomain_id_type>>(*subdomain_ids);
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
std::unique_ptr< std::set< subdomain_id_type > > _subdomain_ids
A default set of subdomain ids in which to search for points.
This is the EquationSystems class.
const EquationSystems & _eqn_systems
The equation systems handler, from which the data are gathered.
virtual std::unique_ptr< FunctionBase< Number > > clone() const override
void add_scaled(const TypeVector< T2 > &, const T &)
Add a scaled value to this vector without creating a temporary.
BuildType
enum defining how to build the tree.
const NumericVector< Number > & _vector
A reference to the vector that holds the data that is to be interpolated.
bool _out_of_mesh_mode
true if out-of-mesh mode is enabled.
const unsigned int invalid_uint
A number which is used quite often to represent an invalid or uninitialized value for an unsigned int...
std::unique_ptr< PointLocatorBase > sub_point_locator() const
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
static Point inverse_map(const unsigned int dim, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true, const bool extra_checks=true)
class FEComputeData hides arbitrary data to be passed to and from children of FEBase through the FEIn...
const DofMap & _dof_map
Need access to the DofMap of the other system.
void resize(const unsigned int n)
Resize the vector.
const Elem * find_element(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
Helper function to reduce code duplication.
const FEType & variable_type(const unsigned int c) const
This is the base class from which all geometric element types are derived.
const FunctionBase * _master
Const pointer to our master, initialized to nullptr.
bool _initialized
When init() was called so that everything is ready for calls to operator() (...), then this bool is t...
Gradient gradient(const Point &p, const Real time=0.)
std::vector< Gradient > dshape
Storage for the computed shape derivative values.
This class defines a vector in LIBMESH_DIM dimensional Real or Complex space.
~MeshFunction()
Destructor.
std::vector< std::vector< Real > > local_transform
Storage for local to global mapping at p.
The libMesh namespace provides an interface to certain functionality in the library.
void disable_out_of_mesh_mode()
Disables out-of-mesh mode.
Real get_close_to_point_tol() const
Get the close-to-point tolerance.
This is the MeshBase class.
void _gradient_on_elem(const Point &p, const Elem *element, std::vector< Gradient > &output)
Helper function for finding a gradient as evaluated from a specific element.
This class handles the numbering of degrees of freedom on a mesh.
DenseVector< Number > _out_of_mesh_value
Value to return outside the mesh if out-of-mesh mode is enabled.
void enable_out_of_mesh_mode(const DenseVector< Number > &value)
Enables out-of-mesh mode.
virtual void clear() override
Clears the function.
std::set< const Elem * > find_elements(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
std::map< const Elem *, Number > discontinuous_value(const Point &p, const Real time=0.)
void find_point_neighbors(const Point &p, std::set< const Elem *> &neighbor_set) const
This function finds all active elements (including this one) which are in the same manifold as this e...
static void compute_data(const unsigned int dim, const FEType &fe_t, const Elem *elem, FEComputeData &data)
Lets the appropriate child of FEBase compute the requested data for the input specified in data...
NumberVectorValue Gradient
void unset_point_locator_tolerance()
Turn off the user-specified PointLocator tolerance.
virtual void init() override
Override the FunctionBase::init() member function.
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
const PointLocatorBase & get_point_locator() const
This is the base class for point locators.
An object whose state is distributed along a set of processors.
void add_scaled(const TypeTensor< T2 > &, const T &)
Add a scaled tensor to this tensor without creating a temporary.
std::vector< Number > shape
Storage for the computed shape function values.
ParallelType type() const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual unsigned short dim() const =0
virtual Number operator()(const Point &p, const Real time=0.) override
MeshFunction(const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, std::vector< unsigned int > vars, const FunctionBase< Number > *master=nullptr)
Constructor for mesh based functions with vectors as return value.
const MeshBase & get_mesh() const
virtual unsigned int size() const override final
virtual bool infinite() const =0
This class provides function-like objects for data distributed over a mesh.
void set_point_locator_tolerance(Real tol)
We may want to specify a tolerance for the PointLocator to use, since in some cases the point we want...
const Elem * check_found_elem(const Elem *element, const Point &p) const
Helper function that is called by MeshFunction::find_element() and MeshFunction::find_elements() to e...
processor_id_type processor_id() const
void enable_derivative()
Enable the computation of shape gradients (dshape).
processor_id_type processor_id() const
std::map< const Elem *, Gradient > discontinuous_gradient(const Point &p, const Real time=0.)
A Point defines a location in LIBMESH_DIM dimensional Real space.
Tensor hessian(const Point &p, const Real time=0.)
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
void set_subdomain_ids(const std::set< subdomain_id_type > *subdomain_ids)
Choose a default list of subdomain ids to be searched for points.
std::unique_ptr< PointLocatorBase > _point_locator
A point locator is needed to locate the points in the mesh.
bool is_evaluable(const DofObjectSubclass &obj, unsigned int var_num=libMesh::invalid_uint) const
This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space.
const std::vector< unsigned int > _system_vars
The indices of the variables within the other system for which data are to be gathered.