libMesh
mesh_function.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 // C++ includes
21 
22 
23 // Local Includes
24 #include "libmesh/mesh_function.h"
25 #include "libmesh/dense_vector.h"
26 #include "libmesh/equation_systems.h"
27 #include "libmesh/numeric_vector.h"
28 #include "libmesh/dof_map.h"
29 #include "libmesh/point_locator_tree.h"
30 #include "libmesh/fe_base.h"
31 #include "libmesh/fe_interface.h"
32 #include "libmesh/fe_compute_data.h"
33 #include "libmesh/mesh_base.h"
34 #include "libmesh/point.h"
35 #include "libmesh/elem.h"
36 
37 namespace libMesh
38 {
39 
40 
41 //------------------------------------------------------------------
42 // MeshFunction methods
44  const NumericVector<Number> & vec,
45  const DofMap & dof_map,
46  const std::vector<unsigned int> & vars,
47  const FunctionBase<Number> * master) :
48  FunctionBase<Number> (master),
49  ParallelObject (eqn_systems),
50  _eqn_systems (eqn_systems),
51  _vector (vec),
52  _dof_map (dof_map),
53  _system_vars (vars),
54  _point_locator (libmesh_nullptr),
55  _out_of_mesh_mode (false),
56  _out_of_mesh_value ()
57 {
58 }
59 
60 
61 
63  const NumericVector<Number> & vec,
64  const DofMap & dof_map,
65  const unsigned int var,
66  const FunctionBase<Number> * master) :
67  FunctionBase<Number> (master),
68  ParallelObject (eqn_systems),
69  _eqn_systems (eqn_systems),
70  _vector (vec),
71  _dof_map (dof_map),
72  _system_vars (1,var),
74  _out_of_mesh_mode (false),
76 {
77  // std::vector<unsigned int> buf (1);
78  // buf[0] = var;
79  // _system_vars (buf);
80 }
81 
82 
83 
84 
85 
86 
87 
89 {
90  // only delete the point locator when we are the master
91  if (this->_master == libmesh_nullptr)
92  delete this->_point_locator;
93 }
94 
95 
96 
97 
98 void MeshFunction::init (const Trees::BuildType /*point_locator_build_type*/)
99 {
100  // are indices of the desired variable(s) provided?
101  libmesh_assert_greater (this->_system_vars.size(), 0);
102 
103  // Don't do twice...
104  if (this->_initialized)
105  {
107  return;
108  }
109 
110  /*
111  * set up the PointLocator: either someone else
112  * is the master (go and get the address of his
113  * point locator) or this object is the master
114  * (build the point locator on our own).
115  */
116  if (this->_master != libmesh_nullptr)
117  {
118  // we aren't the master
119  const MeshFunction * master =
120  cast_ptr<const MeshFunction *>(this->_master);
121 
122  if (master->_point_locator == libmesh_nullptr)
123  libmesh_error_msg("ERROR: When the master-servant concept is used, the master has to be initialized first!");
124 
125  else
126  {
127  this->_point_locator = master->_point_locator;
128  }
129  }
130  else
131  {
132  // we are the master: build the point locator
133 
134  // constant reference to the other mesh
135  const MeshBase & mesh = this->_eqn_systems.get_mesh();
136 
137  // build the point locator. Only \p TREE version available
138  //UniquePtr<PointLocatorBase> ap (PointLocatorBase::build (TREE, mesh));
139  //this->_point_locator = ap.release();
140  // this->_point_locator = new PointLocatorTree (mesh, point_locator_build_type);
141  this->_point_locator = mesh.sub_point_locator().release();
142 
143  // Point locator no longer needs to be initialized.
144  // this->_point_locator->init();
145  }
146 
147 
148  // ready for use
149  this->_initialized = true;
150 }
151 
152 
153 void
155 {
156  // only delete the point locator when we are the master
157  if ((this->_point_locator != libmesh_nullptr) && (this->_master == libmesh_nullptr))
158  {
159  delete this->_point_locator;
161  }
162  this->_initialized = false;
163 }
164 
165 
166 
168 {
169  FunctionBase<Number> * mf_clone =
171 
172  if (this->initialized())
173  mf_clone->init();
174 
175  return UniquePtr<FunctionBase<Number>>(mf_clone);
176 }
177 
178 
179 
181  const Real time)
182 {
183  libmesh_assert (this->initialized());
184 
185  DenseVector<Number> buf (1);
186  this->operator() (p, time, buf);
187  return buf(0);
188 }
189 
190 std::map<const Elem *, Number> MeshFunction::discontinuous_value (const Point & p,
191  const Real time)
192 {
193  libmesh_assert (this->initialized());
194 
195  std::map<const Elem *, DenseVector<Number>> buffer;
196  this->discontinuous_value (p, time, buffer);
197  std::map<const Elem *, Number> return_value;
198  for (std::map<const Elem *, DenseVector<Number>>::const_iterator it = buffer.begin(); it != buffer.end(); ++it)
199  return_value[it->first] = it->second(0);
200  // NOTE: If no suitable element is found, then the map return_value is empty. This
201  // puts burden on the user of this function but I don't really see a better way.
202  return return_value;
203 }
204 
206  const Real time)
207 {
208  libmesh_assert (this->initialized());
209 
210  std::vector<Gradient> buf (1);
211  this->gradient(p, time, buf);
212  return buf[0];
213 }
214 
215 std::map<const Elem *, Gradient> MeshFunction::discontinuous_gradient (const Point & p,
216  const Real time)
217 {
218  libmesh_assert (this->initialized());
219 
220  std::map<const Elem *, std::vector<Gradient>> buffer;
221  this->discontinuous_gradient (p, time, buffer);
222  std::map<const Elem *, Gradient> return_value;
223  for (std::map<const Elem *, std::vector<Gradient>>::const_iterator it = buffer.begin(); it != buffer.end(); ++it)
224  return_value[it->first] = it->second[0];
225  // NOTE: If no suitable element is found, then the map return_value is empty. This
226  // puts burden on the user of this function but I don't really see a better way.
227  return return_value;
228 }
229 
230 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
232  const Real time)
233 {
234  libmesh_assert (this->initialized());
235 
236  std::vector<Tensor> buf (1);
237  this->hessian(p, time, buf);
238  return buf[0];
239 }
240 #endif
241 
243  const Real time,
244  DenseVector<Number> & output)
245 {
246  this->operator() (p, time, output, libmesh_nullptr);
247 }
248 
250  const Real,
251  DenseVector<Number> & output,
252  const std::set<subdomain_id_type> * subdomain_ids)
253 {
254  libmesh_assert (this->initialized());
255 
256  const Elem * element = this->find_element(p,subdomain_ids);
257 
258  if (!element)
259  {
260  // We'd better be in out_of_mesh_mode if we couldn't find an
261  // element in the mesh
263  output = _out_of_mesh_value;
264  }
265  else
266  {
267  // resize the output vector to the number of output values
268  // that the user told us
269  output.resize (cast_int<unsigned int>
270  (this->_system_vars.size()));
271 
272 
273  {
274  const unsigned int dim = element->dim();
275 
276 
277  /*
278  * Get local coordinates to feed these into compute_data().
279  * Note that the fe_type can safely be used from the 0-variable,
280  * since the inverse mapping is the same for all FEFamilies
281  */
282  const Point mapped_point (FEInterface::inverse_map (dim,
283  this->_dof_map.variable_type(0),
284  element,
285  p));
286 
287 
288  // loop over all vars
289  for (std::size_t index=0; index < this->_system_vars.size(); index++)
290  {
291  /*
292  * the data for this variable
293  */
294  const unsigned int var = _system_vars[index];
295 
296  if (var == libMesh::invalid_uint)
297  {
299  index < _out_of_mesh_value.size());
300  output(index) = _out_of_mesh_value(index);
301  continue;
302  }
303 
304  const FEType & fe_type = this->_dof_map.variable_type(var);
305 
310  {
311  FEComputeData data (this->_eqn_systems, mapped_point);
312 
313  FEInterface::compute_data (dim, fe_type, element, data);
314 
315  // where the solution values for the var-th variable are stored
316  std::vector<dof_id_type> dof_indices;
317  this->_dof_map.dof_indices (element, dof_indices, var);
318 
319  // interpolate the solution
320  {
321  Number value = 0.;
322 
323  for (std::size_t i=0; i<dof_indices.size(); i++)
324  value += this->_vector(dof_indices[i]) * data.shape[i];
325 
326  output(index) = value;
327  }
328 
329  }
330 
331  // next variable
332  }
333  }
334  }
335 }
336 
337 
339  const Real time,
340  std::map<const Elem *, DenseVector<Number>> & output)
341 {
342  this->discontinuous_value (p, time, output, libmesh_nullptr);
343 }
344 
345 
346 
348  const Real,
349  std::map<const Elem *, DenseVector<Number>> & output,
350  const std::set<subdomain_id_type> * subdomain_ids)
351 {
352  libmesh_assert (this->initialized());
353 
354  // clear the output map
355  output.clear();
356 
357  // get the candidate elements
358  std::set<const Elem *> candidate_element = this->find_elements(p,subdomain_ids);
359 
360  // loop through all candidates, if the set is empty this function will return an
361  // empty map
362  for (std::set<const Elem *>::const_iterator it = candidate_element.begin(); it != candidate_element.end(); ++it)
363  {
364  const Elem * element = (*it);
365 
366  const unsigned int dim = element->dim();
367 
368  // define a temporary vector to store all values
369  DenseVector<Number> temp_output (cast_int<unsigned int>(this->_system_vars.size()));
370 
371  /*
372  * Get local coordinates to feed these into compute_data().
373  * Note that the fe_type can safely be used from the 0-variable,
374  * since the inverse mapping is the same for all FEFamilies
375  */
376  const Point mapped_point (FEInterface::inverse_map (dim,
377  this->_dof_map.variable_type(0),
378  element,
379  p));
380 
381 
382  // loop over all vars
383  for (std::size_t index=0; index < this->_system_vars.size(); index++)
384  {
385  /*
386  * the data for this variable
387  */
388  const unsigned int var = _system_vars[index];
389 
390  if (var == libMesh::invalid_uint)
391  {
393  index < _out_of_mesh_value.size());
394  temp_output(index) = _out_of_mesh_value(index);
395  continue;
396  }
397 
398  const FEType & fe_type = this->_dof_map.variable_type(var);
399 
404  {
405  FEComputeData data (this->_eqn_systems, mapped_point);
406 
407  FEInterface::compute_data (dim, fe_type, element, data);
408 
409  // where the solution values for the var-th variable are stored
410  std::vector<dof_id_type> dof_indices;
411  this->_dof_map.dof_indices (element, dof_indices, var);
412 
413  // interpolate the solution
414  {
415  Number value = 0.;
416 
417  for (std::size_t i=0; i<dof_indices.size(); i++)
418  value += this->_vector(dof_indices[i]) * data.shape[i];
419 
420  temp_output(index) = value;
421  }
422 
423  }
424 
425  // next variable
426  }
427 
428  // Insert temp_output into output
429  output[element] = temp_output;
430  }
431 }
432 
433 
434 
436  const Real,
437  std::vector<Gradient> & output,
438  const std::set<subdomain_id_type> * subdomain_ids)
439 {
440  libmesh_assert (this->initialized());
441 
442  const Elem * element = this->find_element(p,subdomain_ids);
443 
444  if (!element)
445  {
446  output.resize(0);
447  }
448  else
449  {
450  // resize the output vector to the number of output values
451  // that the user told us
452  output.resize (this->_system_vars.size());
453 
454 
455  {
456  const unsigned int dim = element->dim();
457 
458 
459  /*
460  * Get local coordinates to feed these into compute_data().
461  * Note that the fe_type can safely be used from the 0-variable,
462  * since the inverse mapping is the same for all FEFamilies
463  */
464  const Point mapped_point (FEInterface::inverse_map (dim,
465  this->_dof_map.variable_type(0),
466  element,
467  p));
468 
469  std::vector<Point> point_list (1, mapped_point);
470 
471  // loop over all vars
472  for (std::size_t index=0; index < this->_system_vars.size(); index++)
473  {
474  /*
475  * the data for this variable
476  */
477  const unsigned int var = _system_vars[index];
478 
479  if (var == libMesh::invalid_uint)
480  {
482  index < _out_of_mesh_value.size());
483  output[index] = Gradient(_out_of_mesh_value(index));
484  continue;
485  }
486 
487  const FEType & fe_type = this->_dof_map.variable_type(var);
488 
489  UniquePtr<FEBase> point_fe (FEBase::build(dim, fe_type));
490  const std::vector<std::vector<RealGradient>> & dphi = point_fe->get_dphi();
491  point_fe->reinit(element, &point_list);
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 (std::size_t i=0; i<dof_indices.size(); i++)
501  grad.add_scaled(dphi[i][0], this->_vector(dof_indices[i]));
502 
503  output[index] = grad;
504  }
505  }
506  }
507 }
508 
509 
511  const Real time,
512  std::map<const Elem *, std::vector<Gradient>> & output)
513 {
514  this->discontinuous_gradient (p, time, output, libmesh_nullptr);
515 }
516 
517 
518 
520  const Real,
521  std::map<const Elem *, std::vector<Gradient>> & output,
522  const std::set<subdomain_id_type> * subdomain_ids)
523 {
524  libmesh_assert (this->initialized());
525 
526  // clear the output map
527  output.clear();
528 
529  // get the candidate elements
530  std::set<const Elem *> candidate_element = this->find_elements(p,subdomain_ids);
531 
532  // loop through all candidates, if the set is empty this function will return an
533  // empty map
534  for (std::set<const Elem *>::const_iterator it = candidate_element.begin(); it != candidate_element.end(); ++it)
535  {
536  const Elem * element = (*it);
537 
538  const unsigned int dim = element->dim();
539 
540  // define a temporary vector to store all values
541  std::vector<Gradient> temp_output (cast_int<unsigned int>(this->_system_vars.size()));
542 
543  /*
544  * Get local coordinates to feed these into compute_data().
545  * Note that the fe_type can safely be used from the 0-variable,
546  * since the inverse mapping is the same for all FEFamilies
547  */
548  const Point mapped_point (FEInterface::inverse_map (dim,
549  this->_dof_map.variable_type(0),
550  element,
551  p));
552 
553 
554  // loop over all vars
555  std::vector<Point> point_list (1, mapped_point);
556  for (std::size_t index = 0 ; index < this->_system_vars.size(); ++index)
557  {
558  /*
559  * the data for this variable
560  */
561  const unsigned int var = _system_vars[index];
562 
563  if (var == libMesh::invalid_uint)
564  {
566  index < _out_of_mesh_value.size());
567  temp_output[index] = Gradient(_out_of_mesh_value(index));
568  continue;
569  }
570 
571  const FEType & fe_type = this->_dof_map.variable_type(var);
572 
573  UniquePtr<FEBase> point_fe (FEBase::build(dim, fe_type));
574  const std::vector<std::vector<RealGradient>> & dphi = point_fe->get_dphi();
575  point_fe->reinit(element, &point_list);
576 
577  // where the solution values for the var-th variable are stored
578  std::vector<dof_id_type> dof_indices;
579  this->_dof_map.dof_indices (element, dof_indices, var);
580 
581  Gradient grad(0.);
582 
583  for (std::size_t i = 0; i < dof_indices.size(); ++i)
584  grad.add_scaled(dphi[i][0], this->_vector(dof_indices[i]));
585 
586  temp_output[index] = grad;
587 
588  // next variable
589  }
590 
591  // Insert temp_output into output
592  output[element] = temp_output;
593  }
594 }
595 
596 
597 
598 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
599 void MeshFunction::hessian (const Point & p,
600  const Real,
601  std::vector<Tensor> & output,
602  const std::set<subdomain_id_type> * subdomain_ids)
603 {
604  libmesh_assert (this->initialized());
605 
606  const Elem * element = this->find_element(p,subdomain_ids);
607 
608  if (!element)
609  {
610  output.resize(0);
611  }
612  else
613  {
614  // resize the output vector to the number of output values
615  // that the user told us
616  output.resize (this->_system_vars.size());
617 
618 
619  {
620  const unsigned int dim = element->dim();
621 
622 
623  /*
624  * Get local coordinates to feed these into compute_data().
625  * Note that the fe_type can safely be used from the 0-variable,
626  * since the inverse mapping is the same for all FEFamilies
627  */
628  const Point mapped_point (FEInterface::inverse_map (dim,
629  this->_dof_map.variable_type(0),
630  element,
631  p));
632 
633  std::vector<Point> point_list (1, mapped_point);
634 
635  // loop over all vars
636  for (std::size_t index=0; index < this->_system_vars.size(); index++)
637  {
638  /*
639  * the data for this variable
640  */
641  const unsigned int var = _system_vars[index];
642 
643  if (var == libMesh::invalid_uint)
644  {
646  index < _out_of_mesh_value.size());
647  output[index] = Tensor(_out_of_mesh_value(index));
648  continue;
649  }
650  const FEType & fe_type = this->_dof_map.variable_type(var);
651 
652  UniquePtr<FEBase> point_fe (FEBase::build(dim, fe_type));
653  const std::vector<std::vector<RealTensor>> & d2phi =
654  point_fe->get_d2phi();
655  point_fe->reinit(element, &point_list);
656 
657  // where the solution values for the var-th variable are stored
658  std::vector<dof_id_type> dof_indices;
659  this->_dof_map.dof_indices (element, dof_indices, var);
660 
661  // interpolate the solution
662  Tensor hess;
663 
664  for (std::size_t i=0; i<dof_indices.size(); i++)
665  hess.add_scaled(d2phi[i][0], this->_vector(dof_indices[i]));
666 
667  output[index] = hess;
668  }
669  }
670  }
671 }
672 #endif
673 
675  const std::set<subdomain_id_type> * subdomain_ids) const
676 {
677  /* Ensure that in the case of a master mesh function, the
678  out-of-mesh mode is enabled either for both or for none. This is
679  important because the out-of-mesh mode is also communicated to
680  the point locator. Since this is time consuming, enable it only
681  in debug mode. */
682 #ifdef DEBUG
683  if (this->_master != libmesh_nullptr)
684  {
685  const MeshFunction * master =
686  cast_ptr<const MeshFunction *>(this->_master);
688  libmesh_error_msg("ERROR: If you use out-of-mesh-mode in connection with master mesh " \
689  << "functions, you must enable out-of-mesh mode for both the master and the slave mesh function.");
690  }
691 #endif
692 
693  // locate the point in the other mesh
694  const Elem * element = this->_point_locator->operator()(p,subdomain_ids);
695 
696  // If we have an element, but it's not a local element, then we
697  // either need to have a serialized vector or we need to find a
698  // local element sharing the same point.
699  if (element &&
700  (element->processor_id() != this->processor_id()) &&
701  _vector.type() != SERIAL)
702  {
703  // look for a local element containing the point
704  std::set<const Elem *> point_neighbors;
705  element->find_point_neighbors(p, point_neighbors);
706  element = libmesh_nullptr;
707  std::set<const Elem *>::const_iterator it = point_neighbors.begin();
708  const std::set<const Elem *>::const_iterator end = point_neighbors.end();
709  for (; it != end; ++it)
710  {
711  const Elem * elem = *it;
712  if (elem->processor_id() == this->processor_id())
713  {
714  element = elem;
715  break;
716  }
717  }
718  }
719 
720  return element;
721 }
722 
723 std::set<const Elem *> MeshFunction::find_elements(const Point & p,
724  const std::set<subdomain_id_type> * subdomain_ids) const
725 {
726  /* Ensure that in the case of a master mesh function, the
727  out-of-mesh mode is enabled either for both or for none. This is
728  important because the out-of-mesh mode is also communicated to
729  the point locator. Since this is time consuming, enable it only
730  in debug mode. */
731 #ifdef DEBUG
732  if (this->_master != libmesh_nullptr)
733  {
734  const MeshFunction * master =
735  cast_ptr<const MeshFunction *>(this->_master);
737  libmesh_error_msg("ERROR: If you use out-of-mesh-mode in connection with master mesh " \
738  << "functions, you must enable out-of-mesh mode for both the master and the slave mesh function.");
739  }
740 #endif
741 
742  // locate the point in the other mesh
743  std::set<const Elem *> candidate_elements;
744  std::set<const Elem *> final_candidate_elements;
745  this->_point_locator->operator()(p,candidate_elements,subdomain_ids);
746  for (std::set<const Elem *>::const_iterator elem_it = candidate_elements.begin(); elem_it != candidate_elements.end(); ++elem_it)
747  {
748  const Elem * element = (*elem_it);
749  // If we have an element, but it's not a local element, then we
750  // either need to have a serialized vector or we need to find a
751  // local element sharing the same point.
752  if (element &&
753  (element->processor_id() != this->processor_id()) &&
754  _vector.type() != SERIAL)
755  {
756  // look for a local element containing the point
757  std::set<const Elem *> point_neighbors;
758  element->find_point_neighbors(p, point_neighbors);
759  std::set<const Elem *>::const_iterator it = point_neighbors.begin();
760  const std::set<const Elem *>::const_iterator end = point_neighbors.end();
761  for (; it != end; ++it)
762  {
763  const Elem * elem = *it;
764  if (elem->processor_id() == this->processor_id())
765  {
766  final_candidate_elements.insert(elem);
767  break;
768  }
769  }
770  }
771  else
772  final_candidate_elements.insert(element);
773  }
774 
775  return final_candidate_elements;
776 }
777 
779 {
780  libmesh_assert (this->initialized());
781  return *_point_locator;
782 }
783 
785 {
786  libmesh_assert (this->initialized());
788  _out_of_mesh_mode = true;
790 }
791 
793 {
794  DenseVector<Number> v(1);
795  v(0) = value;
796  this->enable_out_of_mesh_mode(v);
797 }
798 
800 {
801  libmesh_assert (this->initialized());
803  _out_of_mesh_mode = false;
804 }
805 
807 {
809 }
810 
812 {
814 }
815 
816 } // namespace libMesh
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
virtual void disable_out_of_mesh_mode()=0
Disables out-of-mesh mode (default).
const Elem * find_element(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=libmesh_nullptr) const
Helper function to reduce code duplication.
This is the EquationSystems class.
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.
Definition: libmesh.h:184
virtual void set_close_to_point_tol(Real close_to_point_tol)
Set a tolerance to use when determining if a point is contained within the mesh.
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1697
PointLocatorBase * _point_locator
A point locator is needed to locate the points in the mesh.
void add_scaled(const TypeTensor< T2 > &, const T)
Add a scaled tensor to this tensor without creating a temporary.
Definition: type_tensor.h:777
void add_scaled(const TypeVector< T2 > &, const T)
Add a scaled value to this vector without creating a temporary.
Definition: type_vector.h:624
unsigned int dim
class FEComputeData hides arbitrary data to be passed to and from children of FEBase through the FEIn...
const DofMap & _dof_map
Need access to the DofMap of the other system.
void resize(const unsigned int n)
Resize the vector.
Definition: dense_vector.h:350
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
MeshBase & mesh
const FunctionBase * _master
Const pointer to our master, initialized to NULL.
bool _initialized
When init() was called so that everything is ready for calls to operator() (...), then this bool is t...
const class libmesh_nullptr_t libmesh_nullptr
Gradient gradient(const Point &p, const Real time=0.)
const PointLocatorBase & get_point_locator(void) const
IterBase * end
Also have a polymorphic pointer to the end object, this prevents iterating past the end...
~MeshFunction()
Destructor.
Definition: mesh_function.C:88
The libMesh namespace provides an interface to certain functionality in the library.
Number operator()(const Point &p, const Real time=0.) libmesh_override
void find_point_neighbors(const Point &p, std::set< const Elem * > &neighbor_set) const
This function finds all active elements (including this one) which are in the same manifold as this e...
Definition: elem.C:605
virtual void init()
The actual initialization process.
Definition: function_base.h:81
virtual unsigned int size() const libmesh_override
Definition: dense_vector.h:87
void disable_out_of_mesh_mode(void)
Disables out-of-mesh mode.
This is the MeshBase class.
Definition: mesh_base.h:68
libmesh_assert(j)
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:167
DenseVector< Number > _out_of_mesh_value
Value to return outside the mesh if out-of-mesh mode is enabled.
virtual UniquePtr< FunctionBase< Number > > clone() const libmesh_override
void enable_out_of_mesh_mode(const DenseVector< Number > &value)
Enables out-of-mesh mode.
std::map< const Elem *, Number > discontinuous_value(const Point &p, const Real time=0.)
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...
Definition: fe_interface.C:827
NumberVectorValue Gradient
void unset_point_locator_tolerance()
Turn off the user-specified PointLocator tolerance.
virtual void clear() libmesh_override
Clears the function.
static Point inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: fe_interface.C:569
This is the base class for point locators.
This class forms the base class for all other classes that are expected to be implemented in parallel...
std::vector< Number > shape
Storage for the computed shape function values.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::set< const Elem * > find_elements(const Point &p, const std::set< subdomain_id_type > *subdomain_ids=libmesh_nullptr) const
NumberTensorValue Tensor
static const bool value
Definition: xdr_io.C:108
ParallelType type() const
virtual unsigned int dim() const =0
BuildType
enum defining how to build the tree.
Definition: tree_base.h:55
const MeshBase & get_mesh() const
virtual void unset_close_to_point_tol()
Specify that we do not want to use a user-specified tolerance to determine if a point is contained wi...
This class provides function-like objects for data distributed over a mesh.
Definition: mesh_function.h:53
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...
IterBase * data
Ideally this private member data should have protected access.
virtual void enable_out_of_mesh_mode()=0
Enables out-of-mesh mode.
std::map< const Elem *, Gradient > discontinuous_gradient(const Point &p, const Real time=0.)
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
Tensor hessian(const Point &p, const Real time=0.)
UniquePtr< PointLocatorBase > sub_point_locator() const
Definition: mesh_base.C:534
virtual void init() libmesh_override
Override the FunctionBase::init() member function by calling our own and specifying the Trees::NODES ...
Definition: mesh_function.h:96
MeshFunction(const EquationSystems &eqn_systems, const NumericVector< Number > &vec, const DofMap &dof_map, const std::vector< unsigned int > &vars, const FunctionBase< Number > *master=libmesh_nullptr)
Constructor for mesh based functions with vectors as return value.
Definition: mesh_function.C:43
This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space.
processor_id_type processor_id() const
const std::vector< unsigned int > _system_vars
The indices of the variables within the other system for which data are to be gathered.
processor_id_type processor_id() const
Definition: dof_object.h:694
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:1917
static UniquePtr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.