libMesh::MeshFunction Class Reference

#include <mesh_function.h>

Inheritance diagram for libMesh::MeshFunction:

Public Member Functions

 MeshFunction (const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, const std::vector< unsigned int > &vars, const FunctionBase< Number > *master=nullptr)
 
 MeshFunction (const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, const unsigned int var, const FunctionBase< Number > *master=nullptr)
 
 MeshFunction (MeshFunction &&)=delete
 
 MeshFunction (const MeshFunction &)=delete
 
MeshFunctionoperator= (const MeshFunction &)=delete
 
MeshFunctionoperator= (MeshFunction &&)=delete
 
 ~MeshFunction ()
 
virtual void init () override
 
void init (const Trees::BuildType point_locator_build_type)
 
virtual void clear () override
 
virtual std::unique_ptr< FunctionBase< Number > > clone () const override
 
Number operator() (const Point &p, const Real time=0.) override
 
std::map< const Elem *, Numberdiscontinuous_value (const Point &p, const Real time=0.)
 
Gradient gradient (const Point &p, const Real time=0.)
 
std::map< const Elem *, Gradientdiscontinuous_gradient (const Point &p, const Real time=0.)
 
Tensor hessian (const Point &p, const Real time=0.)
 
void operator() (const Point &p, const Real time, DenseVector< Number > &output) override
 
void operator() (const Point &p, const Real time, DenseVector< Number > &output, const std::set< subdomain_id_type > *subdomain_ids)
 
void discontinuous_value (const Point &p, const Real time, std::map< const Elem *, DenseVector< Number >> &output)
 
void discontinuous_value (const Point &p, const Real time, std::map< const Elem *, DenseVector< Number >> &output, const std::set< subdomain_id_type > *subdomain_ids)
 
void gradient (const Point &p, const Real time, std::vector< Gradient > &output, const std::set< subdomain_id_type > *subdomain_ids=nullptr)
 
void discontinuous_gradient (const Point &p, const Real time, std::map< const Elem *, std::vector< Gradient >> &output)
 
void discontinuous_gradient (const Point &p, const Real time, std::map< const Elem *, std::vector< Gradient >> &output, const std::set< subdomain_id_type > *subdomain_ids)
 
void hessian (const Point &p, const Real time, std::vector< Tensor > &output, const std::set< subdomain_id_type > *subdomain_ids=nullptr)
 
const PointLocatorBaseget_point_locator (void) const
 
void enable_out_of_mesh_mode (const DenseVector< Number > &value)
 
void enable_out_of_mesh_mode (const Number &value)
 
void disable_out_of_mesh_mode (void)
 
void set_point_locator_tolerance (Real tol)
 
void unset_point_locator_tolerance ()
 
void operator() (const Point &p, DenseVector< Number > &output)
 
virtual Number component (unsigned int i, const Point &p, Real time=0.)
 
bool initialized () const
 
void set_is_time_dependent (bool is_time_dependent)
 
bool is_time_dependent () const
 
const Parallel::Communicatorcomm () const
 
processor_id_type n_processors () const
 
processor_id_type processor_id () const
 

Protected Member Functions

const Elemfind_element (const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
 
std::set< const Elem * > find_elements (const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
 

Protected Attributes

const EquationSystems_eqn_systems
 
const NumericVector< Number > & _vector
 
const DofMap_dof_map
 
const std::vector< unsigned int > _system_vars
 
PointLocatorBase_point_locator
 
bool _out_of_mesh_mode
 
DenseVector< Number_out_of_mesh_value
 
const FunctionBase_master
 
bool _initialized
 
bool _is_time_dependent
 
const Parallel::Communicator_communicator
 

Detailed Description

This class provides function-like objects for data distributed over a mesh.

Author
Daniel Dreyer
Date
2003

Definition at line 53 of file mesh_function.h.

Constructor & Destructor Documentation

◆ MeshFunction() [1/4]

libMesh::MeshFunction::MeshFunction ( const EquationSystems eqn_systems,
const NumericVector< Number > &  vec,
const DofMap dof_map,
const std::vector< unsigned int > &  vars,
const FunctionBase< Number > *  master = nullptr 
)

Constructor for mesh based functions with vectors as return value. Optionally takes a master function. If the MeshFunction is to be evaluated outside of the local partition of the mesh, then both the mesh in eqn_systems and the coefficient vector vec should be serialized.

Definition at line 43 of file mesh_function.C.

Referenced by clone().

47  :
48  FunctionBase<Number> (master),
49  ParallelObject (eqn_systems),
50  _eqn_systems (eqn_systems),
51  _vector (vec),
52  _dof_map (dof_map),
53  _system_vars (vars),
54  _point_locator (nullptr),
55  _out_of_mesh_mode (false),
57 {
58 }
ParallelObject(const Parallel::Communicator &comm_in)
const EquationSystems & _eqn_systems
const NumericVector< Number > & _vector
PointLocatorBase * _point_locator
const DofMap & _dof_map
DenseVector< Number > _out_of_mesh_value
const std::vector< unsigned int > _system_vars

◆ MeshFunction() [2/4]

libMesh::MeshFunction::MeshFunction ( const EquationSystems eqn_systems,
const NumericVector< Number > &  vec,
const DofMap dof_map,
const unsigned int  var,
const FunctionBase< Number > *  master = nullptr 
)

Constructor for mesh based functions with a number as return value. Optionally takes a master function. If the MeshFunction is to be evaluated outside of the local partition of the mesh, then both the mesh in eqn_systems and the coefficient vector vec should be serialized.

Definition at line 62 of file mesh_function.C.

66  :
67  FunctionBase<Number> (master),
68  ParallelObject (eqn_systems),
69  _eqn_systems (eqn_systems),
70  _vector (vec),
71  _dof_map (dof_map),
72  _system_vars (1,var),
73  _point_locator (nullptr),
74  _out_of_mesh_mode (false),
76 {
77  // std::vector<unsigned int> buf (1);
78  // buf[0] = var;
79  // _system_vars (buf);
80 }
ParallelObject(const Parallel::Communicator &comm_in)
const EquationSystems & _eqn_systems
const NumericVector< Number > & _vector
PointLocatorBase * _point_locator
const DofMap & _dof_map
DenseVector< Number > _out_of_mesh_value
const std::vector< unsigned int > _system_vars

◆ MeshFunction() [3/4]

libMesh::MeshFunction::MeshFunction ( MeshFunction &&  )
delete

This class is sometimes responsible for cleaning up the _point_locator, so it can't be default (shallow) copy constructed or move constructed.

◆ MeshFunction() [4/4]

libMesh::MeshFunction::MeshFunction ( const MeshFunction )
delete

◆ ~MeshFunction()

libMesh::MeshFunction::~MeshFunction ( )

Destructor.

Definition at line 88 of file mesh_function.C.

References libMesh::FunctionBase< Number >::_master, and _point_locator.

89 {
90  // only delete the point locator when we are the master
91  if (this->_master == nullptr)
92  delete this->_point_locator;
93 }
PointLocatorBase * _point_locator
const FunctionBase * _master

Member Function Documentation

◆ clear()

void libMesh::MeshFunction::clear ( )
overridevirtual

Clears the function.

Reimplemented from libMesh::FunctionBase< Number >.

Definition at line 148 of file mesh_function.C.

References libMesh::FunctionBase< Number >::_initialized, libMesh::FunctionBase< Number >::_master, and _point_locator.

149 {
150  // only delete the point locator when we are the master
151  if ((this->_point_locator != nullptr) && (this->_master == nullptr))
152  {
153  delete this->_point_locator;
154  this->_point_locator = nullptr;
155  }
156  this->_initialized = false;
157 }
PointLocatorBase * _point_locator
const FunctionBase * _master

◆ clone()

std::unique_ptr< FunctionBase< Number > > libMesh::MeshFunction::clone ( ) const
overridevirtual
Returns
A new copy of the function.

The new copy uses the original as a master function to enable simultaneous evaluations of the copies in different threads.

Note
This implies the copy should not be used after the original is destroyed.

Implements libMesh::FunctionBase< Number >.

Definition at line 161 of file mesh_function.C.

References _dof_map, _eqn_systems, _system_vars, _vector, libMesh::FunctionBase< Output >::init(), libMesh::FunctionBase< Number >::initialized(), and MeshFunction().

162 {
163  FunctionBase<Number> * mf_clone =
165 
166  if (this->initialized())
167  mf_clone->init();
168 
169  return std::unique_ptr<FunctionBase<Number>>(mf_clone);
170 }
const EquationSystems & _eqn_systems
const NumericVector< Number > & _vector
MeshFunction(const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, const std::vector< unsigned int > &vars, const FunctionBase< Number > *master=nullptr)
Definition: mesh_function.C:43
const DofMap & _dof_map
const std::vector< unsigned int > _system_vars

◆ comm()

const Parallel::Communicator& libMesh::ParallelObject::comm ( ) const
inlineinherited
Returns
A reference to the Parallel::Communicator object used by this mesh.

Definition at line 89 of file parallel_object.h.

References libMesh::ParallelObject::_communicator.

Referenced by libMesh::__libmesh_petsc_diff_solver_jacobian(), libMesh::__libmesh_petsc_diff_solver_monitor(), libMesh::__libmesh_petsc_diff_solver_residual(), libMesh::__libmesh_tao_equality_constraints(), libMesh::__libmesh_tao_equality_constraints_jacobian(), libMesh::__libmesh_tao_gradient(), libMesh::__libmesh_tao_hessian(), libMesh::__libmesh_tao_inequality_constraints(), libMesh::__libmesh_tao_inequality_constraints_jacobian(), libMesh::__libmesh_tao_objective(), libMesh::MeshRefinement::_coarsen_elements(), libMesh::ExactSolution::_compute_error(), libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult_add(), libMesh::EquationSystems::_read_impl(), libMesh::MeshRefinement::_refine_elements(), libMesh::MeshRefinement::_smooth_flags(), libMesh::PetscDMWrapper::add_dofs_helper(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::ImplicitSystem::add_matrix(), libMesh::System::add_vector(), libMesh::UnstructuredMesh::all_second_order(), libMesh::MeshTools::Modification::all_tri(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assemble_qoi(), libMesh::MeshCommunication::assign_global_indices(), libMesh::DofMap::attach_matrix(), libMesh::MeshTools::Generation::build_extrusion(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::EquationSystems::build_parallel_elemental_solution_vector(), libMesh::EquationSystems::build_parallel_solution_vector(), libMesh::PetscDMWrapper::build_section(), libMesh::PetscDMWrapper::build_sf(), libMesh::MeshBase::cache_elem_dims(), libMesh::System::calculate_norm(), libMesh::DofMap::check_dirichlet_bcid_consistency(), libMesh::PetscDMWrapper::check_section_n_dofs(), libMesh::Nemesis_IO_Helper::compute_num_global_elem_blocks(), libMesh::Nemesis_IO_Helper::compute_num_global_nodesets(), libMesh::Nemesis_IO_Helper::compute_num_global_sidesets(), libMesh::Problem_Interface::computeF(), libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), libMesh::ExodusII_IO::copy_elemental_solution(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::MeshTools::create_bounding_box(), libMesh::MeshTools::create_nodal_bounding_box(), libMesh::MeshRefinement::create_parent_error_vector(), libMesh::MeshTools::create_processor_bounding_box(), libMesh::MeshTools::create_subdomain_bounding_box(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), DMlibMeshFunction(), DMlibMeshJacobian(), DMlibMeshSetSystem_libMesh(), DMVariableBounds_libMesh(), libMesh::MeshRefinement::eliminate_unrefined_patches(), libMesh::EpetraVector< T >::EpetraVector(), libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::MeshRefinement::flag_elements_by_elem_fraction(), libMesh::MeshRefinement::flag_elements_by_error_fraction(), libMesh::MeshRefinement::flag_elements_by_nelem_target(), libMesh::CondensedEigenSystem::get_eigenpair(), libMesh::DofMap::get_info(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::LocationMap< T >::init(), libMesh::TimeSolver::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::PetscDMWrapper::init_and_attach_petscdm(), libMesh::EigenSystem::init_data(), libMesh::EigenSystem::init_matrices(), libMesh::OptimizationSystem::initialize_equality_constraints_storage(), libMesh::OptimizationSystem::initialize_inequality_constraints_storage(), libMesh::MeshTools::libmesh_assert_consistent_distributed(), libMesh::MeshTools::libmesh_assert_consistent_distributed_nodes(), libMesh::MeshTools::libmesh_assert_contiguous_dof_ids(), libMesh::MeshTools::libmesh_assert_parallel_consistent_new_node_procids(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_topology_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_boundary_ids(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_flags(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_object_ids(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_p_levels(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::MeshTools::libmesh_assert_valid_unique_ids(), libMesh::libmesh_petsc_snes_fd_residual(), libMesh::libmesh_petsc_snes_jacobian(), libMesh::libmesh_petsc_snes_mffd_residual(), libMesh::libmesh_petsc_snes_postcheck(), libMesh::libmesh_petsc_snes_residual(), libMesh::libmesh_petsc_snes_residual_helper(), libMesh::MeshRefinement::limit_level_mismatch_at_edge(), libMesh::MeshRefinement::limit_level_mismatch_at_node(), libMesh::MeshRefinement::limit_overrefined_boundary(), libMesh::MeshRefinement::limit_underrefined_boundary(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshCommunication::make_elems_parallel_consistent(), libMesh::MeshRefinement::make_flags_parallel_consistent(), libMesh::MeshCommunication::make_new_node_proc_ids_parallel_consistent(), libMesh::MeshCommunication::make_new_nodes_parallel_consistent(), libMesh::MeshCommunication::make_node_ids_parallel_consistent(), libMesh::MeshCommunication::make_node_proc_ids_parallel_consistent(), libMesh::MeshCommunication::make_node_unique_ids_parallel_consistent(), libMesh::MeshCommunication::make_nodes_parallel_consistent(), libMesh::MeshCommunication::make_p_levels_parallel_consistent(), libMesh::MeshRefinement::make_refinement_compatible(), libMesh::FEMSystem::mesh_position_set(), libMesh::DistributedMesh::n_active_elem(), libMesh::MeshTools::n_active_levels(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::CondensedEigenSystem::n_global_non_condensed_dofs(), libMesh::MeshTools::n_levels(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::MeshTools::n_p_levels(), libMesh::BoundaryInfo::n_shellface_conds(), libMesh::DistributedMesh::parallel_max_elem_id(), libMesh::DistributedMesh::parallel_max_node_id(), libMesh::ReplicatedMesh::parallel_max_unique_id(), libMesh::DistributedMesh::parallel_max_unique_id(), libMesh::DistributedMesh::parallel_n_elem(), libMesh::DistributedMesh::parallel_n_nodes(), libMesh::SparsityPattern::Build::parallel_sync(), libMesh::MeshTools::paranoid_n_levels(), libMesh::petsc_auto_fieldsplit(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::MeshBase::prepare_for_use(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::System::read_header(), libMesh::System::read_legacy_data(), libMesh::System::read_SCALAR_dofs(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::System::read_serialized_vector(), libMesh::MeshBase::recalculate_n_partitions(), libMesh::MeshRefinement::refine_and_coarsen_elements(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::CheckpointIO::select_split_config(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::PetscDiffSolver::setup_petsc_data(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::split_mesh(), libMesh::MeshBase::subdomain_ids(), libMesh::BoundaryInfo::sync(), libMesh::MeshRefinement::test_level_one(), libMesh::MeshRefinement::test_unflagged(), libMesh::MeshTools::total_weight(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), and libMesh::XdrIO::write_serialized_nodesets().

90  { return _communicator; }
const Parallel::Communicator & _communicator

◆ component()

Number libMesh::FunctionBase< Number >::component ( unsigned int  i,
const Point p,
Real  time = 0. 
)
inlinevirtualinherited
Returns
The vector component i at coordinate p and time time.
Note
Subclasses aren't required to override this, since the default implementation is based on the full vector evaluation, which is often correct.
Subclasses are recommended to override this, since the default implementation is based on a vector evaluation, which is usually unnecessarily inefficient.

Definition at line 228 of file function_base.h.

231 {
232  DenseVector<Output> outvec(i+1);
233  (*this)(p, time, outvec);
234  return outvec(i);
235 }

◆ disable_out_of_mesh_mode()

void libMesh::MeshFunction::disable_out_of_mesh_mode ( void  )

Disables out-of-mesh mode. This is also the default.

Definition at line 788 of file mesh_function.C.

References _out_of_mesh_mode, _point_locator, libMesh::PointLocatorBase::disable_out_of_mesh_mode(), and libMesh::FunctionBase< Number >::initialized().

789 {
790  libmesh_assert (this->initialized());
792  _out_of_mesh_mode = false;
793 }
virtual void disable_out_of_mesh_mode()=0
PointLocatorBase * _point_locator

◆ discontinuous_gradient() [1/3]

std::map< const Elem *, Gradient > libMesh::MeshFunction::discontinuous_gradient ( const Point p,
const Real  time = 0. 
)
Returns
A map of first derivatives (gradients) of variable 0 at point p and for time. map is from element to Gradient and accounts for double defined values on faces if the gradient is discontinuous

Definition at line 209 of file mesh_function.C.

References libMesh::FunctionBase< Number >::initialized().

Referenced by discontinuous_gradient().

211 {
212  libmesh_assert (this->initialized());
213 
214  std::map<const Elem *, std::vector<Gradient>> buffer;
215  this->discontinuous_gradient (p, time, buffer);
216  std::map<const Elem *, Gradient> return_value;
217  for (const auto & pr : buffer)
218  return_value[pr.first] = pr.second[0];
219  // NOTE: If no suitable element is found, then the map return_value is empty. This
220  // puts burden on the user of this function but I don't really see a better way.
221  return return_value;
222 }
std::map< const Elem *, Gradient > discontinuous_gradient(const Point &p, const Real time=0.)

◆ discontinuous_gradient() [2/3]

void libMesh::MeshFunction::discontinuous_gradient ( const Point p,
const Real  time,
std::map< const Elem *, std::vector< Gradient >> &  output 
)

Similar to gradient, but with the difference that multiple values on faces are explicitly permitted. This is useful for evaluating gradients on faces where the values to the left and right are different.

Definition at line 508 of file mesh_function.C.

References discontinuous_gradient().

511 {
512  this->discontinuous_gradient (p, time, output, nullptr);
513 }
std::map< const Elem *, Gradient > discontinuous_gradient(const Point &p, const Real time=0.)

◆ discontinuous_gradient() [3/3]

void libMesh::MeshFunction::discontinuous_gradient ( const Point p,
const Real  time,
std::map< const Elem *, std::vector< Gradient >> &  output,
const std::set< subdomain_id_type > *  subdomain_ids 
)

Similar to gradient, but with the difference that multiple values on faces are explicitly permitted. This is useful for evaluating gradients on faces where the values to the left and right are different.

Definition at line 517 of file mesh_function.C.

References _dof_map, _out_of_mesh_mode, _out_of_mesh_value, _system_vars, libMesh::TypeVector< T >::add_scaled(), libMesh::FEGenericBase< OutputType >::build(), libMesh::DofMap::dof_indices(), find_elements(), libMesh::FunctionBase< Number >::initialized(), libMesh::invalid_uint, libMesh::FEInterface::inverse_map(), libMesh::DenseVector< T >::size(), and libMesh::DofMap::variable_type().

521 {
522  libmesh_assert (this->initialized());
523 
524  // clear the output map
525  output.clear();
526 
527  // get the candidate elements
528  std::set<const Elem *> candidate_element = this->find_elements(p,subdomain_ids);
529 
530  // loop through all candidates, if the set is empty this function will return an
531  // empty map
532  for (const auto & element : candidate_element)
533  {
534  const unsigned int dim = element->dim();
535 
536  // define a temporary vector to store all values
537  std::vector<Gradient> temp_output (cast_int<unsigned int>(this->_system_vars.size()));
538 
539  /*
540  * Get local coordinates to feed these into compute_data().
541  * Note that the fe_type can safely be used from the 0-variable,
542  * since the inverse mapping is the same for all FEFamilies
543  */
544  const Point mapped_point (FEInterface::inverse_map (dim,
545  this->_dof_map.variable_type(0),
546  element,
547  p));
548 
549 
550  // loop over all vars
551  std::vector<Point> point_list (1, mapped_point);
552  for (unsigned int index=0,
553  sz = cast_int<unsigned int>(this->_system_vars.size());
554  index != sz; ++index)
555  {
556  /*
557  * the data for this variable
558  */
559  const unsigned int var = _system_vars[index];
560 
561  if (var == libMesh::invalid_uint)
562  {
563  libmesh_assert (_out_of_mesh_mode &&
564  index < _out_of_mesh_value.size());
565  temp_output[index] = Gradient(_out_of_mesh_value(index));
566  continue;
567  }
568 
569  const FEType & fe_type = this->_dof_map.variable_type(var);
570 
571  std::unique_ptr<FEBase> point_fe (FEBase::build(dim, fe_type));
572  const std::vector<std::vector<RealGradient>> & dphi = point_fe->get_dphi();
573  point_fe->reinit(element, &point_list);
574 
575  // where the solution values for the var-th variable are stored
576  std::vector<dof_id_type> dof_indices;
577  this->_dof_map.dof_indices (element, dof_indices, var);
578 
579  Gradient grad(0.);
580 
581  for (std::size_t i = 0; i < dof_indices.size(); ++i)
582  grad.add_scaled(dphi[i][0], this->_vector(dof_indices[i]));
583 
584  temp_output[index] = grad;
585 
586  // next variable
587  }
588 
589  // Insert temp_output into output
590  output[element] = temp_output;
591  }
592 }
virtual unsigned int size() const override
Definition: dense_vector.h:92
const unsigned int invalid_uint
Definition: libmesh.h:245
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Definition: dof_map.C:1930
const DofMap & _dof_map
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1792
DenseVector< Number > _out_of_mesh_value
std::set< const Elem * > find_elements(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
NumberVectorValue Gradient
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
static Point inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: fe_interface.C:590
const std::vector< unsigned int > _system_vars

◆ discontinuous_value() [1/3]

std::map< const Elem *, Number > libMesh::MeshFunction::discontinuous_value ( const Point p,
const Real  time = 0. 
)
Returns
A map of values of variable 0 at point p and for time.

The std::map is from element to Number and accounts for doubly-defined values on faces if discontinuous variables are used.

Definition at line 184 of file mesh_function.C.

References libMesh::FunctionBase< Number >::initialized().

Referenced by discontinuous_value().

186 {
187  libmesh_assert (this->initialized());
188 
189  std::map<const Elem *, DenseVector<Number>> buffer;
190  this->discontinuous_value (p, time, buffer);
191  std::map<const Elem *, Number> return_value;
192  for (const auto & pr : buffer)
193  return_value[pr.first] = pr.second(0);
194  // NOTE: If no suitable element is found, then the map return_value is empty. This
195  // puts burden on the user of this function but I don't really see a better way.
196  return return_value;
197 }
std::map< const Elem *, Number > discontinuous_value(const Point &p, const Real time=0.)

◆ discontinuous_value() [2/3]

void libMesh::MeshFunction::discontinuous_value ( const Point p,
const Real  time,
std::map< const Elem *, DenseVector< Number >> &  output 
)

Similar to operator() with the same parameter list, but with the difference that multiple values on faces are explicitly permitted. This is useful for discontinuous shape functions that are evaluated on faces.

Definition at line 334 of file mesh_function.C.

References discontinuous_value().

337 {
338  this->discontinuous_value (p, time, output, nullptr);
339 }
std::map< const Elem *, Number > discontinuous_value(const Point &p, const Real time=0.)

◆ discontinuous_value() [3/3]

void libMesh::MeshFunction::discontinuous_value ( const Point p,
const Real  time,
std::map< const Elem *, DenseVector< Number >> &  output,
const std::set< subdomain_id_type > *  subdomain_ids 
)

Similar to operator() with the same parameter list, but with the difference that multiple values on faces are explicitly permitted. This is useful for discontinuous shape functions that are evaluated on faces.

Build an FEComputeData that contains both input and output data for the specific compute_data method.

Definition at line 343 of file mesh_function.C.

References _dof_map, _eqn_systems, _out_of_mesh_mode, _out_of_mesh_value, _system_vars, _vector, libMesh::FEInterface::compute_data(), data, libMesh::DofMap::dof_indices(), find_elements(), libMesh::FunctionBase< Number >::initialized(), libMesh::invalid_uint, libMesh::FEInterface::inverse_map(), libMesh::DenseVector< T >::size(), value, and libMesh::DofMap::variable_type().

347 {
348  libmesh_assert (this->initialized());
349 
350  // clear the output map
351  output.clear();
352 
353  // get the candidate elements
354  std::set<const Elem *> candidate_element = this->find_elements(p,subdomain_ids);
355 
356  // loop through all candidates, if the set is empty this function will return an
357  // empty map
358  for (const auto & element : candidate_element)
359  {
360  const unsigned int dim = element->dim();
361 
362  // define a temporary vector to store all values
363  DenseVector<Number> temp_output (cast_int<unsigned int>(this->_system_vars.size()));
364 
365  /*
366  * Get local coordinates to feed these into compute_data().
367  * Note that the fe_type can safely be used from the 0-variable,
368  * since the inverse mapping is the same for all FEFamilies
369  */
370  const Point mapped_point (FEInterface::inverse_map (dim,
371  this->_dof_map.variable_type(0),
372  element,
373  p));
374 
375 
376  // loop over all vars
377  for (unsigned int index=0,
378  sz = cast_int<unsigned int>(this->_system_vars.size());
379  index != sz; ++index)
380  {
381  /*
382  * the data for this variable
383  */
384  const unsigned int var = _system_vars[index];
385 
386  if (var == libMesh::invalid_uint)
387  {
388  libmesh_assert (_out_of_mesh_mode &&
389  index < _out_of_mesh_value.size());
390  temp_output(index) = _out_of_mesh_value(index);
391  continue;
392  }
393 
394  const FEType & fe_type = this->_dof_map.variable_type(var);
395 
400  {
401  FEComputeData data (this->_eqn_systems, mapped_point);
402 
403  FEInterface::compute_data (dim, fe_type, element, data);
404 
405  // where the solution values for the var-th variable are stored
406  std::vector<dof_id_type> dof_indices;
407  this->_dof_map.dof_indices (element, dof_indices, var);
408 
409  // interpolate the solution
410  {
411  Number value = 0.;
412 
413  for (std::size_t i=0; i<dof_indices.size(); i++)
414  value += this->_vector(dof_indices[i]) * data.shape[i];
415 
416  temp_output(index) = value;
417  }
418 
419  }
420 
421  // next variable
422  }
423 
424  // Insert temp_output into output
425  output[element] = temp_output;
426  }
427 }
virtual unsigned int size() const override
Definition: dense_vector.h:92
const EquationSystems & _eqn_systems
const NumericVector< Number > & _vector
const unsigned int invalid_uint
Definition: libmesh.h:245
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Definition: dof_map.C:1930
const DofMap & _dof_map
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1792
DenseVector< Number > _out_of_mesh_value
std::set< const Elem * > find_elements(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
static void compute_data(const unsigned int dim, const FEType &fe_t, const Elem *elem, FEComputeData &data)
Definition: fe_interface.C:837
static Point inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: fe_interface.C:590
static const bool value
Definition: xdr_io.C:109
IterBase * data
const std::vector< unsigned int > _system_vars

◆ enable_out_of_mesh_mode() [1/2]

void libMesh::MeshFunction::enable_out_of_mesh_mode ( const DenseVector< Number > &  value)

Enables out-of-mesh mode. In this mode, if asked for a point that is not contained in any element, the MeshFunction will return the given value instead of crashing. This mode is off per default. If you use a master mesh function and you want to enable this mode, you will have to enable it for the master mesh function as well and for all mesh functions that have the same master mesh function. You may, however, specify different values.

Definition at line 773 of file mesh_function.C.

References _out_of_mesh_mode, _out_of_mesh_value, _point_locator, libMesh::PointLocatorBase::enable_out_of_mesh_mode(), libMesh::FunctionBase< Number >::initialized(), and value.

Referenced by enable_out_of_mesh_mode().

774 {
775  libmesh_assert (this->initialized());
777  _out_of_mesh_mode = true;
779 }
PointLocatorBase * _point_locator
DenseVector< Number > _out_of_mesh_value
static const bool value
Definition: xdr_io.C:109
virtual void enable_out_of_mesh_mode()=0

◆ enable_out_of_mesh_mode() [2/2]

void libMesh::MeshFunction::enable_out_of_mesh_mode ( const Number value)

Enables out-of-mesh mode. In this mode, if asked for a point that is not contained in any element, the MeshFunction will return the given value instead of crashing. This mode is off per default. If you use a master mesh function and you want to enable this mode, you will have to enable it for the master mesh function as well and for all mesh functions that have the same master mesh function. You may, however, specify different values.

Definition at line 781 of file mesh_function.C.

References enable_out_of_mesh_mode(), and value.

782 {
783  DenseVector<Number> v(1);
784  v(0) = value;
785  this->enable_out_of_mesh_mode(v);
786 }
void enable_out_of_mesh_mode(const DenseVector< Number > &value)
static const bool value
Definition: xdr_io.C:109

◆ find_element()

const Elem * libMesh::MeshFunction::find_element ( const Point p,
const std::set< subdomain_id_type > *  subdomain_ids = nullptr 
) const
protected

Helper function to reduce code duplication

Definition at line 674 of file mesh_function.C.

References libMesh::FunctionBase< Number >::_master, _out_of_mesh_mode, _point_locator, _vector, libMesh::Elem::find_point_neighbors(), libMesh::ParallelObject::processor_id(), libMesh::DofObject::processor_id(), libMesh::SERIAL, and libMesh::NumericVector< T >::type().

Referenced by gradient(), hessian(), and operator()().

676 {
677  /* Ensure that in the case of a master mesh function, the
678  out-of-mesh mode is enabled either for both or for none. This is
679  important because the out-of-mesh mode is also communicated to
680  the point locator. Since this is time consuming, enable it only
681  in debug mode. */
682 #ifdef DEBUG
683  if (this->_master != nullptr)
684  {
685  const MeshFunction * master =
686  cast_ptr<const MeshFunction *>(this->_master);
687  if (_out_of_mesh_mode!=master->_out_of_mesh_mode)
688  libmesh_error_msg("ERROR: If you use out-of-mesh-mode in connection with master mesh " \
689  << "functions, you must enable out-of-mesh mode for both the master and the slave mesh function.");
690  }
691 #endif
692 
693  // locate the point in the other mesh
694  const Elem * element = this->_point_locator->operator()(p,subdomain_ids);
695 
696  // If we have an element, but it's not a local element, then we
697  // either need to have a serialized vector or we need to find a
698  // local element sharing the same point.
699  if (element &&
700  (element->processor_id() != this->processor_id()) &&
701  _vector.type() != SERIAL)
702  {
703  // look for a local element containing the point
704  std::set<const Elem *> point_neighbors;
705  element->find_point_neighbors(p, point_neighbors);
706  element = nullptr;
707  for (const auto & elem : point_neighbors)
708  if (elem->processor_id() == this->processor_id())
709  {
710  element = elem;
711  break;
712  }
713  }
714 
715  return element;
716 }
const NumericVector< Number > & _vector
PointLocatorBase * _point_locator
MeshFunction(const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, const std::vector< unsigned int > &vars, const FunctionBase< Number > *master=nullptr)
Definition: mesh_function.C:43
const FunctionBase * _master
ParallelType type() const
processor_id_type processor_id() const

◆ find_elements()

std::set< const Elem * > libMesh::MeshFunction::find_elements ( const Point p,
const std::set< subdomain_id_type > *  subdomain_ids = nullptr 
) const
protected
Returns
All elements that are close to a point p.

This is similar to find_element() but covers cases where p is on the boundary.

Definition at line 718 of file mesh_function.C.

References libMesh::FunctionBase< Number >::_master, _out_of_mesh_mode, _point_locator, _vector, libMesh::ParallelObject::processor_id(), libMesh::SERIAL, and libMesh::NumericVector< T >::type().

Referenced by discontinuous_gradient(), and discontinuous_value().

720 {
721  /* Ensure that in the case of a master mesh function, the
722  out-of-mesh mode is enabled either for both or for none. This is
723  important because the out-of-mesh mode is also communicated to
724  the point locator. Since this is time consuming, enable it only
725  in debug mode. */
726 #ifdef DEBUG
727  if (this->_master != nullptr)
728  {
729  const MeshFunction * master =
730  cast_ptr<const MeshFunction *>(this->_master);
731  if (_out_of_mesh_mode!=master->_out_of_mesh_mode)
732  libmesh_error_msg("ERROR: If you use out-of-mesh-mode in connection with master mesh " \
733  << "functions, you must enable out-of-mesh mode for both the master and the slave mesh function.");
734  }
735 #endif
736 
737  // locate the point in the other mesh
738  std::set<const Elem *> candidate_elements;
739  std::set<const Elem *> final_candidate_elements;
740  this->_point_locator->operator()(p,candidate_elements,subdomain_ids);
741  for (const auto & element : candidate_elements)
742  {
743  // If we have an element, but it's not a local element, then we
744  // either need to have a serialized vector or we need to find a
745  // local element sharing the same point.
746  if (element &&
747  (element->processor_id() != this->processor_id()) &&
748  _vector.type() != SERIAL)
749  {
750  // look for a local element containing the point
751  std::set<const Elem *> point_neighbors;
752  element->find_point_neighbors(p, point_neighbors);
753  for (const auto & elem : point_neighbors)
754  if (elem->processor_id() == this->processor_id())
755  {
756  final_candidate_elements.insert(elem);
757  break;
758  }
759  }
760  else
761  final_candidate_elements.insert(element);
762  }
763 
764  return final_candidate_elements;
765 }
const NumericVector< Number > & _vector
PointLocatorBase * _point_locator
MeshFunction(const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, const std::vector< unsigned int > &vars, const FunctionBase< Number > *master=nullptr)
Definition: mesh_function.C:43
const FunctionBase * _master
ParallelType type() const
processor_id_type processor_id() const

◆ get_point_locator()

const PointLocatorBase & libMesh::MeshFunction::get_point_locator ( void  ) const
Returns
The current PointLocator object, for use elsewhere.
Note
The MeshFunction object must be initialized before this is called.

Definition at line 767 of file mesh_function.C.

References _point_locator, and libMesh::FunctionBase< Number >::initialized().

768 {
769  libmesh_assert (this->initialized());
770  return *_point_locator;
771 }
PointLocatorBase * _point_locator

◆ gradient() [1/2]

Gradient libMesh::MeshFunction::gradient ( const Point p,
const Real  time = 0. 
)
Returns
The first derivatives of variable 0 at point p and for time, which defaults to zero.

Definition at line 199 of file mesh_function.C.

References libMesh::FunctionBase< Number >::initialized().

Referenced by libMesh::ExactErrorEstimator::find_squared_element_error().

201 {
202  libmesh_assert (this->initialized());
203 
204  std::vector<Gradient> buf (1);
205  this->gradient(p, time, buf);
206  return buf[0];
207 }
Gradient gradient(const Point &p, const Real time=0.)

◆ gradient() [2/2]

void libMesh::MeshFunction::gradient ( const Point p,
const Real  time,
std::vector< Gradient > &  output,
const std::set< subdomain_id_type > *  subdomain_ids = nullptr 
)

Computes gradients at coordinate p and for time time, which defaults to zero, optionally restricting the point to the passed subdomain_ids. This is useful in cases where there are multiple dimensioned elements, for example.

Definition at line 431 of file mesh_function.C.

References _dof_map, _out_of_mesh_mode, _out_of_mesh_value, _system_vars, libMesh::TypeVector< T >::add_scaled(), libMesh::FEGenericBase< OutputType >::build(), libMesh::Elem::dim(), libMesh::DofMap::dof_indices(), find_element(), libMesh::FunctionBase< Number >::initialized(), libMesh::invalid_uint, libMesh::FEInterface::inverse_map(), libMesh::DenseVector< T >::size(), and libMesh::DofMap::variable_type().

435 {
436  libmesh_assert (this->initialized());
437 
438  const Elem * element = this->find_element(p,subdomain_ids);
439 
440  if (!element)
441  {
442  output.resize(0);
443  }
444  else
445  {
446  // resize the output vector to the number of output values
447  // that the user told us
448  output.resize (this->_system_vars.size());
449 
450 
451  {
452  const unsigned int dim = element->dim();
453 
454 
455  /*
456  * Get local coordinates to feed these into compute_data().
457  * Note that the fe_type can safely be used from the 0-variable,
458  * since the inverse mapping is the same for all FEFamilies
459  */
460  const Point mapped_point (FEInterface::inverse_map (dim,
461  this->_dof_map.variable_type(0),
462  element,
463  p));
464 
465  std::vector<Point> point_list (1, mapped_point);
466 
467  // loop over all vars
468  for (unsigned int index=0,
469  sz = cast_int<unsigned int>(this->_system_vars.size());
470  index != sz; ++index)
471  {
472  /*
473  * the data for this variable
474  */
475  const unsigned int var = _system_vars[index];
476 
477  if (var == libMesh::invalid_uint)
478  {
479  libmesh_assert (_out_of_mesh_mode &&
480  index < _out_of_mesh_value.size());
481  output[index] = Gradient(_out_of_mesh_value(index));
482  continue;
483  }
484 
485  const FEType & fe_type = this->_dof_map.variable_type(var);
486 
487  std::unique_ptr<FEBase> point_fe (FEBase::build(dim, fe_type));
488  const std::vector<std::vector<RealGradient>> & dphi = point_fe->get_dphi();
489  point_fe->reinit(element, &point_list);
490 
491  // where the solution values for the var-th variable are stored
492  std::vector<dof_id_type> dof_indices;
493  this->_dof_map.dof_indices (element, dof_indices, var);
494 
495  // interpolate the solution
496  Gradient grad(0.);
497 
498  for (std::size_t i=0; i<dof_indices.size(); i++)
499  grad.add_scaled(dphi[i][0], this->_vector(dof_indices[i]));
500 
501  output[index] = grad;
502  }
503  }
504  }
505 }
virtual unsigned int size() const override
Definition: dense_vector.h:92
const unsigned int invalid_uint
Definition: libmesh.h:245
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Definition: dof_map.C:1930
const DofMap & _dof_map
const Elem * find_element(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1792
DenseVector< Number > _out_of_mesh_value
NumberVectorValue Gradient
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
static Point inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: fe_interface.C:590
const std::vector< unsigned int > _system_vars

◆ hessian() [1/2]

Tensor libMesh::MeshFunction::hessian ( const Point p,
const Real  time = 0. 
)
Returns
The second derivatives of variable 0 at point p and for time, which defaults to zero.

Definition at line 225 of file mesh_function.C.

References libMesh::FunctionBase< Number >::initialized().

Referenced by libMesh::ExactErrorEstimator::find_squared_element_error().

227 {
228  libmesh_assert (this->initialized());
229 
230  std::vector<Tensor> buf (1);
231  this->hessian(p, time, buf);
232  return buf[0];
233 }
Tensor hessian(const Point &p, const Real time=0.)

◆ hessian() [2/2]

void libMesh::MeshFunction::hessian ( const Point p,
const Real  time,
std::vector< Tensor > &  output,
const std::set< subdomain_id_type > *  subdomain_ids = nullptr 
)

Computes gradients at coordinate p and for time time, which defaults to zero, optionally restricting the point to the passed subdomain_ids. This is useful in cases where there are multiple dimensioned elements, for example.

Definition at line 597 of file mesh_function.C.

References _dof_map, _out_of_mesh_mode, _out_of_mesh_value, _system_vars, libMesh::TypeTensor< T >::add_scaled(), libMesh::FEGenericBase< OutputType >::build(), libMesh::Elem::dim(), libMesh::DofMap::dof_indices(), find_element(), libMesh::FunctionBase< Number >::initialized(), libMesh::invalid_uint, libMesh::FEInterface::inverse_map(), libMesh::DenseVector< T >::size(), and libMesh::DofMap::variable_type().

601 {
602  libmesh_assert (this->initialized());
603 
604  const Elem * element = this->find_element(p,subdomain_ids);
605 
606  if (!element)
607  {
608  output.resize(0);
609  }
610  else
611  {
612  // resize the output vector to the number of output values
613  // that the user told us
614  output.resize (this->_system_vars.size());
615 
616 
617  {
618  const unsigned int dim = element->dim();
619 
620 
621  /*
622  * Get local coordinates to feed these into compute_data().
623  * Note that the fe_type can safely be used from the 0-variable,
624  * since the inverse mapping is the same for all FEFamilies
625  */
626  const Point mapped_point (FEInterface::inverse_map (dim,
627  this->_dof_map.variable_type(0),
628  element,
629  p));
630 
631  std::vector<Point> point_list (1, mapped_point);
632 
633  // loop over all vars
634  for (unsigned int index=0,
635  sz = cast_int<unsigned int>(this->_system_vars.size());
636  index != sz; ++index)
637  {
638  /*
639  * the data for this variable
640  */
641  const unsigned int var = _system_vars[index];
642 
643  if (var == libMesh::invalid_uint)
644  {
645  libmesh_assert (_out_of_mesh_mode &&
646  index < _out_of_mesh_value.size());
647  output[index] = Tensor(_out_of_mesh_value(index));
648  continue;
649  }
650  const FEType & fe_type = this->_dof_map.variable_type(var);
651 
652  std::unique_ptr<FEBase> point_fe (FEBase::build(dim, fe_type));
653  const std::vector<std::vector<RealTensor>> & d2phi =
654  point_fe->get_d2phi();
655  point_fe->reinit(element, &point_list);
656 
657  // where the solution values for the var-th variable are stored
658  std::vector<dof_id_type> dof_indices;
659  this->_dof_map.dof_indices (element, dof_indices, var);
660 
661  // interpolate the solution
662  Tensor hess;
663 
664  for (std::size_t i=0; i<dof_indices.size(); i++)
665  hess.add_scaled(d2phi[i][0], this->_vector(dof_indices[i]));
666 
667  output[index] = hess;
668  }
669  }
670  }
671 }
virtual unsigned int size() const override
Definition: dense_vector.h:92
const unsigned int invalid_uint
Definition: libmesh.h:245
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Definition: dof_map.C:1930
void add_scaled(const TypeTensor< T2 > &, const T)
Definition: type_tensor.h:808
const DofMap & _dof_map
const Elem * find_element(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1792
DenseVector< Number > _out_of_mesh_value
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
static Point inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: fe_interface.C:590
NumberTensorValue Tensor
const std::vector< unsigned int > _system_vars

◆ init() [1/2]

virtual void libMesh::MeshFunction::init ( )
inlineoverridevirtual

Override the FunctionBase::init() member function by calling our own and specifying the Trees::NODES method. specifies the method to use when building a PointLocator

Reimplemented from libMesh::FunctionBase< Number >.

Definition at line 110 of file mesh_function.h.

References init(), and libMesh::Trees::NODES.

Referenced by init().

110 { this->init(Trees::NODES); }
virtual void init() override

◆ init() [2/2]

void libMesh::MeshFunction::init ( const Trees::BuildType  point_locator_build_type)

The actual initialization process. Takes an optional argument which specifies the method to use when building a PointLocator

Definition at line 98 of file mesh_function.C.

References _eqn_systems, libMesh::FunctionBase< Number >::_initialized, libMesh::FunctionBase< Number >::_master, _point_locator, _system_vars, libMesh::EquationSystems::get_mesh(), mesh, and libMesh::MeshBase::sub_point_locator().

99 {
100  // are indices of the desired variable(s) provided?
101  libmesh_assert_greater (this->_system_vars.size(), 0);
102 
103  // Don't do twice...
104  if (this->_initialized)
105  {
106  libmesh_assert(this->_point_locator);
107  return;
108  }
109 
110  /*
111  * set up the PointLocator: either someone else
112  * is the master (go and get the address of his
113  * point locator) or this object is the master
114  * (build the point locator on our own).
115  */
116  if (this->_master != nullptr)
117  {
118  // we aren't the master
119  const MeshFunction * master =
120  cast_ptr<const MeshFunction *>(this->_master);
121 
122  if (master->_point_locator == nullptr)
123  libmesh_error_msg("ERROR: When the master-servant concept is used, the master has to be initialized first!");
124 
125  else
126  {
127  this->_point_locator = master->_point_locator;
128  }
129  }
130  else
131  {
132  // we are the master: build the point locator
133 
134  // constant reference to the other mesh
135  const MeshBase & mesh = this->_eqn_systems.get_mesh();
136 
137  // Get PointLocator object from the Mesh. We are responsible for
138  // deleting this only if we are the master.
139  this->_point_locator = mesh.sub_point_locator().release();
140  }
141 
142  // ready for use
143  this->_initialized = true;
144 }
const EquationSystems & _eqn_systems
PointLocatorBase * _point_locator
MeshFunction(const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, const std::vector< unsigned int > &vars, const FunctionBase< Number > *master=nullptr)
Definition: mesh_function.C:43
MeshBase & mesh
const FunctionBase * _master
const MeshBase & get_mesh() const
const std::vector< unsigned int > _system_vars

◆ initialized()

bool libMesh::FunctionBase< Number >::initialized ( ) const
inlineinherited
Returns
true when this object is properly initialized and ready for use, false otherwise.

Definition at line 206 of file function_base.h.

Referenced by clone(), disable_out_of_mesh_mode(), discontinuous_gradient(), discontinuous_value(), enable_out_of_mesh_mode(), get_point_locator(), gradient(), hessian(), and operator()().

207 {
208  return (this->_initialized);
209 }

◆ is_time_dependent()

bool libMesh::FunctionBase< Number >::is_time_dependent ( ) const
inlineinherited
Returns
true when the function this object represents is actually time-dependent, false otherwise.

Definition at line 220 of file function_base.h.

221 {
222  return (this->_is_time_dependent);
223 }

◆ n_processors()

processor_id_type libMesh::ParallelObject::n_processors ( ) const
inlineinherited
Returns
The number of processors in the group.

Definition at line 95 of file parallel_object.h.

References libMesh::ParallelObject::_communicator, and libMesh::Parallel::Communicator::size().

Referenced by libMesh::BoundaryInfo::_find_id_maps(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::DistributedMesh::add_elem(), libMesh::DistributedMesh::add_node(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assembly(), libMesh::AztecLinearSolver< T >::AztecLinearSolver(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::EquationSystems::build_parallel_elemental_solution_vector(), libMesh::DistributedMesh::clear(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::UnstructuredMesh::create_pid_mesh(), libMesh::MeshTools::create_processor_bounding_box(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::EnsightIO::EnsightIO(), libMesh::MeshBase::get_info(), libMesh::SystemSubsetBySubdomain::init(), libMesh::PetscDMWrapper::init_and_attach_petscdm(), libMesh::Nemesis_IO_Helper::initialize(), libMesh::DistributedMesh::insert_elem(), libMesh::MeshTools::libmesh_assert_contiguous_dof_ids(), libMesh::MeshTools::libmesh_assert_parallel_consistent_new_node_procids(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_topology_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_boundary_ids(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::DofMap::local_variable_indices(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshBase::n_active_elem_on_proc(), libMesh::MeshBase::n_elem_on_proc(), libMesh::MeshBase::n_nodes_on_proc(), libMesh::MeshBase::partition(), libMesh::PetscLinearSolver< T >::PetscLinearSolver(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::NameBasedIO::read(), libMesh::Nemesis_IO::read(), libMesh::CheckpointIO::read(), libMesh::CheckpointIO::read_connectivity(), libMesh::XdrIO::read_header(), libMesh::CheckpointIO::read_nodes(), libMesh::System::read_parallel_data(), libMesh::System::read_SCALAR_dofs(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::System::read_serialized_vector(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::GMVIO::write_binary(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::System::write_parallel_data(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), and libMesh::XdrIO::write_serialized_nodesets().

96  { return cast_int<processor_id_type>(_communicator.size()); }
processor_id_type size() const
Definition: communicator.h:175
const Parallel::Communicator & _communicator

◆ operator()() [1/4]

void libMesh::FunctionBase< Number >::operator() ( const Point p,
DenseVector< Number > &  output 
)
inlineinherited

Evaluation function for time-independent vector-valued functions. Sets output values in the passed-in output DenseVector.

Definition at line 241 of file function_base.h.

243 {
244  // Call the time-dependent function with t=0.
245  this->operator()(p, 0., output);
246 }
virtual Number operator()(const Point &p, const Real time=0.)=0

◆ operator()() [2/4]

Number libMesh::MeshFunction::operator() ( const Point p,
const Real  time = 0. 
)
overridevirtual
Returns
The value of variable 0 at point p and for time, which defaults to zero.

Implements libMesh::FunctionBase< Number >.

Definition at line 174 of file mesh_function.C.

References libMesh::FunctionBase< Number >::initialized().

Referenced by operator()().

176 {
177  libmesh_assert (this->initialized());
178 
179  DenseVector<Number> buf (1);
180  this->operator() (p, time, buf);
181  return buf(0);
182 }
Number operator()(const Point &p, const Real time=0.) override

◆ operator()() [3/4]

void libMesh::MeshFunction::operator() ( const Point p,
const Real  time,
DenseVector< Number > &  output 
)
overridevirtual

Computes values at coordinate p and for time time, which defaults to zero, optionally restricting the point to the passed subdomain_ids. This is useful in cases where there are multiple dimensioned elements, for example.

Implements libMesh::FunctionBase< Number >.

Definition at line 236 of file mesh_function.C.

References operator()().

239 {
240  this->operator() (p, time, output, nullptr);
241 }
Number operator()(const Point &p, const Real time=0.) override

◆ operator()() [4/4]

void libMesh::MeshFunction::operator() ( const Point p,
const Real  time,
DenseVector< Number > &  output,
const std::set< subdomain_id_type > *  subdomain_ids 
)

Computes values at coordinate p and for time time, restricting the point to the passed subdomain_ids. This is useful in cases where there are multiple dimensioned elements, for example.

Build an FEComputeData that contains both input and output data for the specific compute_data method.

Definition at line 243 of file mesh_function.C.

References _dof_map, _eqn_systems, _out_of_mesh_mode, _out_of_mesh_value, _system_vars, _vector, libMesh::FEInterface::compute_data(), data, libMesh::Elem::dim(), libMesh::DofMap::dof_indices(), find_element(), libMesh::FunctionBase< Number >::initialized(), libMesh::invalid_uint, libMesh::FEInterface::inverse_map(), libMesh::DenseVector< T >::resize(), libMesh::DenseVector< T >::size(), value, and libMesh::DofMap::variable_type().

247 {
248  libmesh_assert (this->initialized());
249 
250  const Elem * element = this->find_element(p,subdomain_ids);
251 
252  if (!element)
253  {
254  // We'd better be in out_of_mesh_mode if we couldn't find an
255  // element in the mesh
256  libmesh_assert (_out_of_mesh_mode);
257  output = _out_of_mesh_value;
258  }
259  else
260  {
261  // resize the output vector to the number of output values
262  // that the user told us
263  output.resize (cast_int<unsigned int>
264  (this->_system_vars.size()));
265 
266 
267  {
268  const unsigned int dim = element->dim();
269 
270 
271  /*
272  * Get local coordinates to feed these into compute_data().
273  * Note that the fe_type can safely be used from the 0-variable,
274  * since the inverse mapping is the same for all FEFamilies
275  */
276  const Point mapped_point (FEInterface::inverse_map (dim,
277  this->_dof_map.variable_type(0),
278  element,
279  p));
280 
281 
282  // loop over all vars
283  for (unsigned int index=0,
284  sz = cast_int<unsigned int>(this->_system_vars.size());
285  index != sz; ++index)
286  {
287  /*
288  * the data for this variable
289  */
290  const unsigned int var = _system_vars[index];
291 
292  if (var == libMesh::invalid_uint)
293  {
294  libmesh_assert (_out_of_mesh_mode &&
295  index < _out_of_mesh_value.size());
296  output(index) = _out_of_mesh_value(index);
297  continue;
298  }
299 
300  const FEType & fe_type = this->_dof_map.variable_type(var);
301 
306  {
307  FEComputeData data (this->_eqn_systems, mapped_point);
308 
309  FEInterface::compute_data (dim, fe_type, element, data);
310 
311  // where the solution values for the var-th variable are stored
312  std::vector<dof_id_type> dof_indices;
313  this->_dof_map.dof_indices (element, dof_indices, var);
314 
315  // interpolate the solution
316  {
317  Number value = 0.;
318 
319  for (std::size_t i=0; i<dof_indices.size(); i++)
320  value += this->_vector(dof_indices[i]) * data.shape[i];
321 
322  output(index) = value;
323  }
324 
325  }
326 
327  // next variable
328  }
329  }
330  }
331 }
virtual unsigned int size() const override
Definition: dense_vector.h:92
const EquationSystems & _eqn_systems
const NumericVector< Number > & _vector
const unsigned int invalid_uint
Definition: libmesh.h:245
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Definition: dof_map.C:1930
const DofMap & _dof_map
void resize(const unsigned int n)
Definition: dense_vector.h:355
const Elem * find_element(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1792
DenseVector< Number > _out_of_mesh_value
static void compute_data(const unsigned int dim, const FEType &fe_t, const Elem *elem, FEComputeData &data)
Definition: fe_interface.C:837
static Point inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: fe_interface.C:590
static const bool value
Definition: xdr_io.C:109
IterBase * data
const std::vector< unsigned int > _system_vars

◆ operator=() [1/2]

MeshFunction& libMesh::MeshFunction::operator= ( const MeshFunction )
delete

This class contains const references so it can't be assigned.

◆ operator=() [2/2]

MeshFunction& libMesh::MeshFunction::operator= ( MeshFunction &&  )
delete

◆ processor_id()

processor_id_type libMesh::ParallelObject::processor_id ( ) const
inlineinherited
Returns
The rank of this processor in the group.

Definition at line 101 of file parallel_object.h.

References libMesh::ParallelObject::_communicator, and libMesh::Parallel::Communicator::rank().

Referenced by libMesh::BoundaryInfo::_find_id_maps(), libMesh::EquationSystems::_read_impl(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::DistributedMesh::add_elem(), libMesh::BoundaryInfo::add_elements(), libMesh::DofMap::add_neighbors_to_send_list(), libMesh::DistributedMesh::add_node(), libMesh::UnstructuredMesh::all_second_order(), libMesh::MeshTools::Modification::all_tri(), libMesh::FEMSystem::assembly(), libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::Nemesis_IO_Helper::build_element_and_node_maps(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::EquationSystems::build_parallel_elemental_solution_vector(), libMesh::DistributedMesh::clear(), libMesh::ExodusII_IO_Helper::close(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::compute_communication_map_parameters(), libMesh::Nemesis_IO_Helper::compute_internal_and_border_elems_and_internal_nodes(), libMesh::Nemesis_IO_Helper::compute_node_communication_maps(), libMesh::Nemesis_IO_Helper::compute_num_global_elem_blocks(), libMesh::Nemesis_IO_Helper::compute_num_global_nodesets(), libMesh::Nemesis_IO_Helper::compute_num_global_sidesets(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::ExodusII_IO_Helper::create(), libMesh::DistributedMesh::delete_elem(), libMesh::DistributedMesh::delete_node(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DistributedMesh::DistributedMesh(), libMesh::DofMap::end_dof(), libMesh::DofMap::end_old_dof(), libMesh::EnsightIO::EnsightIO(), find_element(), find_elements(), libMesh::UnstructuredMesh::find_neighbors(), libMesh::DofMap::first_dof(), libMesh::DofMap::first_old_dof(), libMesh::Nemesis_IO_Helper::get_cmap_params(), libMesh::Nemesis_IO_Helper::get_eb_info_global(), libMesh::Nemesis_IO_Helper::get_elem_cmap(), libMesh::Nemesis_IO_Helper::get_elem_map(), libMesh::MeshBase::get_info(), libMesh::DofMap::get_info(), libMesh::Nemesis_IO_Helper::get_init_global(), libMesh::Nemesis_IO_Helper::get_init_info(), libMesh::Nemesis_IO_Helper::get_loadbal_param(), libMesh::Nemesis_IO_Helper::get_node_cmap(), libMesh::Nemesis_IO_Helper::get_node_map(), libMesh::Nemesis_IO_Helper::get_ns_param_global(), libMesh::Nemesis_IO_Helper::get_ss_param_global(), libMesh::SparsityPattern::Build::handle_vi_vj(), libMesh::SystemSubsetBySubdomain::init(), libMesh::ExodusII_IO_Helper::initialize(), libMesh::ExodusII_IO_Helper::initialize_element_variables(), libMesh::ExodusII_IO_Helper::initialize_global_variables(), libMesh::ExodusII_IO_Helper::initialize_nodal_variables(), libMesh::DistributedMesh::insert_elem(), libMesh::DofMap::is_evaluable(), libMesh::SparsityPattern::Build::join(), libMesh::DofMap::last_dof(), libMesh::MeshTools::libmesh_assert_consistent_distributed(), libMesh::MeshTools::libmesh_assert_consistent_distributed_nodes(), libMesh::MeshTools::libmesh_assert_contiguous_dof_ids(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_object_ids(), libMesh::DofMap::local_variable_indices(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshBase::n_active_local_elem(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::DofMap::n_local_dofs(), libMesh::System::n_local_dofs(), libMesh::MeshBase::n_local_elem(), libMesh::MeshBase::n_local_nodes(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::BoundaryInfo::n_shellface_conds(), libMesh::SparsityPattern::Build::operator()(), libMesh::DistributedMesh::own_node(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::Nemesis_IO_Helper::put_cmap_params(), libMesh::Nemesis_IO_Helper::put_elem_cmap(), libMesh::Nemesis_IO_Helper::put_elem_map(), libMesh::Nemesis_IO_Helper::put_loadbal_param(), libMesh::Nemesis_IO_Helper::put_node_cmap(), libMesh::Nemesis_IO_Helper::put_node_map(), libMesh::NameBasedIO::read(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::CheckpointIO::read(), libMesh::ExodusII_IO_Helper::read_elem_num_map(), libMesh::ExodusII_IO_Helper::read_global_values(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::System::read_header(), libMesh::System::read_legacy_data(), libMesh::ExodusII_IO_Helper::read_node_num_map(), libMesh::System::read_parallel_data(), libMesh::System::read_SCALAR_dofs(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::System::read_serialized_data(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::System::read_serialized_vector(), libMesh::System::read_serialized_vectors(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::CheckpointIO::select_split_config(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::MeshTools::total_weight(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::Parallel::Packing< T >::unpack(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), libMesh::CheckpointIO::write(), libMesh::EquationSystems::write(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::ExodusII_IO::write_element_data(), libMesh::ExodusII_IO_Helper::write_element_values(), libMesh::ExodusII_IO_Helper::write_elements(), libMesh::ExodusII_IO::write_global_data(), libMesh::ExodusII_IO_Helper::write_global_values(), libMesh::System::write_header(), libMesh::ExodusII_IO::write_information_records(), libMesh::ExodusII_IO_Helper::write_information_records(), libMesh::ExodusII_IO_Helper::write_nodal_coordinates(), libMesh::UCDIO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), libMesh::ExodusII_IO_Helper::write_nodal_values(), libMesh::Nemesis_IO_Helper::write_nodesets(), libMesh::ExodusII_IO_Helper::write_nodesets(), libMesh::System::write_parallel_data(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bc_names(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::System::write_serialized_data(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::XdrIO::write_serialized_subdomain_names(), libMesh::System::write_serialized_vector(), libMesh::System::write_serialized_vectors(), libMesh::Nemesis_IO_Helper::write_sidesets(), libMesh::ExodusII_IO_Helper::write_sidesets(), libMesh::ExodusII_IO::write_timestep(), libMesh::ExodusII_IO_Helper::write_timestep(), and libMesh::ExodusII_IO::write_timestep_discontinuous().

102  { return cast_int<processor_id_type>(_communicator.rank()); }
const Parallel::Communicator & _communicator
processor_id_type rank() const
Definition: communicator.h:173

◆ set_is_time_dependent()

void libMesh::FunctionBase< Number >::set_is_time_dependent ( bool  is_time_dependent)
inlineinherited

Function to set whether this is a time-dependent function or not. This is intended to be only used by subclasses who cannot natively determine time-dependence. In such a case, this function should be used immediately following construction.

Definition at line 213 of file function_base.h.

◆ set_point_locator_tolerance()

void libMesh::MeshFunction::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 to evaluate at might be slightly outside the mesh (due to numerical rounding issues, for example).

Definition at line 795 of file mesh_function.C.

References _point_locator, and libMesh::PointLocatorBase::set_close_to_point_tol().

796 {
798 }
virtual void set_close_to_point_tol(Real close_to_point_tol)
PointLocatorBase * _point_locator

◆ unset_point_locator_tolerance()

void libMesh::MeshFunction::unset_point_locator_tolerance ( )

Turn off the user-specified PointLocator tolerance.

Definition at line 800 of file mesh_function.C.

References _point_locator, and libMesh::PointLocatorBase::unset_close_to_point_tol().

801 {
803 }
PointLocatorBase * _point_locator
virtual void unset_close_to_point_tol()

Member Data Documentation

◆ _communicator

◆ _dof_map

const DofMap& libMesh::MeshFunction::_dof_map
protected

Need access to the DofMap of the other system.

Definition at line 339 of file mesh_function.h.

Referenced by clone(), discontinuous_gradient(), discontinuous_value(), gradient(), hessian(), and operator()().

◆ _eqn_systems

const EquationSystems& libMesh::MeshFunction::_eqn_systems
protected

The equation systems handler, from which the data are gathered.

Definition at line 328 of file mesh_function.h.

Referenced by clone(), discontinuous_value(), init(), and operator()().

◆ _initialized

bool libMesh::FunctionBase< Number >::_initialized
protectedinherited

When init() was called so that everything is ready for calls to operator() (...), then this bool is true.

Definition at line 180 of file function_base.h.

Referenced by clear(), and init().

◆ _is_time_dependent

bool libMesh::FunctionBase< Number >::_is_time_dependent
protectedinherited

Cache whether or not this function is actually time-dependent.

Definition at line 185 of file function_base.h.

◆ _master

const FunctionBase* libMesh::FunctionBase< Number >::_master
protectedinherited

Const pointer to our master, initialized to nullptr. There may be cases where multiple functions are required, but to save memory, one master handles some centralized data.

Definition at line 174 of file function_base.h.

Referenced by clear(), find_element(), find_elements(), init(), and ~MeshFunction().

◆ _out_of_mesh_mode

bool libMesh::MeshFunction::_out_of_mesh_mode
protected

true if out-of-mesh mode is enabled. See enable_out_of_mesh_mode() for more details. Default is false.

Definition at line 357 of file mesh_function.h.

Referenced by disable_out_of_mesh_mode(), discontinuous_gradient(), discontinuous_value(), enable_out_of_mesh_mode(), find_element(), find_elements(), gradient(), hessian(), and operator()().

◆ _out_of_mesh_value

DenseVector<Number> libMesh::MeshFunction::_out_of_mesh_value
protected

Value to return outside the mesh if out-of-mesh mode is enabled. See enable_out_of_mesh_mode() for more details.

Definition at line 363 of file mesh_function.h.

Referenced by discontinuous_gradient(), discontinuous_value(), enable_out_of_mesh_mode(), gradient(), hessian(), and operator()().

◆ _point_locator

PointLocatorBase* libMesh::MeshFunction::_point_locator
protected

◆ _system_vars

const std::vector<unsigned int> libMesh::MeshFunction::_system_vars
protected

The indices of the variables within the other system for which data are to be gathered.

Definition at line 345 of file mesh_function.h.

Referenced by clone(), discontinuous_gradient(), discontinuous_value(), gradient(), hessian(), init(), and operator()().

◆ _vector

const NumericVector<Number>& libMesh::MeshFunction::_vector
protected

A reference to the vector that holds the data that is to be interpolated.

Definition at line 334 of file mesh_function.h.

Referenced by clone(), discontinuous_value(), find_element(), find_elements(), and operator()().


The documentation for this class was generated from the following files: