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=libmesh_nullptr)
 
 MeshFunction (const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, const unsigned int var, const FunctionBase< Number > *master=libmesh_nullptr)
 
 ~MeshFunction ()
 
virtual void init () libmesh_override
 
void init (const Trees::BuildType point_locator_build_type)
 
virtual void clear () libmesh_override
 
virtual std::unique_ptr< FunctionBase< Number > > clone () const libmesh_override
 
Number operator() (const Point &p, const Real time=0.) libmesh_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) libmesh_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=libmesh_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=libmesh_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=libmesh_nullptr) const
 
std::set< const Elem * > find_elements (const Point &p, const std::set< subdomain_id_type > *subdomain_ids=libmesh_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

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 = libmesh_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),
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
const class libmesh_nullptr_t libmesh_nullptr
DenseVector< Number > _out_of_mesh_value
const std::vector< unsigned int > _system_vars
libMesh::MeshFunction::MeshFunction ( const EquationSystems eqn_systems,
const NumericVector< Number > &  vec,
const DofMap dof_map,
const unsigned int  var,
const FunctionBase< Number > *  master = libmesh_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),
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
const class libmesh_nullptr_t libmesh_nullptr
DenseVector< Number > _out_of_mesh_value
const std::vector< unsigned int > _system_vars
libMesh::MeshFunction::~MeshFunction ( )

Destructor.

Definition at line 88 of file mesh_function.C.

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

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

Member Function Documentation

void libMesh::MeshFunction::clear ( )
virtual

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, _point_locator, and libmesh_nullptr.

Referenced by init().

149 {
150  // only delete the point locator when we are the master
151  if ((this->_point_locator != libmesh_nullptr) && (this->_master == libmesh_nullptr))
152  {
153  delete this->_point_locator;
155  }
156  this->_initialized = false;
157 }
PointLocatorBase * _point_locator
const FunctionBase * _master
const class libmesh_nullptr_t libmesh_nullptr
std::unique_ptr< FunctionBase< Number > > libMesh::MeshFunction::clone ( ) const
virtual
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().

Referenced by init().

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
const DofMap & _dof_map
MeshFunction(const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, const std::vector< unsigned int > &vars, const FunctionBase< Number > *master=libmesh_nullptr)
Definition: mesh_function.C:43
const std::vector< unsigned int > _system_vars
const Parallel::Communicator& libMesh::ParallelObject::comm ( ) const
inlineinherited
Returns
A reference to the Parallel::Communicator object used by this mesh.

Definition at line 87 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_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_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::ParmetisPartitioner::_do_repartition(), libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult(), libMesh::SlepcEigenSolver< 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::ImplicitSystem::add_matrix(), libMesh::System::add_vector(), libMesh::EigenSparseLinearSolver< T >::adjoint_solve(), libMesh::UnstructuredMesh::all_second_order(), libMesh::MeshTools::Modification::all_tri(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assemble_qoi(), libMesh::MeshCommunication::assign_global_indices(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::DofMap::attach_matrix(), libMesh::Parallel::BinSorter< KeyType, IdxType >::binsort(), libMesh::Parallel::Sort< KeyType, IdxType >::binsort(), libMesh::MeshCommunication::broadcast(), libMesh::SparseMatrix< T >::build(), libMesh::MeshTools::Generation::build_extrusion(), libMesh::Parallel::Histogram< KeyType, IdxType >::build_histogram(), libMesh::PetscNonlinearSolver< T >::build_mat_null_space(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::EquationSystems::build_parallel_solution_vector(), libMesh::MeshBase::cache_elem_dims(), libMesh::System::calculate_norm(), libMesh::DofMap::check_dirichlet_bcid_consistency(), libMesh::DistributedVector< T >::clone(), libMesh::EigenSparseVector< T >::clone(), libMesh::LaspackVector< T >::clone(), libMesh::EpetraVector< T >::clone(), libMesh::PetscVector< T >::clone(), libMesh::EpetraVector< T >::close(), libMesh::Parallel::Sort< KeyType, IdxType >::communicate_bins(), 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::MeshCommunication::gather(), libMesh::MeshCommunication::gather_neighboring_elements(), libMesh::CondensedEigenSystem::get_eigenpair(), libMesh::DofMap::get_info(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::EquationSystems::get_solution(), libMesh::LocationMap< T >::init(), libMesh::TopologyMap::init(), libMesh::TimeSolver::init(), libMesh::TaoOptimizationSolver< T >::init(), libMesh::PetscNonlinearSolver< T >::init(), libMesh::DistributedVector< T >::init(), libMesh::EpetraVector< T >::init(), libMesh::PetscVector< T >::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::EigenSystem::init_data(), libMesh::EigenSystem::init_matrices(), libMesh::ParmetisPartitioner::initialize(), libMesh::OptimizationSystem::initialize_equality_constraints_storage(), libMesh::OptimizationSystem::initialize_inequality_constraints_storage(), 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_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::DistributedVector< T >::max(), libMesh::FEMSystem::mesh_position_set(), libMesh::MeshSerializer::MeshSerializer(), libMesh::DistributedVector< T >::min(), 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::Partitioner::partition(), libMesh::LinearPartitioner::partition_range(), libMesh::MetisPartitioner::partition_range(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::petsc_auto_fieldsplit(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::MeshBase::prepare_for_use(), libMesh::SparseMatrix< T >::print(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::MeshBase::recalculate_n_partitions(), libMesh::MeshCommunication::redistribute(), libMesh::MeshRefinement::refine_and_coarsen_elements(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::CheckpointIO::select_split_config(), libMesh::MeshCommunication::send_coarse_ghosts(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::Partitioner::set_parent_processor_ids(), libMesh::PetscDiffSolver::setup_petsc_data(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::Parallel::Sort< KeyType, IdxType >::sort(), libMesh::split_mesh(), libMesh::MeshBase::subdomain_ids(), libMesh::BoundaryInfo::sync(), libMesh::Parallel::sync_element_data_by_parent_id(), libMesh::Parallel::sync_node_data_by_element_id(), libMesh::MeshRefinement::test_level_one(), libMesh::MeshRefinement::test_unflagged(), libMesh::MeshTools::total_weight(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::DistributedVector< T >::zero_clone(), libMesh::EigenSparseVector< T >::zero_clone(), libMesh::LaspackVector< T >::zero_clone(), libMesh::EpetraVector< T >::zero_clone(), and libMesh::PetscVector< T >::zero_clone().

88  { return _communicator; }
const Parallel::Communicator & _communicator
virtual Number libMesh::FunctionBase< Number >::component ( unsigned int  i,
const Point p,
Real  time = 0. 
)
virtualinherited
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.
void libMesh::MeshFunction::disable_out_of_mesh_mode ( void  )

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

Definition at line 793 of file mesh_function.C.

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

Referenced by init().

794 {
795  libmesh_assert (this->initialized());
797  _out_of_mesh_mode = false;
798 }
virtual void disable_out_of_mesh_mode()=0
PointLocatorBase * _point_locator
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(), and init().

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 (std::map<const Elem *, std::vector<Gradient>>::const_iterator it = buffer.begin(); it != buffer.end(); ++it)
218  return_value[it->first] = it->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.)
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 504 of file mesh_function.C.

References discontinuous_gradient(), and libmesh_nullptr.

507 {
508  this->discontinuous_gradient (p, time, output, libmesh_nullptr);
509 }
const class libmesh_nullptr_t libmesh_nullptr
std::map< const Elem *, Gradient > discontinuous_gradient(const Point &p, const Real time=0.)
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 513 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_elements(), libMesh::FunctionBase< Number >::initialized(), libMesh::invalid_uint, libMesh::FEInterface::inverse_map(), libMesh::DenseVector< T >::size(), and libMesh::DofMap::variable_type().

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

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 (std::map<const Elem *, DenseVector<Number>>::const_iterator it = buffer.begin(); it != buffer.end(); ++it)
193  return_value[it->first] = it->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.)
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 332 of file mesh_function.C.

References discontinuous_value(), and libmesh_nullptr.

335 {
336  this->discontinuous_value (p, time, output, libmesh_nullptr);
337 }
const class libmesh_nullptr_t libmesh_nullptr
std::map< const Elem *, Number > discontinuous_value(const Point &p, const Real time=0.)
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 341 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_elements(), libMesh::FunctionBase< Number >::initialized(), libMesh::invalid_uint, libMesh::FEInterface::inverse_map(), libMesh::FEComputeData::shape, libMesh::DenseVector< T >::size(), value, and libMesh::DofMap::variable_type().

345 {
346  libmesh_assert (this->initialized());
347 
348  // clear the output map
349  output.clear();
350 
351  // get the candidate elements
352  std::set<const Elem *> candidate_element = this->find_elements(p,subdomain_ids);
353 
354  // loop through all candidates, if the set is empty this function will return an
355  // empty map
356  for (std::set<const Elem *>::const_iterator it = candidate_element.begin(); it != candidate_element.end(); ++it)
357  {
358  const Elem * element = (*it);
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 (std::size_t index=0; index < this->_system_vars.size(); index++)
378  {
379  /*
380  * the data for this variable
381  */
382  const unsigned int var = _system_vars[index];
383 
384  if (var == libMesh::invalid_uint)
385  {
386  libmesh_assert (_out_of_mesh_mode &&
387  index < _out_of_mesh_value.size());
388  temp_output(index) = _out_of_mesh_value(index);
389  continue;
390  }
391 
392  const FEType & fe_type = this->_dof_map.variable_type(var);
393 
398  {
399  FEComputeData data (this->_eqn_systems, mapped_point);
400 
401  FEInterface::compute_data (dim, fe_type, element, data);
402 
403  // where the solution values for the var-th variable are stored
404  std::vector<dof_id_type> dof_indices;
405  this->_dof_map.dof_indices (element, dof_indices, var);
406 
407  // interpolate the solution
408  {
409  Number value = 0.;
410 
411  for (std::size_t i=0; i<dof_indices.size(); i++)
412  value += this->_vector(dof_indices[i]) * data.shape[i];
413 
414  temp_output(index) = value;
415  }
416 
417  }
418 
419  // next variable
420  }
421 
422  // Insert temp_output into output
423  output[element] = temp_output;
424  }
425 }
const EquationSystems & _eqn_systems
const NumericVector< Number > & _vector
const unsigned int invalid_uint
Definition: libmesh.h:184
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1719
const DofMap & _dof_map
virtual unsigned int size() const libmesh_override
Definition: dense_vector.h:87
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:794
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:547
std::set< const Elem * > find_elements(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=libmesh_nullptr) const
static const bool value
Definition: xdr_io.C:108
IterBase * data
const std::vector< unsigned int > _system_vars
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Definition: dof_map.C:1924
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 778 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(), and init().

779 {
780  libmesh_assert (this->initialized());
782  _out_of_mesh_mode = true;
784 }
PointLocatorBase * _point_locator
DenseVector< Number > _out_of_mesh_value
static const bool value
Definition: xdr_io.C:108
virtual void enable_out_of_mesh_mode()=0
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 786 of file mesh_function.C.

References enable_out_of_mesh_mode(), and value.

787 {
788  DenseVector<Number> v(1);
789  v(0) = value;
790  this->enable_out_of_mesh_mode(v);
791 }
void enable_out_of_mesh_mode(const DenseVector< Number > &value)
static const bool value
Definition: xdr_io.C:108
const Elem * libMesh::MeshFunction::find_element ( const Point p,
const std::set< subdomain_id_type > *  subdomain_ids = libmesh_nullptr 
) const
protected

Helper function to reduce code duplication

Definition at line 668 of file mesh_function.C.

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

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

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

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

Referenced by discontinuous_gradient(), discontinuous_value(), and init().

719 {
720  /* Ensure that in the case of a master mesh function, the
721  out-of-mesh mode is enabled either for both or for none. This is
722  important because the out-of-mesh mode is also communicated to
723  the point locator. Since this is time consuming, enable it only
724  in debug mode. */
725 #ifdef DEBUG
726  if (this->_master != libmesh_nullptr)
727  {
728  const MeshFunction * master =
729  cast_ptr<const MeshFunction *>(this->_master);
730  if (_out_of_mesh_mode!=master->_out_of_mesh_mode)
731  libmesh_error_msg("ERROR: If you use out-of-mesh-mode in connection with master mesh " \
732  << "functions, you must enable out-of-mesh mode for both the master and the slave mesh function.");
733  }
734 #endif
735 
736  // locate the point in the other mesh
737  std::set<const Elem *> candidate_elements;
738  std::set<const Elem *> final_candidate_elements;
739  this->_point_locator->operator()(p,candidate_elements,subdomain_ids);
740  for (std::set<const Elem *>::const_iterator elem_it = candidate_elements.begin(); elem_it != candidate_elements.end(); ++elem_it)
741  {
742  const Elem * element = (*elem_it);
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  std::set<const Elem *>::const_iterator it = point_neighbors.begin();
754  const std::set<const Elem *>::const_iterator end = point_neighbors.end();
755  for (; it != end; ++it)
756  {
757  const Elem * elem = *it;
758  if (elem->processor_id() == this->processor_id())
759  {
760  final_candidate_elements.insert(elem);
761  break;
762  }
763  }
764  }
765  else
766  final_candidate_elements.insert(element);
767  }
768 
769  return final_candidate_elements;
770 }
const NumericVector< Number > & _vector
PointLocatorBase * _point_locator
const FunctionBase * _master
const class libmesh_nullptr_t libmesh_nullptr
IterBase * end
ParallelType type() const
MeshFunction(const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, const std::vector< unsigned int > &vars, const FunctionBase< Number > *master=libmesh_nullptr)
Definition: mesh_function.C:43
processor_id_type processor_id() const
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 772 of file mesh_function.C.

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

Referenced by init().

773 {
774  libmesh_assert (this->initialized());
775  return *_point_locator;
776 }
PointLocatorBase * _point_locator
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(), and init().

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.)
void libMesh::MeshFunction::gradient ( const Point p,
const Real  time,
std::vector< Gradient > &  output,
const std::set< subdomain_id_type > *  subdomain_ids = libmesh_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 429 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().

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

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.)
void libMesh::MeshFunction::hessian ( const Point p,
const Real  time,
std::vector< Tensor > &  output,
const std::set< subdomain_id_type > *  subdomain_ids = libmesh_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 593 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().

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

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 96 of file mesh_function.h.

References clear(), clone(), disable_out_of_mesh_mode(), discontinuous_gradient(), discontinuous_value(), enable_out_of_mesh_mode(), find_element(), find_elements(), get_point_locator(), gradient(), hessian(), init(), libmesh_nullptr, libMesh::Trees::NODES, operator()(), libMesh::Real, set_point_locator_tolerance(), unset_point_locator_tolerance(), and value.

Referenced by init().

96 { this->init(Trees::NODES); }
virtual void init() libmesh_override
Definition: mesh_function.h:96
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(), libmesh_nullptr, 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 != libmesh_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 == libmesh_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
MeshBase & mesh
const FunctionBase * _master
const class libmesh_nullptr_t libmesh_nullptr
const MeshBase & get_mesh() const
MeshFunction(const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, const std::vector< unsigned int > &vars, const FunctionBase< Number > *master=libmesh_nullptr)
Definition: mesh_function.C:43
const std::vector< unsigned int > _system_vars
bool libMesh::FunctionBase< Number >::initialized ( ) const
inherited
Returns
true when this object is properly initialized and ready for use, false otherwise.

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

bool libMesh::FunctionBase< Number >::is_time_dependent ( ) const
inherited
Returns
true when the function this object represents is actually time-dependent, false otherwise.
processor_id_type libMesh::ParallelObject::n_processors ( ) const
inlineinherited
Returns
The number of processors in the group.

Definition at line 93 of file parallel_object.h.

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

Referenced by libMesh::ParmetisPartitioner::_do_repartition(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::DistributedMesh::add_elem(), libMesh::DistributedMesh::add_node(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assembly(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::AztecLinearSolver< T >::AztecLinearSolver(), libMesh::MeshCommunication::broadcast(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::DistributedMesh::clear(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::MeshTools::correct_node_proc_ids(), 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::DistributedMesh::DistributedMesh(), libMesh::EnsightIO::EnsightIO(), libMesh::MeshCommunication::gather(), libMesh::MeshCommunication::gather_neighboring_elements(), libMesh::MeshBase::get_info(), libMesh::EquationSystems::get_solution(), libMesh::DistributedVector< T >::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::ParmetisPartitioner::initialize(), libMesh::Nemesis_IO_Helper::initialize(), libMesh::DistributedMesh::insert_elem(), 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::SparsityPattern::Build::parallel_sync(), libMesh::Partitioner::partition(), libMesh::MeshBase::partition(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::PetscLinearSolver< T >::PetscLinearSolver(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::SparseMatrix< T >::print(), libMesh::NameBasedIO::read(), libMesh::Nemesis_IO::read(), libMesh::CheckpointIO::read(), libMesh::CheckpointIO::read_connectivity(), libMesh::XdrIO::read_header(), libMesh::CheckpointIO::read_nodes(), libMesh::MeshCommunication::redistribute(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::Partitioner::repartition(), libMesh::MeshCommunication::send_coarse_ghosts(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::Parallel::Sort< KeyType, IdxType >::sort(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::GMVIO::write_binary(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), and libMesh::XdrIO::write_serialized_nodesets().

94  { return cast_int<processor_id_type>(_communicator.size()); }
unsigned int size() const
Definition: parallel.h:725
const Parallel::Communicator & _communicator
void libMesh::FunctionBase< Number >::operator() ( const Point p,
DenseVector< Number > &  output 
)
inherited

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

Number libMesh::MeshFunction::operator() ( const Point p,
const Real  time = 0. 
)
virtual
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 init(), and 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.) libmesh_override
void libMesh::MeshFunction::operator() ( const Point p,
const Real  time,
DenseVector< Number > &  output 
)
virtual

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 libmesh_nullptr, and operator()().

239 {
240  this->operator() (p, time, output, libmesh_nullptr);
241 }
const class libmesh_nullptr_t libmesh_nullptr
Number operator()(const Point &p, const Real time=0.) libmesh_override
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::FEComputeData::shape, 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 (std::size_t index=0; index < this->_system_vars.size(); index++)
284  {
285  /*
286  * the data for this variable
287  */
288  const unsigned int var = _system_vars[index];
289 
290  if (var == libMesh::invalid_uint)
291  {
292  libmesh_assert (_out_of_mesh_mode &&
293  index < _out_of_mesh_value.size());
294  output(index) = _out_of_mesh_value(index);
295  continue;
296  }
297 
298  const FEType & fe_type = this->_dof_map.variable_type(var);
299 
304  {
305  FEComputeData data (this->_eqn_systems, mapped_point);
306 
307  FEInterface::compute_data (dim, fe_type, element, data);
308 
309  // where the solution values for the var-th variable are stored
310  std::vector<dof_id_type> dof_indices;
311  this->_dof_map.dof_indices (element, dof_indices, var);
312 
313  // interpolate the solution
314  {
315  Number value = 0.;
316 
317  for (std::size_t i=0; i<dof_indices.size(); i++)
318  value += this->_vector(dof_indices[i]) * data.shape[i];
319 
320  output(index) = value;
321  }
322 
323  }
324 
325  // next variable
326  }
327  }
328  }
329 }
const Elem * find_element(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=libmesh_nullptr) const
const EquationSystems & _eqn_systems
const NumericVector< Number > & _vector
const unsigned int invalid_uint
Definition: libmesh.h:184
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1719
const DofMap & _dof_map
void resize(const unsigned int n)
Definition: dense_vector.h:350
virtual unsigned int size() const libmesh_override
Definition: dense_vector.h:87
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:794
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:547
static const bool value
Definition: xdr_io.C:108
IterBase * data
const std::vector< unsigned int > _system_vars
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Definition: dof_map.C:1924
processor_id_type libMesh::ParallelObject::processor_id ( ) const
inlineinherited
Returns
The rank of this processor in the group.

Definition at line 99 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::DistributedMesh::add_elem(), libMesh::BoundaryInfo::add_elements(), libMesh::DofMap::add_neighbors_to_send_list(), libMesh::DistributedMesh::add_node(), libMesh::MeshRefinement::add_node(), libMesh::MeshTools::Modification::all_tri(), libMesh::FEMSystem::assembly(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::MeshCommunication::broadcast(), libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::Nemesis_IO_Helper::build_element_and_node_maps(), libMesh::ParmetisPartitioner::build_graph(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::DofMap::build_sparsity(), 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::EnsightIO::EnsightIO(), find_element(), find_elements(), libMesh::UnstructuredMesh::find_neighbors(), libMesh::MeshCommunication::gather(), libMesh::MeshCommunication::gather_neighboring_elements(), 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::EquationSystems::get_solution(), libMesh::Nemesis_IO_Helper::get_ss_param_global(), libMesh::SparsityPattern::Build::handle_vi_vj(), libMesh::DistributedVector< T >::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::ParmetisPartitioner::initialize(), 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::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), 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::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::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::SparsityPattern::Build::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::SparsityPattern::Build::parallel_sync(), libMesh::MetisPartitioner::partition_range(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::SparseMatrix< T >::print(), libMesh::NumericVector< T >::print_global(), 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::ExodusII_IO_Helper::read_node_num_map(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::MeshCommunication::redistribute(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::CheckpointIO::select_split_config(), libMesh::MeshCommunication::send_coarse_ghosts(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::MeshTools::total_weight(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::Parallel::Packing< Node * >::unpack(), libMesh::Parallel::Packing< Elem * >::unpack(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::MeshTools::weight(), 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::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::XdrIO::write_serialized_bc_names(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::XdrIO::write_serialized_subdomain_names(), 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().

100  { return cast_int<processor_id_type>(_communicator.rank()); }
const Parallel::Communicator & _communicator
unsigned int rank() const
Definition: parallel.h:723
void libMesh::FunctionBase< Number >::set_is_time_dependent ( bool  is_time_dependent)
inherited

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.

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 800 of file mesh_function.C.

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

Referenced by init().

801 {
803 }
virtual void set_close_to_point_tol(Real close_to_point_tol)
PointLocatorBase * _point_locator
void libMesh::MeshFunction::unset_point_locator_tolerance ( )

Turn off the user-specified PointLocator tolerance.

Definition at line 805 of file mesh_function.C.

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

Referenced by init().

806 {
808 }
PointLocatorBase * _point_locator
virtual void unset_close_to_point_tol()

Member Data Documentation

const DofMap& libMesh::MeshFunction::_dof_map
protected

Need access to the DofMap of the other system.

Definition at line 325 of file mesh_function.h.

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

const EquationSystems& libMesh::MeshFunction::_eqn_systems
protected

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

Definition at line 314 of file mesh_function.h.

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

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 175 of file function_base.h.

Referenced by clear(), and init().

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

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

Definition at line 180 of file function_base.h.

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

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

Definition at line 169 of file function_base.h.

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

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 343 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()().

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 349 of file mesh_function.h.

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

PointLocatorBase* libMesh::MeshFunction::_point_locator
protected
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 331 of file mesh_function.h.

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

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

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

Definition at line 320 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: