libMesh
fem_context.h
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 #ifndef LIBMESH_FEM_CONTEXT_H
21 #define LIBMESH_FEM_CONTEXT_H
22 
23 // Local Includes
24 #include "libmesh/diff_context.h"
25 #include "libmesh/id_types.h"
26 #include "libmesh/fe_type.h"
27 #include "libmesh/fe_base.h"
28 #include "libmesh/vector_value.h"
29 
30 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
31 #include "libmesh/tensor_value.h"
32 #endif
33 
34 // C++ includes
35 #include <map>
36 
37 namespace libMesh
38 {
39 
40 // Forward Declarations
41 class BoundaryInfo;
42 class Elem;
43 template <typename T> class FEGenericBase;
44 typedef FEGenericBase<Real> FEBase;
45 class QBase;
46 class Point;
47 template <typename T> class NumericVector;
48 
61 class FEMContext : public DiffContext
62 {
63 public:
64 
68  explicit
69  FEMContext (const System & sys);
70 
75  explicit
76  FEMContext (const System & sys, int extra_quadrature_order);
77 
81  virtual ~FEMContext ();
82 
87 
94 #ifdef LIBMESH_ENABLE_DEPRECATED
95  std::vector<boundary_id_type> side_boundary_ids() const;
96 #endif
97 
101  void side_boundary_ids(std::vector<boundary_id_type> & vec_to_fill) const;
102 
109  Number interior_value(unsigned int var, unsigned int qp) const;
110 
117  Number side_value(unsigned int var, unsigned int qp) const;
118 
125  Number point_value(unsigned int var, const Point & p) const;
126 
133  Gradient interior_gradient(unsigned int var, unsigned int qp) const;
134 
141  Gradient side_gradient(unsigned int var, unsigned int qp) const;
142 
149  Gradient point_gradient(unsigned int var, const Point & p) const;
150 
151 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
152 
158  Tensor interior_hessian(unsigned int var, unsigned int qp) const;
159 
166  Tensor side_hessian(unsigned int var, unsigned int qp) const;
167 
174  Tensor point_hessian(unsigned int var, const Point & p) const;
175 
176 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
177 
184  Number fixed_interior_value(unsigned int var, unsigned int qp) const;
185 
192  Number fixed_side_value(unsigned int var, unsigned int qp) const;
193 
200  Number fixed_point_value(unsigned int var, const Point & p) const;
201 
208  Gradient fixed_interior_gradient(unsigned int var, unsigned int qp) const;
209 
216  Gradient fixed_side_gradient(unsigned int var, unsigned int qp) const;
217 
224  Gradient fixed_point_gradient(unsigned int var, const Point & p) const;
225 
226 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
227 
233  Tensor fixed_interior_hessian(unsigned int var, unsigned int qp) const;
234 
241  Tensor fixed_side_hessian(unsigned int var, unsigned int qp) const;
242 
249  Tensor fixed_point_hessian (unsigned int var, const Point & p) const;
250 
251 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
252 
261  template<typename OutputShape>
262  void get_element_fe( unsigned int var, FEGenericBase<OutputShape> *& fe ) const
263  { this->get_element_fe<OutputShape>(var,fe,this->get_dim()); }
264 
273  FEBase * get_element_fe( unsigned int var ) const
274  { return this->get_element_fe(var,this->get_dim()); }
275 
280  template<typename OutputShape>
281  void get_element_fe( unsigned int var, FEGenericBase<OutputShape> *& fe,
282  unsigned char dim ) const;
283 
288  FEBase * get_element_fe( unsigned int var, unsigned char dim ) const;
289 
298  template<typename OutputShape>
299  void get_side_fe( unsigned int var, FEGenericBase<OutputShape> *& fe ) const
300  { this->get_side_fe<OutputShape>(var,fe,this->get_dim()); }
301 
310  FEBase * get_side_fe( unsigned int var ) const
311  { return this->get_side_fe(var,this->get_dim()); }
312 
317  template<typename OutputShape>
318  void get_side_fe( unsigned int var, FEGenericBase<OutputShape> *& fe,
319  unsigned char dim ) const;
320 
325  FEBase * get_side_fe( unsigned int var, unsigned char dim ) const;
326 
330  template<typename OutputShape>
331  void get_edge_fe( unsigned int var, FEGenericBase<OutputShape> *& fe ) const;
332 
336  FEBase * get_edge_fe( unsigned int var ) const;
337 
344  template<typename OutputType>
345  void interior_value(unsigned int var,
346  unsigned int qp,
347  OutputType & u) const;
348 
353  template<typename OutputType>
354  void interior_values(unsigned int var,
355  const NumericVector<Number> & _system_vector,
356  std::vector<OutputType> & interior_values_vector) const;
357 
364  template<typename OutputType>
365  void side_value(unsigned int var,
366  unsigned int qp,
367  OutputType & u) const;
368 
373  template<typename OutputType>
374  void side_values(unsigned int var,
375  const NumericVector<Number> & _system_vector,
376  std::vector<OutputType> & side_values_vector) const;
377 
387  template<typename OutputType>
388  void point_value(unsigned int var,
389  const Point & p,
390  OutputType & u,
391  const Real tolerance = TOLERANCE) const;
392 
399  template<typename OutputType>
400  void interior_gradient(unsigned int var, unsigned int qp,
401  OutputType & du) const;
402 
409  template<typename OutputType>
410  void interior_gradients(unsigned int var,
411  const NumericVector<Number> & _system_vector,
412  std::vector<OutputType> & interior_gradients_vector) const;
413 
420  template<typename OutputType>
421  void side_gradient(unsigned int var,
422  unsigned int qp,
423  OutputType & du) const;
424 
431  template<typename OutputType>
432  void side_gradients(unsigned int var,
433  const NumericVector<Number> & _system_vector,
434  std::vector<OutputType> & side_gradients_vector) const;
435 
445  template<typename OutputType>
446  void point_gradient(unsigned int var,
447  const Point & p,
448  OutputType & grad_u,
449  const Real tolerance = TOLERANCE) const;
450 
451 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
452 
458  template<typename OutputType>
459  void interior_hessian(unsigned int var,
460  unsigned int qp,
461  OutputType & d2u) const;
462 
468  template<typename OutputType>
469  void interior_hessians(unsigned int var,
470  const NumericVector<Number> & _system_vector,
471  std::vector<OutputType> & d2u_vals) const;
472 
479  template<typename OutputType>
480  void side_hessian(unsigned int var,
481  unsigned int qp,
482  OutputType & d2u) const;
483 
489  template<typename OutputType>
490  void side_hessians(unsigned int var,
491  const NumericVector<Number> & _system_vector,
492  std::vector<OutputType> & d2u_vals) const;
493 
503  template<typename OutputType>
504  void point_hessian(unsigned int var,
505  const Point & p,
506  OutputType & hess_u,
507  const Real tolerance = TOLERANCE) const;
508 
509 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
510 
516  template<typename OutputType>
517  void interior_rate(unsigned int var,
518  unsigned int qp,
519  OutputType & u) const;
520 
525  template<typename OutputType>
526  void side_rate(unsigned int var,
527  unsigned int qp,
528  OutputType & u) const;
529 
534  template<typename OutputType>
535  void point_rate(unsigned int var,
536  const Point & p,
537  OutputType & u) const;
538 
544  template<typename OutputType>
545  void interior_accel(unsigned int var,
546  unsigned int qp,
547  OutputType & u) const;
548 
553  template<typename OutputType>
554  void side_accel(unsigned int var,
555  unsigned int qp,
556  OutputType & u) const;
557 
562  template<typename OutputType>
563  void point_accel(unsigned int var,
564  const Point & p,
565  OutputType & u) const;
566 
573  template<typename OutputType>
574  void fixed_interior_value(unsigned int var,
575  unsigned int qp,
576  OutputType & u) const;
577 
584  template<typename OutputType>
585  void fixed_side_value(unsigned int var,
586  unsigned int qp,
587  OutputType & u) const;
588 
598  template<typename OutputType>
599  void fixed_point_value(unsigned int var,
600  const Point & p,
601  OutputType & u,
602  const Real tolerance = TOLERANCE) const;
603 
610  template<typename OutputType>
611  void fixed_interior_gradient(unsigned int var,
612  unsigned int qp,
613  OutputType & grad_u) const;
614 
621  template<typename OutputType>
622  void fixed_side_gradient(unsigned int var,
623  unsigned int qp,
624  OutputType & grad_u) const;
625 
635  template<typename OutputType>
636  void fixed_point_gradient(unsigned int var,
637  const Point & p,
638  OutputType & grad_u,
639  const Real tolerance = TOLERANCE) const;
640 
641 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
642 
648  template<typename OutputType>
649  void fixed_interior_hessian(unsigned int var,
650  unsigned int qp,
651  OutputType & hess_u) const;
652 
659  template<typename OutputType>
660  void fixed_side_hessian(unsigned int var,
661  unsigned int qp,
662  OutputType & hess_u) const;
663 
673  template<typename OutputType>
674  void fixed_point_hessian(unsigned int var,
675  const Point & p,
676  OutputType & hess_u,
677  const Real tolerance = TOLERANCE) const;
678 
679 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES
680 
685  template<typename OutputType>
686  void interior_curl(unsigned int var,
687  unsigned int qp,
688  OutputType & curl_u) const;
689 
697  template<typename OutputType>
698  void point_curl(unsigned int var,
699  const Point & p,
700  OutputType & curl_u,
701  const Real tolerance = TOLERANCE) const;
702 
707  template<typename OutputType>
708  void interior_div(unsigned int var,
709  unsigned int qp,
710  OutputType & div_u) const;
711 
712  // should be protected:
718  virtual void elem_reinit(Real theta) libmesh_override;
719 
725  virtual void elem_side_reinit(Real theta) libmesh_override;
726 
732  virtual void elem_edge_reinit(Real theta) libmesh_override;
733 
738  virtual void nonlocal_reinit(Real theta) libmesh_override;
739 
743  virtual void pre_fe_reinit(const System &,
744  const Elem * e);
745 
749  virtual void elem_fe_reinit(const std::vector<Point> * const pts = libmesh_nullptr);
750 
754  virtual void side_fe_reinit();
755 
759  virtual void edge_fe_reinit();
760 
765  const QBase & get_element_qrule() const
766  { return this->get_element_qrule(this->get_elem_dim()); }
767 
772  const QBase & get_side_qrule() const
773  { return this->get_side_qrule(this->get_elem_dim()); }
774 
778  const QBase & get_element_qrule( unsigned char dim ) const
780  return *(this->_element_qrule[dim]); }
781 
785  const QBase & get_side_qrule( unsigned char dim ) const
786  {
788  return *(this->_side_qrule[dim]);
789  }
790 
794  const QBase & get_edge_qrule() const
795  { return *(this->_edge_qrule); }
796 
805  virtual void set_mesh_system(System * sys)
806  { this->_mesh_sys = sys; }
807 
811  const System * get_mesh_system() const
812  { return this->_mesh_sys; }
813 
818  { return this->_mesh_sys; }
819 
823  unsigned int get_mesh_x_var() const
824  { return _mesh_x_var; }
825 
831  void set_mesh_x_var(unsigned int x_var)
832  { _mesh_x_var = x_var; }
833 
837  unsigned int get_mesh_y_var() const
838  { return _mesh_y_var; }
839 
845  void set_mesh_y_var(unsigned int y_var)
846  { _mesh_y_var = y_var; }
847 
851  unsigned int get_mesh_z_var() const
852  { return _mesh_z_var; }
853 
859  void set_mesh_z_var(unsigned int z_var)
860  { _mesh_z_var = z_var; }
861 
865  bool has_elem() const
866  { return (this->_elem != libmesh_nullptr); }
867 
871  const Elem & get_elem() const
872  { libmesh_assert(this->_elem);
873  return *(this->_elem); }
874 
879  { libmesh_assert(this->_elem);
880  return *(const_cast<Elem *>(this->_elem)); }
881 
885  unsigned char get_side() const
886  { return side; }
887 
891  unsigned char get_edge() const
892  { return edge; }
893 
899  unsigned char get_dim() const
900  { return this->_dim; }
901 
906  unsigned char get_elem_dim() const
907  { return _elem_dim; }
908 
913  const std::set<unsigned char> & elem_dimensions() const
914  { return _elem_dims; }
915 
921  void elem_position_set(Real theta);
922 
927  void elem_position_get();
928 
933  enum AlgebraicType { NONE = 0, // Do not reinitialize dof_indices
934  DOFS_ONLY, // Reinitialize dof_indices, not
935  // algebraic structures
936  CURRENT, // Use dof_indices, current solution
937  OLD }; // Use old_dof_indices, custom solution
938 
947  { _atype = atype; }
948 
949  /*
950  * Get the current AlgebraicType setting
951  */
952  AlgebraicType algebraic_type() const { return _atype; }
953 
962  { _custom_solution = custom_sol; }
963 
968 
973 
977  unsigned char side;
978 
982  unsigned char edge;
983 
984 protected:
985 
990 
995 
999  template<typename OutputShape>
1001  const Point & p,
1002  const Real tolerance = TOLERANCE) const;
1003 
1006 
1007 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1008  mutable bool _real_fe_is_inf;
1009  mutable bool _real_grad_fe_is_inf;
1010 #endif
1011 
1012  template<typename OutputShape>
1013  FEGenericBase<OutputShape> * cached_fe( const unsigned int elem_dim,
1014  const FEType fe_type ) const;
1015 
1019  void set_elem( const Elem * e );
1020 
1021  // gcc-3.4, oracle 12.3 require this typedef to be public
1022  // in order to use it in a return type
1023 public:
1024 
1028  typedef const DenseSubVector<Number> & (DiffContext::*diff_subsolution_getter)(unsigned int) const;
1029 
1030 protected:
1034  template <typename OutputType>
1035  struct FENeeded
1036  {
1037  // Rank decrementer helper types
1040 
1041  // Typedefs for "Value getter" function pointer
1044  typedef void (FEMContext::*value_getter) (unsigned int, value_base *&, unsigned char) const;
1045 
1046  // Typedefs for "Grad getter" function pointer
1049  typedef void (FEMContext::*grad_getter) (unsigned int, grad_base *&, unsigned char) const;
1050 
1051  // Typedefs for "Hessian getter" function pointer
1054  typedef void (FEMContext::*hess_getter) (unsigned int, hess_base *&, unsigned char) const;
1055  };
1056 
1057 
1058 
1063  template<typename OutputType,
1064  typename FENeeded<OutputType>::value_getter fe_getter,
1065  diff_subsolution_getter subsolution_getter>
1066  void some_value(unsigned int var, unsigned int qp, OutputType & u) const;
1067 
1072  template<typename OutputType,
1073  typename FENeeded<OutputType>::grad_getter fe_getter,
1074  diff_subsolution_getter subsolution_getter>
1075  void some_gradient(unsigned int var, unsigned int qp, OutputType & u) const;
1076 
1077 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
1078 
1082  template<typename OutputType,
1083  typename FENeeded<OutputType>::hess_getter fe_getter,
1084  diff_subsolution_getter subsolution_getter>
1085  void some_hessian(unsigned int var, unsigned int qp, OutputType & u) const;
1086 #endif
1087 
1093  std::vector<std::map<FEType, FEAbstract *>> _element_fe;
1094  std::vector<std::map<FEType, FEAbstract *>> _side_fe;
1095  std::map<FEType, FEAbstract *> _edge_fe;
1096 
1097 
1104  std::vector<std::vector<FEAbstract *>> _element_fe_var;
1105  std::vector<std::vector<FEAbstract *>> _side_fe_var;
1106  std::vector<FEAbstract *> _edge_fe_var;
1107 
1113 
1117  const Elem * _elem;
1118 
1122  unsigned char _dim;
1123 
1127  unsigned char _elem_dim;
1128 
1133  std::set<unsigned char> _elem_dims;
1134 
1141  std::vector<QBase *> _element_qrule;
1142 
1149  std::vector<QBase *> _side_qrule;
1150 
1159 
1164 
1165 private:
1169  void init_internal_data(const System & sys);
1170 
1179  void _do_elem_position_set(Real theta);
1180 
1186  void _update_time_from_system(Real theta);
1187 };
1188 
1189 
1190 
1191 // ------------------------------------------------------------
1192 // FEMContext inline methods
1193 
1194 inline
1196 {
1197  if (_mesh_sys)
1198  this->_do_elem_position_set(theta);
1199 }
1200 
1201 template<typename OutputShape>
1202 inline
1204  unsigned char dim ) const
1205 {
1206  libmesh_assert( !_element_fe_var[dim].empty() );
1207  libmesh_assert_less ( var, (_element_fe_var[dim].size() ) );
1208  fe = cast_ptr<FEGenericBase<OutputShape> *>( (_element_fe_var[dim][var] ) );
1209 }
1210 
1211 inline
1212 FEBase * FEMContext::get_element_fe( unsigned int var, unsigned char dim ) const
1213 {
1214  libmesh_assert( !_element_fe_var[dim].empty() );
1215  libmesh_assert_less ( var, (_element_fe_var[dim].size() ) );
1216  return cast_ptr<FEBase *>( (_element_fe_var[dim][var] ) );
1217 }
1218 
1219 template<typename OutputShape>
1220 inline
1222  unsigned char dim ) const
1223 {
1224  libmesh_assert( !_side_fe_var[dim].empty() );
1225  libmesh_assert_less ( var, (_side_fe_var[dim].size() ) );
1226  fe = cast_ptr<FEGenericBase<OutputShape> *>( (_side_fe_var[dim][var] ) );
1227 }
1228 
1229 inline
1230 FEBase * FEMContext::get_side_fe( unsigned int var, unsigned char dim ) const
1231 {
1232  libmesh_assert( !_side_fe_var[dim].empty() );
1233  libmesh_assert_less ( var, (_side_fe_var[dim].size() ) );
1234  return cast_ptr<FEBase *>( (_side_fe_var[dim][var] ) );
1235 }
1236 
1237 template<typename OutputShape>
1238 inline
1239 void FEMContext::get_edge_fe( unsigned int var, FEGenericBase<OutputShape> *& fe ) const
1240 {
1241  libmesh_assert_less ( var, _edge_fe_var.size() );
1242  fe = cast_ptr<FEGenericBase<OutputShape> *>( _edge_fe_var[var] );
1243 }
1244 
1245 inline
1246 FEBase * FEMContext::get_edge_fe( unsigned int var ) const
1247 {
1248  libmesh_assert_less ( var, _edge_fe_var.size() );
1249  return cast_ptr<FEBase *>( _edge_fe_var[var] );
1250 }
1251 
1252 
1253 } // namespace libMesh
1254 
1255 #endif // LIBMESH_FEM_CONTEXT_H
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
Number fixed_interior_value(unsigned int var, unsigned int qp) const
Definition: fem_context.C:973
unsigned char get_side() const
Accessor for current side of Elem object.
Definition: fem_context.h:885
const QBase & get_edge_qrule() const
Accessor for element edge quadrature rule.
Definition: fem_context.h:794
void interior_curl(unsigned int var, unsigned int qp, OutputType &curl_u) const
Definition: fem_context.C:512
void set_mesh_z_var(unsigned int z_var)
Accessor for z-variable of moving mesh System.
Definition: fem_context.h:859
FEBase * get_element_fe(unsigned int var) const
Accessor for interior finite element object for scalar-valued variable var for the largest dimension ...
Definition: fem_context.h:273
void elem_position_set(Real theta)
Uses the coordinate data specified by mesh_*_position configuration to set the geometry of elem to th...
Definition: fem_context.h:1195
This class provides all data required for a physics package (e.g.
Definition: diff_context.h:54
virtual void pre_fe_reinit(const System &, const Elem *e)
Reinitializes local data vectors/matrices on the current geometric element.
Definition: fem_context.C:1582
unsigned int get_mesh_x_var() const
Accessor for x-variable of moving mesh System.
Definition: fem_context.h:823
void set_mesh_x_var(unsigned int x_var)
Accessor for x-variable of moving mesh System.
Definition: fem_context.h:831
const NumericVector< Number > * _custom_solution
Data with which to do algebra reinitialization.
Definition: fem_context.h:994
Tensor interior_hessian(unsigned int var, unsigned int qp) const
Definition: fem_context.C:446
void get_edge_fe(unsigned int var, FEGenericBase< OutputShape > *&fe) const
Accessor for edge (3D only!) finite element object for variable var.
Definition: fem_context.h:1239
void set_elem(const Elem *e)
Helper function to promote accessor usage.
Definition: fem_context.C:1788
unsigned int _mesh_z_var
Definition: fem_context.h:972
void side_values(unsigned int var, const NumericVector< Number > &_system_vector, std::vector< OutputType > &side_values_vector) const
Fills a vector of values of the _system_vector at the all the quadrature points on the current elemen...
Definition: fem_context.C:598
void point_curl(unsigned int var, const Point &p, OutputType &curl_u, const Real tolerance=TOLERANCE) const
Definition: fem_context.C:937
Number point_value(unsigned int var, const Point &p) const
Definition: fem_context.C:788
unsigned int dim
Tensor fixed_point_hessian(unsigned int var, const Point &p) const
Definition: fem_context.C:1220
void side_hessians(unsigned int var, const NumericVector< Number > &_system_vector, std::vector< OutputType > &d2u_vals) const
Fills a vector of hessians of the _system_vector at the all the quadrature points on the current elem...
Definition: fem_context.C:745
ImplicitSystem & sys
bool has_elem() const
Test for current Elem object.
Definition: fem_context.h:865
std::vector< boundary_id_type > side_boundary_ids() const
Lists the boundary ids found on the current side.
Definition: fem_context.C:216
TensorTools::MakeReal< OutputType >::type value_shape
Definition: fem_context.h:1042
void interior_hessians(unsigned int var, const NumericVector< Number > &_system_vector, std::vector< OutputType > &d2u_vals) const
Fills a vector of hessians of the _system_vector at the all the quadrature points in the current elem...
Definition: fem_context.C:470
const QBase & get_side_qrule() const
Accessor for element side quadrature rule for the dimension of the current _elem. ...
Definition: fem_context.h:772
UniquePtr< FEGenericBase< Real > > _real_fe
Definition: fem_context.h:1004
FEGenericBase< value_shape > value_base
Definition: fem_context.h:1043
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
void some_value(unsigned int var, unsigned int qp, OutputType &u) const
Helper function to reduce some code duplication in the *interior_value methods.
Definition: fem_context.C:234
void point_accel(unsigned int var, const Point &p, OutputType &u) const
const class libmesh_nullptr_t libmesh_nullptr
std::vector< QBase * > _side_qrule
Quadrature rules for element sides The FEM context will try to find a quadrature rule that correctly ...
Definition: fem_context.h:1149
void get_element_fe(unsigned int var, FEGenericBase< OutputShape > *&fe) const
Accessor for interior finite element object for variable var for the largest dimension in the mesh...
Definition: fem_context.h:262
void set_algebraic_type(const AlgebraicType atype)
Setting which determines whether to initialize algebraic structures (elem_*) on each element and set ...
Definition: fem_context.h:946
int _extra_quadrature_order
The extra quadrature order for this context.
Definition: fem_context.h:1163
static const Real TOLERANCE
void elem_position_get()
Uses the geometry of elem to set the coordinate data specified by mesh_*_position configuration...
Definition: fem_context.C:1441
The libMesh namespace provides an interface to certain functionality in the library.
FEGenericBase< hess_shape > hess_base
Definition: fem_context.h:1053
void some_hessian(unsigned int var, unsigned int qp, OutputType &u) const
Helper function to reduce some code duplication in the *interior_hessian methods. ...
Definition: fem_context.C:296
virtual ~FEMContext()
Destructor.
Definition: fem_context.C:175
const System * get_mesh_system() const
Accessor for moving mesh System.
Definition: fem_context.h:811
Number fixed_side_value(unsigned int var, unsigned int qp) const
Definition: fem_context.C:1049
const Elem & get_elem() const
Accessor for current Elem object.
Definition: fem_context.h:871
Number interior_value(unsigned int var, unsigned int qp) const
Definition: fem_context.C:326
void side_gradients(unsigned int var, const NumericVector< Number > &_system_vector, std::vector< OutputType > &side_gradients_vector) const
Fills a vector with the gradient of the solution variable var at all the quadrature points on the cur...
Definition: fem_context.C:678
const QBase & get_side_qrule(unsigned char dim) const
Accessor for element side quadrature rule.
Definition: fem_context.h:785
void interior_rate(unsigned int var, unsigned int qp, OutputType &u) const
Definition: fem_context.C:1273
libmesh_assert(j)
virtual void nonlocal_reinit(Real theta) libmesh_override
Gives derived classes the opportunity to reinitialize data needed for nonlocal calculations at a new ...
Definition: fem_context.C:1372
virtual void elem_side_reinit(Real theta) libmesh_override
Resets the current time in the context.
Definition: fem_context.C:1342
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
Defines a dense subvector for use in finite element computations.
Gradient fixed_point_gradient(unsigned int var, const Point &p) const
Definition: fem_context.C:1169
UniquePtr< FEGenericBase< RealGradient > > _real_grad_fe
Definition: fem_context.h:1005
unsigned char get_dim() const
Accessor for cached mesh dimension.
Definition: fem_context.h:899
virtual void set_mesh_system(System *sys)
Tells the FEMContext that system sys contains the isoparametric Lagrangian variables which correspond...
Definition: fem_context.h:805
unsigned char edge
Current edge for edge_* to examine.
Definition: fem_context.h:982
const std::set< unsigned char > & elem_dimensions() const
Definition: fem_context.h:913
std::set< unsigned char > _elem_dims
Cached dimensions of elements in the mesh, plus dimension 0 if SCALAR variables are in use...
Definition: fem_context.h:1133
void(FEMContext::* value_getter)(unsigned int, value_base *&, unsigned char) const
Definition: fem_context.h:1044
int8_t boundary_id_type
Definition: id_types.h:51
unsigned char _elem_dim
Cached dimension of this->_elem.
Definition: fem_context.h:1127
virtual void side_fe_reinit()
Reinitializes side FE objects on the current geometric element.
Definition: fem_context.C:1405
void _update_time_from_system(Real theta)
Update the time in the context object for the given value of theta, based on the values of "time" and...
Definition: fem_context.C:1796
void init_internal_data(const System &sys)
Helper function used in constructors to set up internal data.
Definition: fem_context.C:79
FEGenericBase< OutputShape > * cached_fe(const unsigned int elem_dim, const FEType fe_type) const
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
System * _mesh_sys
System from which to acquire moving mesh information.
Definition: fem_context.h:967
Tensor fixed_side_hessian(unsigned int var, unsigned int qp) const
Definition: fem_context.C:1098
void interior_accel(unsigned int var, unsigned int qp, OutputType &u) const
Definition: fem_context.C:1295
unsigned int _mesh_x_var
Variables from which to acquire moving mesh information.
Definition: fem_context.h:972
The BoundaryInfo class contains information relevant to boundary conditions including storing faces...
Definition: boundary_info.h:56
TensorTools::MakeReal< Rank1Decrement >::type grad_shape
Definition: fem_context.h:1047
This class provides all data required for a physics package (e.g.
Definition: fem_context.h:61
TensorTools::MakeReal< Rank2Decrement >::type hess_shape
Definition: fem_context.h:1052
AlgebraicType _atype
Keep track of what type of algebra reinitialization is to be done.
Definition: fem_context.h:989
Tensor point_hessian(unsigned int var, const Point &p) const
Definition: fem_context.C:886
std::map< FEType, FEAbstract * > _edge_fe
Definition: fem_context.h:1095
const BoundaryInfo & _boundary_info
Saved reference to BoundaryInfo on the mesh for this System.
Definition: fem_context.h:1112
std::vector< std::vector< FEAbstract * > > _element_fe_var
Pointers to the same finite element objects, but indexed by variable number.
Definition: fem_context.h:1104
FEGenericBase< grad_shape > grad_base
Definition: fem_context.h:1048
Gradient interior_gradient(unsigned int var, unsigned int qp) const
Definition: fem_context.C:381
virtual void elem_fe_reinit(const std::vector< Point > *const pts=libmesh_nullptr)
Reinitializes interior FE objects on the current geometric element.
Definition: fem_context.C:1382
FEGenericBase< Real > FEBase
const QBase & get_element_qrule(unsigned char dim) const
Accessor for element interior quadrature rule.
Definition: fem_context.h:778
TensorTools::DecrementRank< Rank1Decrement >::type Rank2Decrement
Definition: fem_context.h:1039
unsigned char get_edge() const
Accessor for current edge of Elem object.
Definition: fem_context.h:891
void interior_gradients(unsigned int var, const NumericVector< Number > &_system_vector, std::vector< OutputType > &interior_gradients_vector) const
Fills a vector with the gradient of the solution variable var at all the quadrature points in the cur...
Definition: fem_context.C:408
unsigned char _dim
Cached dimension of largest dimension element in this mesh.
Definition: fem_context.h:1122
virtual void edge_fe_reinit()
Reinitializes edge FE objects on the current geometric element.
Definition: fem_context.C:1425
std::vector< std::map< FEType, FEAbstract * > > _side_fe
Definition: fem_context.h:1094
std::vector< FEAbstract * > _edge_fe_var
Definition: fem_context.h:1106
Elem & get_elem()
Accessor for current Elem object.
Definition: fem_context.h:878
AlgebraicType algebraic_type() const
Definition: fem_context.h:952
void side_rate(unsigned int var, unsigned int qp, OutputType &u) const
Definition: fem_context.C:1285
void(FEMContext::* hess_getter)(unsigned int, hess_base *&, unsigned char) const
Definition: fem_context.h:1054
std::vector< std::map< FEType, FEAbstract * > > _element_fe
Finite element objects for each variable&#39;s interior, sides and edges.
Definition: fem_context.h:1093
bool has_side_boundary_id(boundary_id_type id) const
Reports if the boundary id is found on the current side.
Definition: fem_context.C:209
std::vector< QBase * > _element_qrule
Quadrature rule for element interior.
Definition: fem_context.h:1141
virtual void elem_reinit(Real theta) libmesh_override
Resets the current time in the context.
Definition: fem_context.C:1318
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
unsigned char side
Current side for side_* to examine.
Definition: fem_context.h:977
Number side_value(unsigned int var, unsigned int qp) const
Definition: fem_context.C:575
AlgebraicType
Enum describing what data to use when initializing algebraic structures on each element.
Definition: fem_context.h:933
TensorTools::DecrementRank< OutputType >::type Rank1Decrement
Definition: fem_context.h:1038
const QBase & get_element_qrule() const
Accessor for element interior quadrature rule for the dimension of the current _elem.
Definition: fem_context.h:765
Gradient side_gradient(unsigned int var, unsigned int qp) const
Definition: fem_context.C:633
void set_mesh_y_var(unsigned int y_var)
Accessor for y-variable of moving mesh System.
Definition: fem_context.h:845
Gradient point_gradient(unsigned int var, const Point &p) const
Definition: fem_context.C:834
FEBase * get_side_fe(unsigned int var) const
Accessor for side finite element object for scalar-valued variable var for the largest dimension in t...
Definition: fem_context.h:310
FEMContext(const System &sys)
Constructor.
Definition: fem_context.C:37
Helper nested class for C++03-compatible "template typedef".
Definition: fem_context.h:1035
void some_gradient(unsigned int var, unsigned int qp, OutputType &u) const
Helper function to reduce some code duplication in the *interior_gradient methods.
Definition: fem_context.C:263
unsigned int _mesh_y_var
Definition: fem_context.h:972
Tensor side_hessian(unsigned int var, unsigned int qp) const
Definition: fem_context.C:716
void side_accel(unsigned int var, unsigned int qp, OutputType &u) const
Definition: fem_context.C:1307
FEGenericBase< OutputShape > * build_new_fe(const FEGenericBase< OutputShape > *fe, const Point &p, const Real tolerance=TOLERANCE) const
Helper function to reduce some code duplication in the *_point_* methods.
Definition: fem_context.C:1884
unsigned int get_mesh_z_var() const
Accessor for z-variable of moving mesh System.
Definition: fem_context.h:851
void set_custom_solution(const NumericVector< Number > *custom_sol)
Set a NumericVector to be used in place of current_local_solution for calculating elem_solution...
Definition: fem_context.h:961
void point_rate(unsigned int var, const Point &p, OutputType &u) const
std::vector< std::vector< FEAbstract * > > _side_fe_var
Definition: fem_context.h:1105
const DenseSubVector< Number > &(DiffContext::* diff_subsolution_getter)(unsigned int) const
Helper typedef to simplify refactoring.
Definition: fem_context.h:1028
virtual void elem_edge_reinit(Real theta) libmesh_override
Resets the current time in the context.
Definition: fem_context.C:1357
Gradient fixed_interior_gradient(unsigned int var, unsigned int qp) const
Definition: fem_context.C:996
const Elem * _elem
Current element for element_* to examine.
Definition: fem_context.h:1117
void interior_values(unsigned int var, const NumericVector< Number > &_system_vector, std::vector< OutputType > &interior_values_vector) const
Fills a vector of values of the _system_vector at the all the quadrature points in the current elemen...
Definition: fem_context.C:346
void _do_elem_position_set(Real theta)
Uses the coordinate data specified by mesh_*_position configuration to set the geometry of elem to th...
Definition: fem_context.C:1506
void interior_div(unsigned int var, unsigned int qp, OutputType &div_u) const
Definition: fem_context.C:543
void get_side_fe(unsigned int var, FEGenericBase< OutputShape > *&fe) const
Accessor for edge/face (2D/3D) finite element object for variable var for the largest dimension in th...
Definition: fem_context.h:299
Number fixed_point_value(unsigned int var, const Point &p) const
Definition: fem_context.C:1123
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
System * get_mesh_system()
Accessor for moving mesh System.
Definition: fem_context.h:817
Gradient fixed_side_gradient(unsigned int var, unsigned int qp) const
Definition: fem_context.C:1073
Tensor fixed_interior_hessian(unsigned int var, unsigned int qp) const
Definition: fem_context.C:1023
UniquePtr< QBase > _edge_qrule
Quadrature rules for element edges.
Definition: fem_context.h:1158
The QBase class provides the basic functionality from which various quadrature rules can be derived...
Definition: quadrature.h:53
void ErrorVector unsigned int
Definition: adjoints_ex3.C:360
This class forms the foundation from which generic finite elements may be derived.
unsigned char get_elem_dim() const
Definition: fem_context.h:906
unsigned int get_mesh_y_var() const
Accessor for y-variable of moving mesh System.
Definition: fem_context.h:837
void(FEMContext::* grad_getter)(unsigned int, grad_base *&, unsigned char) const
Definition: fem_context.h:1049
This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space.