Loading [MathJax]/extensions/tex2jax.js
libMesh
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
libMesh::MeshFunction Class Reference

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

#include <mesh_function.h>

Inheritance diagram for libMesh::MeshFunction:
[legend]

Public Member Functions

 MeshFunction (const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, std::vector< unsigned int > vars, const FunctionBase< Number > *master=nullptr)
 Constructor for mesh based functions with vectors as return value. More...
 
 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. More...
 
 MeshFunction (const MeshFunction &mf)
 A regular copy constructor. More...
 
 MeshFunction (MeshFunction &&)=default
 Special functions. More...
 
MeshFunctionoperator= (const MeshFunction &)=delete
 
MeshFunctionoperator= (MeshFunction &&)=delete
 
 ~MeshFunction ()
 Destructor. More...
 
virtual void init () override
 Override the FunctionBase::init() member function. More...
 
void init (const Trees::BuildType point_locator_build_type)
 The actual initialization process. More...
 
virtual void clear () override
 Clears the function. More...
 
virtual std::unique_ptr< FunctionBase< Number > > clone () const override
 
virtual 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.)
 
virtual void operator() (const Point &p, const Real time, DenseVector< Number > &output) override
 Computes values at coordinate p and for time time, which defaults to zero, optionally restricting the point to the MeshFunction subdomain_ids. More...
 
void 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, which parameter overrides the internal subdomain_ids. More...
 
void 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. More...
 
void 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. More...
 
void gradient (const Point &p, const Real time, std::vector< Gradient > &output)
 Computes gradients at coordinate p and for time time, which defaults to zero, optionally restricting the point to the MeshFunction subdomain_ids. More...
 
void gradient (const Point &p, const Real time, std::vector< Gradient > &output, const std::set< subdomain_id_type > *subdomain_ids)
 Computes gradients at coordinate p and for time time, which defaults to zero, optionally restricting the point to the passed subdomain_ids, which parameter overrides the internal subdomain_ids. More...
 
void 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. More...
 
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)
 Similar to gradient, but with the difference that multiple values on faces are explicitly permitted. More...
 
void hessian (const Point &p, const Real time, std::vector< Tensor > &output)
 Computes gradients at coordinate p and for time time, which defaults to zero, optionally restricting the point to the MeshFunction subdomain_ids. More...
 
void hessian (const Point &p, const Real time, std::vector< Tensor > &output, const std::set< subdomain_id_type > *subdomain_ids)
 Computes gradients at coordinate p and for time time, which defaults to zero, optionally restricting the point to the passed subdomain_ids. More...
 
const PointLocatorBaseget_point_locator () const
 
PointLocatorBaseget_point_locator ()
 
void enable_out_of_mesh_mode (const DenseVector< Number > &value)
 Enables out-of-mesh mode. More...
 
void enable_out_of_mesh_mode (const Number &value)
 Enables out-of-mesh mode. More...
 
void disable_out_of_mesh_mode ()
 Disables out-of-mesh mode. More...
 
void set_point_locator_tolerance (Real tol)
 We may want to specify a tolerance for the PointLocator to use, since in some cases the point we want to evaluate at might be slightly outside the mesh (due to numerical rounding issues, for example). More...
 
void unset_point_locator_tolerance ()
 Turn off the user-specified PointLocator tolerance. More...
 
void set_subdomain_ids (const std::set< subdomain_id_type > *subdomain_ids)
 Choose a default list of subdomain ids to be searched for points. More...
 
void operator() (const Point &p, DenseVector< Number > &output)
 Evaluation function for time-independent vector-valued functions. More...
 
virtual Number component (unsigned int i, const Point &p, Real time=0.)
 
bool initialized () const
 
void set_is_time_dependent (bool is_time_dependent)
 Function to set whether this is a time-dependent function or not. More...
 
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
 Helper function to reduce code duplication. More...
 
std::set< const Elem * > find_elements (const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
 
const Elemcheck_found_elem (const Elem *element, const Point &p) const
 Helper function that is called by MeshFunction::find_element() and MeshFunction::find_elements() to ensure that Elems found by the PointLocator are actually "evaluable" on the processor where they are found. More...
 
void _gradient_on_elem (const Point &p, const Elem *element, std::vector< Gradient > &output)
 Helper function for finding a gradient as evaluated from a specific element. More...
 

Protected Attributes

const EquationSystems_eqn_systems
 The equation systems handler, from which the data are gathered. More...
 
const NumericVector< Number > & _vector
 A reference to the vector that holds the data that is to be interpolated. More...
 
const DofMap_dof_map
 Need access to the DofMap of the other system. More...
 
const std::vector< unsigned int_system_vars
 The indices of the variables within the other system for which data are to be gathered. More...
 
std::unique_ptr< PointLocatorBase_point_locator
 A point locator is needed to locate the points in the mesh. More...
 
std::unique_ptr< std::set< subdomain_id_type > > _subdomain_ids
 A default set of subdomain ids in which to search for points. More...
 
bool _out_of_mesh_mode
 true if out-of-mesh mode is enabled. More...
 
DenseVector< Number_out_of_mesh_value
 Value to return outside the mesh if out-of-mesh mode is enabled. More...
 
const FunctionBase_master
 Const pointer to our master, initialized to nullptr. More...
 
bool _initialized
 When init() was called so that everything is ready for calls to operator() (...), then this bool is true. More...
 
bool _is_time_dependent
 Cache whether or not this function is actually time-dependent. More...
 
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 54 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,
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 42 of file mesh_function.C.

46  :
47  FunctionBase<Number> (master),
48  ParallelObject (eqn_systems),
49  _eqn_systems (eqn_systems),
50  _vector (vec),
51  _dof_map (dof_map),
52  _system_vars (std::move(vars)),
53  _out_of_mesh_mode (false)
54 {
55 }
ParallelObject(const Parallel::Communicator &comm_in)
Constructor.
const EquationSystems & _eqn_systems
The equation systems handler, from which the data are gathered.
const NumericVector< Number > & _vector
A reference to the vector that holds the data that is to be interpolated.
bool _out_of_mesh_mode
true if out-of-mesh mode is enabled.
const DofMap & _dof_map
Need access to the DofMap of the other system.
const std::vector< unsigned int > _system_vars
The indices of the variables within the other system for which data are to be gathered.

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

63  :
64  FunctionBase<Number> (master),
65  ParallelObject (eqn_systems),
66  _eqn_systems (eqn_systems),
67  _vector (vec),
68  _dof_map (dof_map),
69  _system_vars (1,var),
70  _out_of_mesh_mode (false)
71 {
72 }
ParallelObject(const Parallel::Communicator &comm_in)
Constructor.
const EquationSystems & _eqn_systems
The equation systems handler, from which the data are gathered.
const NumericVector< Number > & _vector
A reference to the vector that holds the data that is to be interpolated.
bool _out_of_mesh_mode
true if out-of-mesh mode is enabled.
const DofMap & _dof_map
Need access to the DofMap of the other system.
const std::vector< unsigned int > _system_vars
The indices of the variables within the other system for which data are to be gathered.

◆ MeshFunction() [3/4]

libMesh::MeshFunction::MeshFunction ( const MeshFunction mf)

A regular copy constructor.

Definition at line 74 of file mesh_function.C.

References _subdomain_ids, libMesh::PointLocatorBase::get_close_to_point_tol(), get_point_locator(), init(), libMesh::PointLocatorBase::initialized(), libMesh::FunctionBase< Output >::initialized(), and set_point_locator_tolerance().

74  :
75  FunctionBase<Number> (mf._master),
76  ParallelObject (mf._eqn_systems),
77  _eqn_systems (mf._eqn_systems),
78  _vector (mf._vector),
79  _dof_map (mf._dof_map),
80  _system_vars (mf._system_vars),
81  _out_of_mesh_mode (mf._out_of_mesh_mode)
82 {
83  // Initialize the mf and set the point locator if the
84  // input mf had done so.
85  if(mf.initialized())
86  {
87  this->MeshFunction::init();
88 
89  if(mf.get_point_locator().initialized())
90  this->set_point_locator_tolerance(mf.get_point_locator().get_close_to_point_tol());
91 
92  }
93 
94  if (mf._subdomain_ids)
96  std::make_unique<std::set<subdomain_id_type>>
97  (*mf._subdomain_ids);
98 }
ParallelObject(const Parallel::Communicator &comm_in)
Constructor.
std::unique_ptr< std::set< subdomain_id_type > > _subdomain_ids
A default set of subdomain ids in which to search for points.
const EquationSystems & _eqn_systems
The equation systems handler, from which the data are gathered.
const NumericVector< Number > & _vector
A reference to the vector that holds the data that is to be interpolated.
bool _out_of_mesh_mode
true if out-of-mesh mode is enabled.
const DofMap & _dof_map
Need access to the DofMap of the other system.
virtual void init() override
Override the FunctionBase::init() member function.
void set_point_locator_tolerance(Real tol)
We may want to specify a tolerance for the PointLocator to use, since in some cases the point we want...
const std::vector< unsigned int > _system_vars
The indices of the variables within the other system for which data are to be gathered.

◆ MeshFunction() [4/4]

libMesh::MeshFunction::MeshFunction ( MeshFunction &&  )
default

Special functions.

  • This class contains a unique_ptr so it can't be default copy constructed.
  • This class contains const references so it can't be default copy/move assigned.
  • The destructor is defaulted out-of-line.

◆ ~MeshFunction()

libMesh::MeshFunction::~MeshFunction ( )
default

Destructor.

Member Function Documentation

◆ _gradient_on_elem()

void libMesh::MeshFunction::_gradient_on_elem ( const Point p,
const Elem element,
std::vector< Gradient > &  output 
)
protected

Helper function for finding a gradient as evaluated from a specific element.

Definition at line 457 of file mesh_function.C.

References _dof_map, _eqn_systems, _out_of_mesh_mode, _out_of_mesh_value, _system_vars, _vector, libMesh::TypeVector< T >::add_scaled(), libMesh::FEGenericBase< OutputType >::build(), libMesh::FEInterface::compute_data(), dim, libMesh::Elem::dim(), libMesh::DofMap::dof_indices(), libMesh::FEComputeData::dshape, libMesh::FEComputeData::enable_derivative(), libMesh::index_range(), libMesh::Elem::infinite(), libMesh::invalid_uint, libMesh::FEMap::inverse_map(), libMesh::libmesh_assert(), libMesh::FEComputeData::local_transform, libMesh::DenseVector< T >::size(), and libMesh::DofMap::variable_type().

Referenced by discontinuous_gradient(), and gradient().

460 {
461  libmesh_assert(element);
462 
463  // resize the output vector to the number of output values
464  // that the user told us
465  output.resize (this->_system_vars.size());
466 
467  const unsigned int dim = element->dim();
468 
469  // Get local coordinates to feed these into compute_data().
470  // Note that the fe_type can safely be used from the 0-variable,
471  // since the inverse mapping is the same for all FEFamilies
472  const Point mapped_point (FEMap::inverse_map (dim, element,
473  p));
474 
475  std::vector<Point> point_list (1, mapped_point);
476 
477  // loop over all vars
478  for (auto index : index_range(this->_system_vars))
479  {
480  // the data for this variable
481  const unsigned int var = _system_vars[index];
482 
483  if (var == libMesh::invalid_uint)
484  {
486  index < _out_of_mesh_value.size());
487  output[index] = Gradient(_out_of_mesh_value(index));
488  continue;
489  }
490 
491  const FEType & fe_type = this->_dof_map.variable_type(var);
492 
493  // where the solution values for the var-th variable are stored
494  std::vector<dof_id_type> dof_indices;
495  this->_dof_map.dof_indices (element, dof_indices, var);
496 
497  // interpolate the solution
498  Gradient grad(0.);
499 
500  // for performance-reasons, we use different algorithms now.
501  // TODO: Check that both give the same result for finite elements.
502  // Otherwive it is wrong...
503  if (!element->infinite())
504  {
505  std::unique_ptr<FEBase> point_fe (FEBase::build(dim, fe_type));
506  const std::vector<std::vector<RealGradient>> & dphi = point_fe->get_dphi();
507  point_fe->reinit(element, &point_list);
508 
509  for (auto i : index_range(dof_indices))
510  grad.add_scaled(dphi[i][0], this->_vector(dof_indices[i]));
511 
512  }
513 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
514  else
515  {
516  // Build an FEComputeData that contains both input and output data
517  // for the specific compute_data method.
518  //
519  // TODO: enable this for a vector of points as well...
520  FEComputeData data (this->_eqn_systems, mapped_point);
521  data.enable_derivative();
522  FEInterface::compute_data (dim, fe_type, element, data);
523 
524  // grad [x] = data.dshape[i](v) * dv/dx * dof_index [i]
525  // sum over all indices
526  for (auto i : index_range(dof_indices))
527  {
528  // local coordinates
529  for (std::size_t v=0; v<dim; v++)
530  for (std::size_t xyz=0; xyz<LIBMESH_DIM; xyz++)
531  {
532  // FIXME: this needs better syntax: It is matrix-vector multiplication.
533  grad(xyz) += data.local_transform[v][xyz]
534  * data.dshape[i](v)
535  * this->_vector(dof_indices[i]);
536  }
537  }
538  }
539 #endif
540  output[index] = grad;
541  }
542 }
const EquationSystems & _eqn_systems
The equation systems handler, from which the data are gathered.
const NumericVector< Number > & _vector
A reference to the vector that holds the data that is to be interpolated.
bool _out_of_mesh_mode
true if out-of-mesh mode is enabled.
const unsigned int invalid_uint
A number which is used quite often to represent an invalid or uninitialized value for an unsigned int...
Definition: libmesh.h:293
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Fills the vector di with the global degree of freedom indices for the element.
Definition: dof_map.C:2184
unsigned int dim
static Point inverse_map(const unsigned int dim, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true, const bool extra_checks=true)
Definition: fe_map.C:1628
const DofMap & _dof_map
Need access to the DofMap of the other system.
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:2314
DenseVector< Number > _out_of_mesh_value
Value to return outside the mesh if out-of-mesh mode is enabled.
static void compute_data(const unsigned int dim, const FEType &fe_t, const Elem *elem, FEComputeData &data)
Lets the appropriate child of FEBase compute the requested data for the input specified in data...
NumberVectorValue Gradient
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
libmesh_assert(ctx)
virtual unsigned int size() const override final
Definition: dense_vector.h:104
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:117
const std::vector< unsigned int > _system_vars
The indices of the variables within the other system for which data are to be gathered.

◆ check_found_elem()

const Elem * libMesh::MeshFunction::check_found_elem ( const Elem element,
const Point p 
) const
protected

Helper function that is called by MeshFunction::find_element() and MeshFunction::find_elements() to ensure that Elems found by the PointLocator are actually "evaluable" on the processor where they are found.

Definition at line 689 of file mesh_function.C.

References _dof_map, _vector, libMesh::Elem::find_point_neighbors(), libMesh::GHOSTED, libMesh::DofMap::is_evaluable(), libMesh::ParallelObject::processor_id(), libMesh::DofObject::processor_id(), libMesh::SERIAL, and libMesh::NumericVector< T >::type().

Referenced by find_element(), and find_elements().

690 {
691  // If the PointLocator did not find an Element containing Point p,
692  // OR if the PointLocator found a local element containting Point p,
693  // we can simply return it.
694  if (!element || (element->processor_id() == this->processor_id()))
695  return element;
696 
697  // Otherwise, the PointLocator returned a valid but non-local
698  // element. Therefore, we have to do some further checks:
699  //
700  // 1.) If we have a SERIAL _vector, then we can just return the
701  // non-local element, because all the DOFs are available.
702  if (_vector.type() == SERIAL)
703  return element;
704 
705  // 2.) If we have a GHOSTED _vector, then we can return a non-local
706  // Element provided that all of its DOFs are ghosted on this
707  // processor. That should be faster than the other option, which
708  // is to search through all the Elem's point_neighbors() for a
709  // suitable local Elem.
710  if (_vector.type() == GHOSTED && _dof_map.is_evaluable(*element))
711  return element;
712 
713  // 3.) If we have a PARALLEL vector, then we search the non-local
714  // Elem's point neighbors for a local one to return instead. If
715  // we don't eventually find one, we just return nullptr.
716  std::set<const Elem *> point_neighbors;
717  element->find_point_neighbors(p, point_neighbors);
718  element = nullptr;
719  for (const auto & elem : point_neighbors)
720  if (elem->processor_id() == this->processor_id())
721  {
722  element = elem;
723  break;
724  }
725 
726  return element;
727 }
const NumericVector< Number > & _vector
A reference to the vector that holds the data that is to be interpolated.
const DofMap & _dof_map
Need access to the DofMap of the other system.
ParallelType type() const
processor_id_type processor_id() const
bool is_evaluable(const DofObjectSubclass &obj, unsigned int var_num=libMesh::invalid_uint) const
Definition: dof_map.C:2617

◆ clear()

void libMesh::MeshFunction::clear ( )
overridevirtual

Clears the function.

Reimplemented from libMesh::FunctionBase< Number >.

Definition at line 142 of file mesh_function.C.

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

143 {
144  // only delete the point locator when we are the master
145  if (_point_locator && !_master)
146  _point_locator.reset();
147 
148  this->_initialized = false;
149 }
const FunctionBase * _master
Const pointer to our master, initialized to nullptr.
bool _initialized
When init() was called so that everything is ready for calls to operator() (...), then this bool is t...
std::unique_ptr< PointLocatorBase > _point_locator
A point locator is needed to locate the points in the mesh.

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

154 {
155  return std::make_unique<MeshFunction>(*this);
156 }

◆ comm()

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

Definition at line 97 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::Partitioner::_find_global_index_by_pid_map(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::PetscLinearSolver< Number >::_petsc_shell_matrix_get_diagonal(), libMesh::SlepcEigenSolver< libMesh::Number >::_petsc_shell_matrix_get_diagonal(), libMesh::PetscLinearSolver< Number >::_petsc_shell_matrix_mult(), libMesh::SlepcEigenSolver< libMesh::Number >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< Number >::_petsc_shell_matrix_mult_add(), libMesh::MeshRefinement::_refine_elements(), libMesh::MeshRefinement::_smooth_flags(), libMesh::DofMap::add_constraints_to_send_list(), add_cube_convex_hull_to_mesh(), libMesh::PetscDMWrapper::add_dofs_helper(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::TransientRBConstruction::add_IC_to_RB_space(), libMesh::RBEIMEvaluation::add_interpolation_data(), libMesh::CondensedEigenSystem::add_matrices(), libMesh::EigenSystem::add_matrices(), libMesh::System::add_matrix(), libMesh::RBConstruction::add_scaled_matrix_and_vector(), libMesh::System::add_variable(), libMesh::System::add_variables(), libMesh::System::add_vector(), libMesh::MeshTools::Modification::all_tri(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::DofMap::allgather_recursive_constraints(), libMesh::TransientRBConstruction::allocate_data_structures(), libMesh::RBConstruction::allocate_data_structures(), libMesh::TransientRBConstruction::assemble_affine_expansion(), libMesh::AdvectionSystem::assemble_claw_rhs(), libMesh::FEMSystem::assemble_qoi(), libMesh::Nemesis_IO::assert_symmetric_cmaps(), libMesh::MeshCommunication::assign_global_indices(), libMesh::Partitioner::assign_partitioning(), libMesh::MeshTools::Generation::build_extrusion(), libMesh::Partitioner::build_graph(), 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_data(), libMesh::System::calculate_norm(), libMesh::DofMap::check_dirichlet_bcid_consistency(), libMesh::RBConstruction::compute_Fq_representor_innerprods(), libMesh::RBConstruction::compute_max_error_bound(), 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::RBConstruction::compute_output_dual_innerprods(), libMesh::RBConstruction::compute_residual_dual_norm_slow(), libMesh::RBSCMConstruction::compute_SCM_bounds_on_training_set(), libMesh::DofMap::computed_sparsity_already(), libMesh::Problem_Interface::computeF(), libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), libMesh::ContinuationSystem::ContinuationSystem(), libMesh::MeshBase::copy_constraint_rows(), libMesh::ExodusII_IO::copy_elemental_solution(), libMesh::ExodusII_IO::copy_nodal_solution(), libMesh::ExodusII_IO::copy_scalar_solution(), libMesh::CondensedEigenSystem::copy_super_to_sub(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::MeshTools::create_bounding_box(), libMesh::DofMap::create_dof_constraints(), 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::PetscMatrix< T >::create_submatrix_nosort(), create_wrapped_function(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::RBEIMEvaluation::distribute_bfs(), libMesh::DofMap::distribute_dofs(), DMlibMeshFunction(), DMlibMeshJacobian(), DMlibMeshSetSystem_libMesh(), DMVariableBounds_libMesh(), libMesh::DTKSolutionTransfer::DTKSolutionTransfer(), libMesh::MeshRefinement::eliminate_unrefined_patches(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_interiors(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_nodes(), libMesh::RBEIMConstruction::enrich_eim_approximation_on_sides(), libMesh::TransientRBConstruction::enrich_RB_space(), libMesh::EpetraVector< T >::EpetraVector(), AssembleOptimization::equality_constraints(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::WeightedPatchRecoveryErrorEstimator::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_error_tolerance(), libMesh::MeshRefinement::flag_elements_by_mean_stddev(), libMesh::MeshRefinement::flag_elements_by_nelem_target(), libMesh::RBEIMEvaluation::gather_bfs(), libMesh::DofMap::gather_constraints(), libMesh::MeshfreeInterpolation::gather_remote_data(), libMesh::CondensedEigenSystem::get_eigenpair(), libMesh::RBEIMEvaluation::get_eim_basis_function_node_value(), libMesh::RBEIMEvaluation::get_eim_basis_function_side_value(), libMesh::RBEIMEvaluation::get_eim_basis_function_value(), libMesh::MeshBase::get_info(), libMesh::System::get_info(), libMesh::DofMap::get_info(), libMesh::RBEIMEvaluation::get_interior_basis_functions_as_vecs(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::RBEIMConstruction::get_max_abs_value(), libMesh::RBEIMConstruction::get_node_max_abs_value(), libMesh::RBEIMEvaluation::get_parametrized_function_node_value(), libMesh::RBEIMEvaluation::get_parametrized_function_side_value(), libMesh::RBEIMEvaluation::get_parametrized_function_value(), libMesh::RBEIMConstruction::get_random_point(), AssembleOptimization::inequality_constraints(), AssembleOptimization::inequality_constraints_jacobian(), libMesh::LocationMap< T >::init(), libMesh::TimeSolver::init(), libMesh::StaticCondensation::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::PetscDMWrapper::init_and_attach_petscdm(), libMesh::AdvectionSystem::init_data(), libMesh::ClawSystem::init_data(), libMesh::PetscDMWrapper::init_petscdm(), libMesh::ExodusII_IO_Helper::initialize(), libMesh::OptimizationSystem::initialize_equality_constraints_storage(), libMesh::OptimizationSystem::initialize_inequality_constraints_storage(), libMesh::RBEIMConstruction::initialize_parametrized_functions_in_training_set(), libMesh::RBEIMConstruction::inner_product(), integrate_function(), libMesh::MeshTools::libmesh_assert_consistent_distributed(), libMesh::MeshTools::libmesh_assert_consistent_distributed_nodes(), libMesh::MeshTools::libmesh_assert_contiguous_dof_ids(), libMesh::MeshTools::libmesh_assert_equal_connectivity(), libMesh::MeshTools::libmesh_assert_equal_points(), 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_constraint_rows(), 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_linesearch_shellfunc(), libMesh::libmesh_petsc_preconditioner_apply(), libMesh::libmesh_petsc_recalculate_monitor(), libMesh::libmesh_petsc_snes_fd_residual(), libMesh::libmesh_petsc_snes_jacobian(), libMesh::libmesh_petsc_snes_mffd_interface(), libMesh::libmesh_petsc_snes_mffd_residual(), libMesh::libmesh_petsc_snes_postcheck(), libMesh::libmesh_petsc_snes_precheck(), 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::LinearImplicitSystem::LinearImplicitSystem(), main(), 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_bcids_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::TransientRBConstruction::mass_matrix_scaled_matvec(), libMesh::FEMSystem::mesh_position_set(), libMesh::TriangulatorInterface::MeshedHole::MeshedHole(), LinearElasticityWithContact::move_mesh(), libMesh::DistributedMesh::n_active_elem(), libMesh::MeshTools::n_active_levels(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::MeshTools::n_connected_components(), libMesh::DofMap::n_constrained_dofs(), libMesh::MeshBase::n_constraint_rows(), libMesh::DofMap::n_dofs(), libMesh::DofMap::n_dofs_per_processor(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::CondensedEigenSystem::n_global_non_condensed_dofs(), libMesh::MeshTools::n_levels(), MixedOrderTest::n_neighbor_links(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::SparsityPattern::Build::n_nonzeros(), libMesh::MeshTools::n_p_levels(), libMesh::BoundaryInfo::n_shellface_conds(), libMesh::RBEIMEvaluation::node_distribute_bfs(), libMesh::RBEIMEvaluation::node_gather_bfs(), libMesh::RBEIMConstruction::node_inner_product(), libMesh::MeshBase::operator==(), 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::BoundaryInfo::parallel_sync_node_ids(), libMesh::BoundaryInfo::parallel_sync_side_ids(), libMesh::MeshTools::paranoid_n_levels(), libMesh::Partitioner::partition(), 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::MeshBase::print_constraint_rows(), libMesh::DofMap::print_dof_constraints(), libMesh::DofMap::process_mesh_constraint_rows(), libMesh::Partitioner::processor_pairs_to_interface_nodes(), libMesh::InterMeshProjection::project_system_vectors(), FEMParameters::read(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::EquationSystems::read(), libMesh::ExodusII_IO::read_header(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::System::read_header(), libMesh::RBEIMEvaluation::read_in_interior_basis_functions(), libMesh::RBEIMEvaluation::read_in_node_basis_functions(), libMesh::RBEIMEvaluation::read_in_side_basis_functions(), libMesh::RBEvaluation::read_in_vectors_from_multiple_files(), libMesh::System::read_legacy_data(), libMesh::TransientRBConstruction::read_riesz_representors_from_files(), libMesh::RBConstruction::read_riesz_representors_from_files(), 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::Nemesis_IO_Helper::read_var_names_impl(), libMesh::MeshBase::recalculate_n_partitions(), libMesh::MeshRefinement::refine_and_coarsen_elements(), libMesh::SimplexRefiner::refine_via_edges(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::DistributedMesh::renumber_nodes_and_elements(), LinearElasticityWithContact::residual_and_jacobian(), OverlappingAlgebraicGhostingTest::run_ghosting_test(), OverlappingCouplingGhostingTest::run_sparsity_pattern_test(), scale_mesh_and_plot(), libMesh::DofMap::scatter_constraints(), libMesh::CheckpointIO::select_split_config(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::send_and_insert_dof_values(), libMesh::TransientRBConstruction::set_error_temporal_data(), libMesh::Partitioner::set_interface_node_processor_ids_BFS(), libMesh::Partitioner::set_interface_node_processor_ids_linear(), libMesh::Partitioner::set_interface_node_processor_ids_petscpartitioner(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::Partitioner::set_parent_processor_ids(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::PetscDiffSolver::setup_petsc_data(), libMesh::RBEIMEvaluation::side_distribute_bfs(), libMesh::RBEIMEvaluation::side_gather_bfs(), libMesh::RBEIMConstruction::side_inner_product(), libMesh::Partitioner::single_partition(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::ClawSystem::solve_conservation_law(), libMesh::split_mesh(), libMesh::RBEIMConstruction::store_eim_solutions_for_training_set(), libMesh::MeshBase::subdomain_ids(), libMesh::BoundaryInfo::sync(), ConstraintOperatorTest::test1DCoarseningNewNodes(), ConstraintOperatorTest::test1DCoarseningOperator(), libMesh::MeshRefinement::test_level_one(), MeshfunctionDFEM::test_mesh_function_dfem(), MeshfunctionDFEM::test_mesh_function_dfem_grad(), MeshFunctionTest::test_p_level(), libMesh::MeshRefinement::test_unflagged(), DofMapTest::testBadElemFECombo(), SystemsTest::testBlockRestrictedVarNDofs(), BoundaryInfoTest::testBoundaryOnChildrenErrors(), ConstraintOperatorTest::testCoreform(), ConnectedComponentsTest::testEdge(), MeshInputTest::testExodusIGASidesets(), MeshTriangulationTest::testFoundCenters(), PointLocatorTest::testLocator(), BoundaryInfoTest::testMesh(), PointLocatorTest::testPlanar(), MeshTriangulationTest::testPoly2TriRefinementBase(), SystemsTest::testProjectCubeWithMeshFunction(), BoundaryInfoTest::testRenumber(), CheckpointIOTest::testSplitter(), MeshInputTest::testTetgenIO(), MeshTriangulationTest::testTriangulatorInterp(), MeshTriangulationTest::testTriangulatorMeshedHoles(), MeshTriangulationTest::testTriangulatorRoundHole(), libMesh::MeshTools::total_weight(), libMesh::RBConstruction::train_reduced_basis_with_POD(), libMesh::MeshFunctionSolutionTransfer::transfer(), libMesh::MeshfreeSolutionTransfer::transfer(), libMesh::Poly2TriTriangulator::triangulate(), libMesh::TransientRBConstruction::truth_assembly(), libMesh::RBConstruction::truth_assembly(), libMesh::MeshRefinement::uniformly_coarsen(), update_current_local_solution(), libMesh::TransientRBConstruction::update_RB_initial_condition_all_N(), libMesh::TransientRBConstruction::update_RB_system_matrices(), libMesh::RBConstruction::update_RB_system_matrices(), libMesh::TransientRBConstruction::update_residual_terms(), libMesh::RBConstruction::update_residual_terms(), libMesh::MeshTools::volume(), libMesh::STLIO::write(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), libMesh::VTKIO::write_nodal_data(), libMesh::RBEIMEvaluation::write_out_interior_basis_functions(), libMesh::RBEIMEvaluation::write_out_node_basis_functions(), libMesh::RBEIMEvaluation::write_out_side_basis_functions(), libMesh::RBEvaluation::write_out_vectors(), libMesh::TransientRBConstruction::write_riesz_representors_to_files(), libMesh::RBConstruction::write_riesz_representors_to_files(), 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(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::RBDataSerialization::RBEvaluationSerialization::write_to_file(), libMesh::RBDataSerialization::TransientRBEvaluationSerialization::write_to_file(), libMesh::RBDataSerialization::RBEIMEvaluationSerialization::write_to_file(), and libMesh::RBDataSerialization::RBSCMEvaluationSerialization::write_to_file().

98  { 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.
The default implementation calls operator() with a DenseVector of size i+1 which will result in unexpected behaviour if operator() makes any access beyond that limit.

Reimplemented in TripleFunction, SolutionFunction< dim >, PeriodicQuadFunction, SolutionFunction< dim >, SolutionFunction< dim >, SolutionFunction< dim >, SolutionFunction< dim >, SolutionFunction< dim >, and SolutionFunction< dim >.

Definition at line 232 of file function_base.h.

235 {
236  DenseVector<Output> outvec(i+1);
237  (*this)(p, time, outvec);
238  return outvec(i);
239 }

◆ disable_out_of_mesh_mode()

void libMesh::MeshFunction::disable_out_of_mesh_mode ( )

Disables out-of-mesh mode.

This is also the default.

Definition at line 756 of file mesh_function.C.

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

757 {
758  libmesh_assert (this->initialized());
759  _point_locator->disable_out_of_mesh_mode();
760  _out_of_mesh_mode = false;
761 }
bool _out_of_mesh_mode
true if out-of-mesh mode is enabled.
libmesh_assert(ctx)
std::unique_ptr< PointLocatorBase > _point_locator
A point locator is needed to locate the points in the mesh.

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

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

Referenced by discontinuous_gradient().

197 {
198  libmesh_assert (this->initialized());
199 
200  std::map<const Elem *, std::vector<Gradient>> buffer;
201  this->discontinuous_gradient (p, time, buffer);
202  std::map<const Elem *, Gradient> return_value;
203  for (const auto & [elem, vec] : buffer)
204  return_value[elem] = vec[0];
205  // NOTE: If no suitable element is found, then the map return_value is empty. This
206  // puts burden on the user of this function but I don't really see a better way.
207  return return_value;
208 }
libmesh_assert(ctx)
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 419 of file mesh_function.C.

References _subdomain_ids, and discontinuous_gradient().

422 {
423  this->discontinuous_gradient (p, time, output, this->_subdomain_ids.get());
424 }
std::unique_ptr< std::set< subdomain_id_type > > _subdomain_ids
A default set of subdomain ids in which to search for points.
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 428 of file mesh_function.C.

References _gradient_on_elem(), _system_vars, find_elements(), libMesh::FunctionBase< Number >::initialized(), and libMesh::libmesh_assert().

432 {
433  libmesh_assert (this->initialized());
434 
435  // clear the output map
436  output.clear();
437 
438  // get the candidate elements
439  std::set<const Elem *> candidate_element = this->find_elements(p,subdomain_ids);
440 
441  // loop through all candidates, if the set is empty this function will return an
442  // empty map
443  for (const auto & element : candidate_element)
444  {
445  // define a temporary vector to store all values
446  std::vector<Gradient> temp_output (cast_int<unsigned int>(this->_system_vars.size()));
447 
448  this->_gradient_on_elem(p, element, temp_output);
449 
450  // Insert temp_output into output
451  output.emplace(element, std::move(temp_output));
452  }
453 }
void _gradient_on_elem(const Point &p, const Elem *element, std::vector< Gradient > &output)
Helper function for finding a gradient as evaluated from a specific element.
std::set< const Elem * > find_elements(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
libmesh_assert(ctx)
const std::vector< unsigned int > _system_vars
The indices of the variables within the other system for which data are to be gathered.

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

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

Referenced by discontinuous_value().

172 {
173  libmesh_assert (this->initialized());
174 
175  std::map<const Elem *, DenseVector<Number>> buffer;
176  this->discontinuous_value (p, time, buffer);
177  std::map<const Elem *, Number> return_value;
178  for (const auto & [elem, vec] : buffer)
179  return_value[elem] = vec(0);
180  // NOTE: If no suitable element is found, then the map return_value is empty. This
181  // puts burden on the user of this function but I don't really see a better way.
182  return return_value;
183 }
std::map< const Elem *, Number > discontinuous_value(const Point &p, const Real time=0.)
libmesh_assert(ctx)

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

References _subdomain_ids, and discontinuous_value().

312 {
313  this->discontinuous_value (p, time, output, this->_subdomain_ids.get());
314 }
std::unique_ptr< std::set< subdomain_id_type > > _subdomain_ids
A default set of subdomain ids in which to search for points.
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.

Definition at line 318 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(), dim, libMesh::DofMap::dof_indices(), find_elements(), libMesh::index_range(), libMesh::FunctionBase< Number >::initialized(), libMesh::invalid_uint, libMesh::FEMap::inverse_map(), libMesh::libmesh_assert(), libMesh::FEComputeData::shape, libMesh::DenseVector< T >::size(), value, and libMesh::DofMap::variable_type().

322 {
323  libmesh_assert (this->initialized());
324 
325  // clear the output map
326  output.clear();
327 
328  // get the candidate elements
329  std::set<const Elem *> candidate_element = this->find_elements(p,subdomain_ids);
330 
331  // loop through all candidates, if the set is empty this function will return an
332  // empty map
333  for (const auto & element : candidate_element)
334  {
335  const unsigned int dim = element->dim();
336 
337  // define a temporary vector to store all values
338  DenseVector<Number> temp_output (cast_int<unsigned int>(this->_system_vars.size()));
339 
340  // Get local coordinates to feed these into compute_data().
341  // Note that the fe_type can safely be used from the 0-variable,
342  // since the inverse mapping is the same for all FEFamilies
343  const Point mapped_point (FEMap::inverse_map (dim, element, p));
344 
345  // loop over all vars
346  for (auto index : index_range(this->_system_vars))
347  {
348  // the data for this variable
349  const unsigned int var = _system_vars[index];
350 
351  if (var == libMesh::invalid_uint)
352  {
354  index < _out_of_mesh_value.size());
355  temp_output(index) = _out_of_mesh_value(index);
356  continue;
357  }
358 
359  const FEType & fe_type = this->_dof_map.variable_type(var);
360 
361  // Build an FEComputeData that contains both input and output data
362  // for the specific compute_data method.
363  {
364  FEComputeData data (this->_eqn_systems, mapped_point);
365 
366  FEInterface::compute_data (dim, fe_type, element, data);
367 
368  // where the solution values for the var-th variable are stored
369  std::vector<dof_id_type> dof_indices;
370  this->_dof_map.dof_indices (element, dof_indices, var);
371 
372  // interpolate the solution
373  {
374  Number value = 0.;
375 
376  for (auto i : index_range(dof_indices))
377  value += this->_vector(dof_indices[i]) * data.shape[i];
378 
379  temp_output(index) = value;
380  }
381 
382  }
383 
384  // next variable
385  }
386 
387  // Insert temp_output into output
388  output[element] = temp_output;
389  }
390 }
const EquationSystems & _eqn_systems
The equation systems handler, from which the data are gathered.
const NumericVector< Number > & _vector
A reference to the vector that holds the data that is to be interpolated.
bool _out_of_mesh_mode
true if out-of-mesh mode is enabled.
const unsigned int invalid_uint
A number which is used quite often to represent an invalid or uninitialized value for an unsigned int...
Definition: libmesh.h:293
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Fills the vector di with the global degree of freedom indices for the element.
Definition: dof_map.C:2184
unsigned int dim
static Point inverse_map(const unsigned int dim, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true, const bool extra_checks=true)
Definition: fe_map.C:1628
const DofMap & _dof_map
Need access to the DofMap of the other system.
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:2314
DenseVector< Number > _out_of_mesh_value
Value to return outside the mesh if out-of-mesh mode is enabled.
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)
Lets the appropriate child of FEBase compute the requested data for the input specified in data...
libmesh_assert(ctx)
static const bool value
Definition: xdr_io.C:54
virtual unsigned int size() const override final
Definition: dense_vector.h:104
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:117
const std::vector< unsigned int > _system_vars
The indices of the variables within the other system for which data are to be gathered.

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

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

Referenced by enable_out_of_mesh_mode().

742 {
743  libmesh_assert (this->initialized());
744  _point_locator->enable_out_of_mesh_mode();
745  _out_of_mesh_mode = true;
747 }
bool _out_of_mesh_mode
true if out-of-mesh mode is enabled.
DenseVector< Number > _out_of_mesh_value
Value to return outside the mesh if out-of-mesh mode is enabled.
libmesh_assert(ctx)
static const bool value
Definition: xdr_io.C:54
std::unique_ptr< PointLocatorBase > _point_locator
A point locator is needed to locate the points in the mesh.

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

References enable_out_of_mesh_mode(), and value.

750 {
751  DenseVector<Number> v(1);
752  v(0) = value;
753  this->enable_out_of_mesh_mode(v);
754 }
void enable_out_of_mesh_mode(const DenseVector< Number > &value)
Enables out-of-mesh mode.
static const bool value
Definition: xdr_io.C:54

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

References libMesh::FunctionBase< Number >::_master, _out_of_mesh_mode, and check_found_elem().

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

632 {
633  // Ensure that in the case of a master mesh function, the
634  // out-of-mesh mode is enabled either for both or for none. This is
635  // important because the out-of-mesh mode is also communicated to
636  // the point locator. Since this is time consuming, enable it only
637  // in debug mode.
638 #ifdef DEBUG
639  if (this->_master != nullptr)
640  {
641  const MeshFunction * master =
642  cast_ptr<const MeshFunction *>(this->_master);
643  libmesh_error_msg_if(_out_of_mesh_mode!=master->_out_of_mesh_mode,
644  "ERROR: If you use out-of-mesh-mode in connection with master mesh "
645  "functions, you must enable out-of-mesh mode for both the master and the slave mesh function.");
646  }
647 #endif
648 
649  // locate the point in the other mesh
650  const Elem * element = (*_point_locator)(p, subdomain_ids);
651 
652  // Make sure that the element found is evaluable
653  return this->check_found_elem(element, p);
654 }
bool _out_of_mesh_mode
true if out-of-mesh mode is enabled.
const FunctionBase * _master
Const pointer to our master, initialized to nullptr.
MeshFunction(const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, std::vector< unsigned int > vars, const FunctionBase< Number > *master=nullptr)
Constructor for mesh based functions with vectors as return value.
Definition: mesh_function.C:42
const Elem * check_found_elem(const Elem *element, const Point &p) const
Helper function that is called by MeshFunction::find_element() and MeshFunction::find_elements() to e...

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

References libMesh::FunctionBase< Number >::_master, _out_of_mesh_mode, and check_found_elem().

Referenced by discontinuous_gradient(), and discontinuous_value().

658 {
659  // Ensure that in the case of a master mesh function, the
660  // out-of-mesh mode is enabled either for both or for none. This is
661  // important because the out-of-mesh mode is also communicated to
662  // the point locator. Since this is time consuming, enable it only
663  // in debug mode.
664 #ifdef DEBUG
665  if (this->_master != nullptr)
666  {
667  const MeshFunction * master =
668  cast_ptr<const MeshFunction *>(this->_master);
669  libmesh_error_msg_if(_out_of_mesh_mode!=master->_out_of_mesh_mode,
670  "ERROR: If you use out-of-mesh-mode in connection with master mesh "
671  "functions, you must enable out-of-mesh mode for both the master and the slave mesh function.");
672  }
673 #endif
674 
675  // locate the point in the other mesh
676  std::set<const Elem *> candidate_elements;
677  std::set<const Elem *> final_candidate_elements;
678  (*_point_locator)(p, candidate_elements, subdomain_ids);
679 
680  // For each candidate Elem, if it is evaluable, add it to the set of
681  // final candidate Elems.
682  for (const auto & element : candidate_elements)
683  final_candidate_elements.insert(this->check_found_elem(element, p));
684 
685  return final_candidate_elements;
686 }
bool _out_of_mesh_mode
true if out-of-mesh mode is enabled.
const FunctionBase * _master
Const pointer to our master, initialized to nullptr.
MeshFunction(const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, std::vector< unsigned int > vars, const FunctionBase< Number > *master=nullptr)
Constructor for mesh based functions with vectors as return value.
Definition: mesh_function.C:42
const Elem * check_found_elem(const Elem *element, const Point &p) const
Helper function that is called by MeshFunction::find_element() and MeshFunction::find_elements() to e...

◆ get_point_locator() [1/2]

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

Definition at line 729 of file mesh_function.C.

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

Referenced by MeshFunction().

730 {
731  libmesh_assert (this->initialized());
732  return *_point_locator;
733 }
libmesh_assert(ctx)
std::unique_ptr< PointLocatorBase > _point_locator
A point locator is needed to locate the points in the mesh.

◆ get_point_locator() [2/2]

PointLocatorBase & libMesh::MeshFunction::get_point_locator ( )

Definition at line 735 of file mesh_function.C.

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

736 {
737  libmesh_assert (this->initialized());
738  return *_point_locator;
739 }
libmesh_assert(ctx)
std::unique_ptr< PointLocatorBase > _point_locator
A point locator is needed to locate the points in the mesh.

◆ gradient() [1/3]

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

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

Referenced by libMesh::ExactErrorEstimator::find_squared_element_error(), gptr(), and gradient().

187 {
188  libmesh_assert (this->initialized());
189 
190  std::vector<Gradient> buf (1);
191  this->gradient(p, time, buf);
192  return buf[0];
193 }
Gradient gradient(const Point &p, const Real time=0.)
libmesh_assert(ctx)

◆ gradient() [2/3]

void libMesh::MeshFunction::gradient ( const Point p,
const Real  time,
std::vector< Gradient > &  output 
)

Computes gradients at coordinate p and for time time, which defaults to zero, optionally restricting the point to the MeshFunction subdomain_ids.

Definition at line 394 of file mesh_function.C.

References _subdomain_ids, and gradient().

397 {
398  this->gradient(p, time, output, this->_subdomain_ids.get());
399 }
std::unique_ptr< std::set< subdomain_id_type > > _subdomain_ids
A default set of subdomain ids in which to search for points.
Gradient gradient(const Point &p, const Real time=0.)

◆ gradient() [3/3]

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

Computes gradients at coordinate p and for time time, which defaults to zero, optionally restricting the point to the passed subdomain_ids, which parameter overrides the internal subdomain_ids.

Definition at line 403 of file mesh_function.C.

References _gradient_on_elem(), find_element(), libMesh::FunctionBase< Number >::initialized(), and libMesh::libmesh_assert().

407 {
408  libmesh_assert (this->initialized());
409 
410  const Elem * element = this->find_element(p,subdomain_ids);
411 
412  if (!element)
413  output.resize(0);
414  else
415  this->_gradient_on_elem(p, element, output);
416 }
const Elem * find_element(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
Helper function to reduce code duplication.
void _gradient_on_elem(const Point &p, const Elem *element, std::vector< Gradient > &output)
Helper function for finding a gradient as evaluated from a specific element.
libmesh_assert(ctx)

◆ hessian() [1/3]

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

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

Referenced by libMesh::ExactErrorEstimator::find_squared_element_error(), and hessian().

213 {
214  libmesh_assert (this->initialized());
215 
216  std::vector<Tensor> buf (1);
217  this->hessian(p, time, buf);
218  return buf[0];
219 }
libmesh_assert(ctx)
Tensor hessian(const Point &p, const Real time=0.)

◆ hessian() [2/3]

void libMesh::MeshFunction::hessian ( const Point p,
const Real  time,
std::vector< Tensor > &  output 
)

Computes gradients at coordinate p and for time time, which defaults to zero, optionally restricting the point to the MeshFunction subdomain_ids.

Definition at line 546 of file mesh_function.C.

References _subdomain_ids, and hessian().

549 {
550  this->hessian(p, time, output, this->_subdomain_ids.get());
551 }
std::unique_ptr< std::set< subdomain_id_type > > _subdomain_ids
A default set of subdomain ids in which to search for points.
Tensor hessian(const Point &p, const Real time=0.)

◆ hessian() [3/3]

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

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

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

559 {
560  libmesh_assert (this->initialized());
561 
562  const Elem * element = this->find_element(p,subdomain_ids);
563 
564  if (!element)
565  {
566  output.resize(0);
567  }
568  else
569  {
570  if(element->infinite())
571  libmesh_warning("Warning: Requested the Hessian of an Infinite element."
572  << "Second derivatives for Infinite elements"
573  << " are not yet implemented!"
574  << std::endl);
575 
576  // resize the output vector to the number of output values
577  // that the user told us
578  output.resize (this->_system_vars.size());
579 
580 
581  {
582  const unsigned int dim = element->dim();
583 
584 
585  // Get local coordinates to feed these into compute_data().
586  // Note that the fe_type can safely be used from the 0-variable,
587  // since the inverse mapping is the same for all FEFamilies
588  const Point mapped_point (FEMap::inverse_map (dim, element,
589  p));
590 
591  std::vector<Point> point_list (1, mapped_point);
592 
593  // loop over all vars
594  for (auto index : index_range(this->_system_vars))
595  {
596  // the data for this variable
597  const unsigned int var = _system_vars[index];
598 
599  if (var == libMesh::invalid_uint)
600  {
602  index < _out_of_mesh_value.size());
603  output[index] = Tensor(_out_of_mesh_value(index));
604  continue;
605  }
606  const FEType & fe_type = this->_dof_map.variable_type(var);
607 
608  std::unique_ptr<FEBase> point_fe (FEBase::build(dim, fe_type));
609  const std::vector<std::vector<RealTensor>> & d2phi =
610  point_fe->get_d2phi();
611  point_fe->reinit(element, &point_list);
612 
613  // where the solution values for the var-th variable are stored
614  std::vector<dof_id_type> dof_indices;
615  this->_dof_map.dof_indices (element, dof_indices, var);
616 
617  // interpolate the solution
618  Tensor hess;
619 
620  for (auto i : index_range(dof_indices))
621  hess.add_scaled(d2phi[i][0], this->_vector(dof_indices[i]));
622 
623  output[index] = hess;
624  }
625  }
626  }
627 }
const NumericVector< Number > & _vector
A reference to the vector that holds the data that is to be interpolated.
bool _out_of_mesh_mode
true if out-of-mesh mode is enabled.
const unsigned int invalid_uint
A number which is used quite often to represent an invalid or uninitialized value for an unsigned int...
Definition: libmesh.h:293
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Fills the vector di with the global degree of freedom indices for the element.
Definition: dof_map.C:2184
unsigned int dim
static Point inverse_map(const unsigned int dim, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true, const bool extra_checks=true)
Definition: fe_map.C:1628
const DofMap & _dof_map
Need access to the DofMap of the other system.
const Elem * find_element(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
Helper function to reduce code duplication.
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:2314
DenseVector< Number > _out_of_mesh_value
Value to return outside the mesh if out-of-mesh mode is enabled.
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
libmesh_assert(ctx)
void add_scaled(const TypeTensor< T2 > &, const T &)
Add a scaled tensor to this tensor without creating a temporary.
Definition: type_tensor.h:851
NumberTensorValue Tensor
virtual unsigned int size() const override final
Definition: dense_vector.h:104
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:117
const std::vector< unsigned int > _system_vars
The indices of the variables within the other system for which data are to be gathered.

◆ init() [1/2]

void libMesh::MeshFunction::init ( )
overridevirtual

Override the FunctionBase::init() member function.

Reimplemented from libMesh::FunctionBase< Number >.

Definition at line 104 of file mesh_function.C.

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

Referenced by init(), main(), MeshFunction(), libMesh::InterMeshProjection::project_system_vectors(), MeshfunctionDFEM::test_mesh_function_dfem(), MeshfunctionDFEM::test_mesh_function_dfem_grad(), MeshFunctionTest::test_subdomain_id_sets(), SystemsTest::testProjectCubeWithMeshFunction(), and libMesh::MeshFunctionSolutionTransfer::transfer().

105 {
106  // are indices of the desired variable(s) provided?
107  libmesh_assert_greater (this->_system_vars.size(), 0);
108 
109  // Don't do twice...
110  if (this->_initialized)
111  {
113  return;
114  }
115 
116  // The Mesh owns the "master" PointLocator, while handing us a
117  // PointLocator "proxy" that forwards all requests to the master.
118  const MeshBase & mesh = this->_eqn_systems.get_mesh();
119  _point_locator = mesh.sub_point_locator();
120 
121  // ready for use
122  this->_initialized = true;
123 }
const EquationSystems & _eqn_systems
The equation systems handler, from which the data are gathered.
MeshBase & mesh
bool _initialized
When init() was called so that everything is ready for calls to operator() (...), then this bool is t...
libmesh_assert(ctx)
const MeshBase & get_mesh() const
std::unique_ptr< PointLocatorBase > _point_locator
A point locator is needed to locate the points in the mesh.
const std::vector< unsigned int > _system_vars
The indices of the variables within the other system for which data are to be gathered.

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

References init().

129 {
130  libmesh_deprecated();
131 
132  // Call the init() taking no args instead. Note: this is backwards
133  // compatible because the argument was not used for anything
134  // previously anyway.
135  this->init();
136 }
virtual void init() override
Override the FunctionBase::init() member function.

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

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

211 {
212  return (this->_initialized);
213 }
bool _initialized
When init() was called so that everything is ready for calls to operator() (...), then this bool is t...

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

225 {
226  return (this->_is_time_dependent);
227 }
bool _is_time_dependent
Cache whether or not this function is actually time-dependent.

◆ n_processors()

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

Definition at line 103 of file parallel_object.h.

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

Referenced by libMesh::Partitioner::_find_global_index_by_pid_map(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::DofMap::add_constraints_to_send_list(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::DistributedMesh::add_elem(), libMesh::DofMap::add_neighbors_to_send_list(), libMesh::DistributedMesh::add_node(), libMesh::System::add_vector(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::DofMap::allgather_recursive_constraints(), libMesh::FEMSystem::assembly(), libMesh::Nemesis_IO::assert_symmetric_cmaps(), libMesh::Partitioner::assign_partitioning(), libMesh::AztecLinearSolver< T >::AztecLinearSolver(), libMesh::Partitioner::build_graph(), libMesh::EquationSystems::build_parallel_elemental_solution_vector(), libMesh::DistributedMesh::clear(), libMesh::DistributedMesh::clear_elems(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::UnstructuredMesh::copy_nodes_and_elements(), libMesh::ExodusII_IO::copy_scalar_solution(), libMesh::Nemesis_IO::copy_scalar_solution(), libMesh::UnstructuredMesh::create_pid_mesh(), libMesh::MeshTools::create_processor_bounding_box(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_scalar_dofs(), libMesh::DistributedMesh::DistributedMesh(), libMesh::EnsightIO::EnsightIO(), libMesh::RBEIMEvaluation::gather_bfs(), libMesh::MeshBase::get_info(), libMesh::SystemSubsetBySubdomain::init(), libMesh::PetscDMWrapper::init_petscdm(), libMesh::Nemesis_IO_Helper::initialize(), libMesh::ExodusII_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::DofMap::n_dofs_per_processor(), libMesh::MeshBase::n_elem_on_proc(), libMesh::MeshBase::n_nodes_on_proc(), libMesh::RBEIMEvaluation::node_gather_bfs(), libMesh::Partitioner::partition(), libMesh::MeshBase::partition(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::DofMap::prepare_send_list(), libMesh::MeshBase::print_constraint_rows(), libMesh::DofMap::print_dof_constraints(), 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::Partitioner::repartition(), OverlappingFunctorTest::run_partitioner_test(), libMesh::DofMap::scatter_constraints(), libMesh::DistributedMesh::set_next_unique_id(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::PetscDMWrapper::set_point_range_in_section(), WriteVecAndScalar::setupTests(), libMesh::RBEIMEvaluation::side_gather_bfs(), DistributedMeshTest::testRemoteElemError(), CheckpointIOTest::testSplitter(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::GMVIO::write_binary(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::ExodusII_IO_Helper::write_nodal_coordinates(), libMesh::VTKIO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), 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().

104  {
105  processor_id_type returnval =
106  cast_int<processor_id_type>(_communicator.size());
107  libmesh_assert(returnval); // We never have an empty comm
108  return returnval;
109  }
const Parallel::Communicator & _communicator
processor_id_type size() const
uint8_t processor_id_type
libmesh_assert(ctx)

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

247 {
248  // Call the time-dependent function with t=0.
249  this->operator()(p, 0., output);
250 }
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 160 of file mesh_function.C.

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

Referenced by operator()().

162 {
163  libmesh_assert (this->initialized());
164 
165  DenseVector<Number> buf (1);
166  this->operator() (p, time, buf);
167  return buf(0);
168 }
libmesh_assert(ctx)
virtual 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 MeshFunction subdomain_ids.

This is useful in cases where there are multiple dimensioned elements, for example.

Implements libMesh::FunctionBase< Number >.

Definition at line 222 of file mesh_function.C.

References _subdomain_ids, and operator()().

225 {
226  this->operator() (p, time, output, this->_subdomain_ids.get());
227 }
std::unique_ptr< std::set< subdomain_id_type > > _subdomain_ids
A default set of subdomain ids in which to search for points.
virtual 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, which parameter overrides the internal subdomain_ids.

Definition at line 229 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(), dim, libMesh::Elem::dim(), libMesh::DofMap::dof_indices(), find_element(), libMesh::index_range(), libMesh::FunctionBase< Number >::initialized(), libMesh::invalid_uint, libMesh::FEMap::inverse_map(), libMesh::libmesh_assert(), libMesh::DenseVector< T >::resize(), libMesh::FEComputeData::shape, libMesh::DenseVector< T >::size(), value, and libMesh::DofMap::variable_type().

233 {
234  libmesh_assert (this->initialized());
235 
236  const Elem * element = this->find_element(p,subdomain_ids);
237 
238  if (!element)
239  {
240  // We'd better be in out_of_mesh_mode if we couldn't find an
241  // element in the mesh
243  output = _out_of_mesh_value;
244  }
245  else
246  {
247  // resize the output vector to the number of output values
248  // that the user told us
249  output.resize (cast_int<unsigned int>
250  (this->_system_vars.size()));
251 
252 
253  {
254  const unsigned int dim = element->dim();
255 
256 
257  // Get local coordinates to feed these into compute_data().
258  // Note that the fe_type can safely be used from the 0-variable,
259  // since the inverse mapping is the same for all FEFamilies
260  const Point mapped_point (FEMap::inverse_map (dim, element,
261  p));
262 
263  // loop over all vars
264  for (auto index : index_range(this->_system_vars))
265  {
266  // the data for this variable
267  const unsigned int var = _system_vars[index];
268 
269  if (var == libMesh::invalid_uint)
270  {
272  index < _out_of_mesh_value.size());
273  output(index) = _out_of_mesh_value(index);
274  continue;
275  }
276 
277  const FEType & fe_type = this->_dof_map.variable_type(var);
278 
279  // Build an FEComputeData that contains both input and output data
280  // for the specific compute_data method.
281  {
282  FEComputeData data (this->_eqn_systems, mapped_point);
283 
284  FEInterface::compute_data (dim, fe_type, element, data);
285 
286  // where the solution values for the var-th variable are stored
287  std::vector<dof_id_type> dof_indices;
288  this->_dof_map.dof_indices (element, dof_indices, var);
289 
290  // interpolate the solution
291  {
292  Number value = 0.;
293 
294  for (auto i : index_range(dof_indices))
295  value += this->_vector(dof_indices[i]) * data.shape[i];
296 
297  output(index) = value;
298  }
299 
300  }
301 
302  // next variable
303  }
304  }
305  }
306 }
const EquationSystems & _eqn_systems
The equation systems handler, from which the data are gathered.
const NumericVector< Number > & _vector
A reference to the vector that holds the data that is to be interpolated.
bool _out_of_mesh_mode
true if out-of-mesh mode is enabled.
const unsigned int invalid_uint
A number which is used quite often to represent an invalid or uninitialized value for an unsigned int...
Definition: libmesh.h:293
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
Fills the vector di with the global degree of freedom indices for the element.
Definition: dof_map.C:2184
unsigned int dim
static Point inverse_map(const unsigned int dim, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true, const bool extra_checks=true)
Definition: fe_map.C:1628
const DofMap & _dof_map
Need access to the DofMap of the other system.
void resize(const unsigned int n)
Resize the vector.
Definition: dense_vector.h:384
const Elem * find_element(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
Helper function to reduce code duplication.
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:2314
DenseVector< Number > _out_of_mesh_value
Value to return outside the mesh if out-of-mesh mode is enabled.
static void compute_data(const unsigned int dim, const FEType &fe_t, const Elem *elem, FEComputeData &data)
Lets the appropriate child of FEBase compute the requested data for the input specified in data...
libmesh_assert(ctx)
static const bool value
Definition: xdr_io.C:54
virtual unsigned int size() const override final
Definition: dense_vector.h:104
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:117
const std::vector< unsigned int > _system_vars
The indices of the variables within the other system for which data are to be gathered.

◆ operator=() [1/2]

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

◆ 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 114 of file parallel_object.h.

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

Referenced by libMesh::BoundaryInfo::_find_id_maps(), 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::MeshTools::Modification::all_tri(), libMesh::DofMap::allgather_recursive_constraints(), libMesh::FEMSystem::assembly(), libMesh::Nemesis_IO::assert_symmetric_cmaps(), libMesh::Partitioner::assign_partitioning(), libMesh::Nemesis_IO_Helper::build_element_and_node_maps(), libMesh::Partitioner::build_graph(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::EquationSystems::build_parallel_elemental_solution_vector(), libMesh::EquationSystems::build_parallel_solution_vector(), check_found_elem(), libMesh::DistributedMesh::clear(), libMesh::DistributedMesh::clear_elems(), 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::RBConstruction::compute_max_error_bound(), 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::ExodusII_IO::copy_elemental_solution(), libMesh::ExodusII_IO::copy_nodal_solution(), libMesh::ExodusII_IO::copy_scalar_solution(), libMesh::Nemesis_IO::copy_scalar_solution(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::ExodusII_IO_Helper::create(), libMesh::DistributedMesh::delete_elem(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_scalar_dofs(), libMesh::DistributedMesh::DistributedMesh(), libMesh::DofMap::end_dof(), libMesh::DofMap::end_old_dof(), libMesh::EnsightIO::EnsightIO(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::find_dofs_to_send(), libMesh::UnstructuredMesh::find_neighbors(), libMesh::DofMap::first_dof(), libMesh::DofMap::first_old_dof(), libMesh::RBEIMEvaluation::gather_bfs(), 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::RBEIMEvaluation::get_interior_basis_functions_as_vecs(), libMesh::Nemesis_IO_Helper::get_loadbal_param(), libMesh::DofMap::get_local_constraints(), libMesh::MeshBase::get_local_constraints(), 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::LaplaceMeshSmoother::init(), libMesh::StaticCondensation::init(), libMesh::SystemSubsetBySubdomain::init(), HeatSystem::init_data(), 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::TransientRBEvaluation::legacy_write_offline_data_to_files(), libMesh::RBSCMEvaluation::legacy_write_offline_data_to_files(), libMesh::RBEvaluation::legacy_write_offline_data_to_files(), 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(), main(), libMesh::MeshRefinement::make_coarsening_compatible(), AugmentSparsityOnInterface::mesh_reinit(), libMesh::TriangulatorInterface::MeshedHole::MeshedHole(), libMesh::MeshBase::n_active_local_elem(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::MeshTools::n_connected_components(), libMesh::MeshBase::n_constraint_rows(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::DofMap::n_local_dofs(), libMesh::MeshBase::n_local_elem(), libMesh::MeshBase::n_local_nodes(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::BoundaryInfo::n_shellface_conds(), libMesh::RBEIMEvaluation::node_gather_bfs(), libMesh::DistributedMesh::own_node(), libMesh::BoundaryInfo::parallel_sync_node_ids(), libMesh::BoundaryInfo::parallel_sync_side_ids(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::MeshBase::print_constraint_rows(), libMesh::DofMap::print_dof_constraints(), libMesh::DofMap::process_mesh_constraint_rows(), 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::EquationSystems::read(), libMesh::ExodusII_IO_Helper::read_elem_num_map(), libMesh::ExodusII_IO_Helper::read_global_values(), libMesh::ExodusII_IO::read_header(), libMesh::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), libMesh::System::read_header(), libMesh::System::read_legacy_data(), libMesh::DynaIO::read_mesh(), libMesh::ExodusII_IO_Helper::read_node_num_map(), libMesh::System::read_parallel_data(), libMesh::TransientRBConstruction::read_riesz_representors_from_files(), libMesh::RBConstruction::read_riesz_representors_from_files(), 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::Nemesis_IO_Helper::read_var_names_impl(), libMesh::SimplexRefiner::refine_via_edges(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::DistributedMesh::renumber_nodes_and_elements(), libMesh::DofMap::scatter_constraints(), libMesh::CheckpointIO::select_split_config(), libMesh::DistributedMesh::set_next_unique_id(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::RBEIMEvaluation::side_gather_bfs(), ExodusTest< elem_type >::test_read_gold(), ExodusTest< elem_type >::test_write(), MeshInputTest::testAbaqusRead(), MeshInputTest::testBadGmsh(), MeshInputTest::testCopyElementSolutionImpl(), MeshInputTest::testCopyElementVectorImpl(), MeshInputTest::testCopyNodalSolutionImpl(), DefaultCouplingTest::testCoupling(), PointNeighborCouplingTest::testCoupling(), MeshInputTest::testDynaFileMappings(), MeshInputTest::testDynaNoSplines(), MeshInputTest::testDynaReadElem(), MeshInputTest::testDynaReadPatch(), MeshInputTest::testExodusFileMappings(), MeshInputTest::testExodusIGASidesets(), MeshInputTest::testExodusWriteElementDataFromDiscontinuousNodalData(), MeshInputTest::testGoodGmsh(), MeshInputTest::testGoodSTL(), MeshInputTest::testGoodSTLBinary(), MeshInputTest::testLowOrderEdgeBlocks(), SystemsTest::testProjectMatrix3D(), BoundaryInfoTest::testShellFaceConstraints(), MeshInputTest::testSingleElementImpl(), WriteVecAndScalar::testSolution(), CheckpointIOTest::testSplitter(), MeshInputTest::testTetgenIO(), libMesh::MeshTools::total_weight(), libMesh::NetGenMeshInterface::triangulate(), libMesh::MeshRefinement::uniformly_coarsen(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::DTKAdapter::update_variable_values(), libMesh::MeshTools::volume(), libMesh::STLIO::write(), 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_element_values_element_major(), libMesh::ExodusII_IO_Helper::write_elements(), libMesh::ExodusII_IO_Helper::write_elemset_data(), libMesh::ExodusII_IO_Helper::write_elemsets(), 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::VTKIO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_common(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), libMesh::ExodusII_IO_Helper::write_nodal_values(), libMesh::ExodusII_IO_Helper::write_nodeset_data(), libMesh::Nemesis_IO_Helper::write_nodesets(), libMesh::ExodusII_IO_Helper::write_nodesets(), libMesh::RBEIMEvaluation::write_out_interior_basis_functions(), libMesh::RBEIMEvaluation::write_out_node_basis_functions(), libMesh::RBEIMEvaluation::write_out_side_basis_functions(), write_output_solvedata(), libMesh::System::write_parallel_data(), libMesh::RBConstruction::write_riesz_representors_to_files(), 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::ExodusII_IO_Helper::write_sideset_data(), 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().

115  { return cast_int<processor_id_type>(_communicator.rank()); }
processor_id_type rank() const
const Parallel::Communicator & _communicator

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

218 {
220 }
bool _is_time_dependent
Cache whether or not this function is actually time-dependent.

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

References _point_locator.

Referenced by MeshFunction().

764 {
765  _point_locator->set_close_to_point_tol(tol);
766  _point_locator->set_contains_point_tol(tol);
767 }
std::unique_ptr< PointLocatorBase > _point_locator
A point locator is needed to locate the points in the mesh.

◆ set_subdomain_ids()

void libMesh::MeshFunction::set_subdomain_ids ( const std::set< subdomain_id_type > *  subdomain_ids)

Choose a default list of subdomain ids to be searched for points.

If the provided list pointer is null or if no list has been provided, then all subdomain ids are searched. This list can be overridden on a per-evaluation basis by using the method overrides with a similar argument.

Definition at line 774 of file mesh_function.C.

References _subdomain_ids.

775 {
776  if (subdomain_ids)
777  _subdomain_ids = std::make_unique<std::set<subdomain_id_type>>(*subdomain_ids);
778  else
779  _subdomain_ids.reset();
780 }
std::unique_ptr< std::set< subdomain_id_type > > _subdomain_ids
A default set of subdomain ids in which to search for points.

◆ unset_point_locator_tolerance()

void libMesh::MeshFunction::unset_point_locator_tolerance ( )

Turn off the user-specified PointLocator tolerance.

Definition at line 769 of file mesh_function.C.

References _point_locator.

770 {
771  _point_locator->unset_close_to_point_tol();
772 }
std::unique_ptr< PointLocatorBase > _point_locator
A point locator is needed to locate the points in the mesh.

Member Data Documentation

◆ _communicator

const Parallel::Communicator& libMesh::ParallelObject::_communicator
protectedinherited

◆ _dof_map

const DofMap& libMesh::MeshFunction::_dof_map
protected

Need access to the DofMap of the other system.

Definition at line 389 of file mesh_function.h.

Referenced by _gradient_on_elem(), check_found_elem(), discontinuous_value(), 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 378 of file mesh_function.h.

Referenced by _gradient_on_elem(), 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 184 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 189 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 178 of file function_base.h.

Referenced by clear(), find_element(), and find_elements().

◆ _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 412 of file mesh_function.h.

Referenced by _gradient_on_elem(), disable_out_of_mesh_mode(), discontinuous_value(), enable_out_of_mesh_mode(), find_element(), find_elements(), 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 418 of file mesh_function.h.

Referenced by _gradient_on_elem(), discontinuous_value(), enable_out_of_mesh_mode(), hessian(), and operator()().

◆ _point_locator

std::unique_ptr<PointLocatorBase> libMesh::MeshFunction::_point_locator
protected

A point locator is needed to locate the points in the mesh.

Definition at line 401 of file mesh_function.h.

Referenced by clear(), disable_out_of_mesh_mode(), enable_out_of_mesh_mode(), get_point_locator(), init(), set_point_locator_tolerance(), and unset_point_locator_tolerance().

◆ _subdomain_ids

std::unique_ptr<std::set<subdomain_id_type> > libMesh::MeshFunction::_subdomain_ids
protected

A default set of subdomain ids in which to search for points.

Definition at line 406 of file mesh_function.h.

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

◆ _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 395 of file mesh_function.h.

Referenced by _gradient_on_elem(), discontinuous_gradient(), discontinuous_value(), 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 384 of file mesh_function.h.

Referenced by _gradient_on_elem(), check_found_elem(), discontinuous_value(), hessian(), and operator()().


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