libMesh
fe_base.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2025 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_FE_BASE_H
21 #define LIBMESH_FE_BASE_H
22 
23 // Local includes
24 #include "libmesh/libmesh_common.h"
25 #include "libmesh/compare_types.h"
26 #include "libmesh/fe_abstract.h"
27 #include "libmesh/fe_transformation_base.h"
28 #include "libmesh/point.h"
29 #include "libmesh/reference_counted_object.h"
30 #include "libmesh/tensor_tools.h"
31 #include "libmesh/type_n_tensor.h"
32 #include "libmesh/vector_value.h"
33 #include "libmesh/dense_matrix.h"
34 
35 // C++ includes
36 #include <cstddef>
37 #include <vector>
38 #include <memory>
39 
40 namespace libMesh
41 {
42 
43 
44 // forward declarations
45 template <typename T> class DenseMatrix;
46 template <typename T> class DenseVector;
47 class BoundaryInfo;
48 class DofConstraints;
49 class DofMap;
50 class Elem;
51 class MeshBase;
52 template <typename T> class NumericVector;
53 class QBase;
54 template <typename T> class FETransformationBase;
55 class FEType;
56 
57 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
58 class NodeConstraints;
59 #endif
60 
61 #ifdef LIBMESH_ENABLE_PERIODIC
62 class PeriodicBoundaries;
63 class PointLocatorBase;
64 #endif
65 
66 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
67 template <unsigned int Dim, FEFamily T_radial, InfMapType T_map>
68 class InfFE;
69 #endif
70 
84 template <typename OutputType>
85 class FEGenericBase : public FEAbstract
86 {
87 protected:
88 
94  FEGenericBase (const unsigned int dim,
95  const FEType & fet);
96 
97 public:
98 
102  virtual ~FEGenericBase();
103 
112  static std::unique_ptr<FEGenericBase> build (const unsigned int dim,
113  const FEType & type);
114 
119  typedef OutputType OutputShape;
127 
128 
129 
130 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
131 
140  static std::unique_ptr<FEGenericBase> build_InfFE (const unsigned int dim,
141  const FEType & type);
142 
143 #endif
144 
145 #ifdef LIBMESH_ENABLE_AMR
146 
153  static void compute_proj_constraints (DofConstraints & constraints,
154  DofMap & dof_map,
155  const unsigned int variable_number,
156  const Elem * elem);
157 
165  static void coarsened_dof_values(const NumericVector<Number> & global_vector,
166  const DofMap & dof_map,
167  const Elem * coarse_elem,
168  DenseVector<Number> & coarse_dofs,
169  const unsigned int var,
170  const bool use_old_dof_indices = false);
171 
178  static void coarsened_dof_values(const NumericVector<Number> & global_vector,
179  const DofMap & dof_map,
180  const Elem * coarse_elem,
181  DenseVector<Number> & coarse_dofs,
182  const bool use_old_dof_indices = false);
183 
184 #endif // #ifdef LIBMESH_ENABLE_AMR
185 
186 #ifdef LIBMESH_ENABLE_PERIODIC
187 
193  static void compute_periodic_constraints (DofConstraints & constraints,
194  DofMap & dof_map,
195  const PeriodicBoundaries & boundaries,
196  const MeshBase & mesh,
197  const PointLocatorBase * point_locator,
198  const unsigned int variable_number,
199  const Elem * elem);
200 
201 #endif // LIBMESH_ENABLE_PERIODIC
202 
207  const std::vector<std::vector<OutputShape>> & get_phi() const
209  calculate_phi = true; return phi; }
210 
211  const std::vector<std::vector<OutputShape>> & get_dual_phi() const
212  {
214  calculate_dual = true;
215  // Dual phi computation relies on primal phi computation
216  this->request_phi();
217  return dual_phi;
218  }
219 
220  virtual void request_phi() const override
221  { get_phi(); }
222 
223  virtual void request_dual_phi() const override
224  { get_dual_phi(); }
225 
230  const std::vector<std::vector<OutputGradient>> & get_dphi() const
232  calculate_dphi = calculate_dphiref = true; return dphi; }
233 
234  const std::vector<std::vector<OutputGradient>> & get_dual_dphi() const
237 
238  virtual void request_dphi() const override
239  { get_dphi(); }
240 
241  virtual void request_dual_dphi() const override
242  { get_dual_dphi(); }
243 
245  { return dual_coeff; }
246 
251  virtual_for_inffe
252  const std::vector<std::vector<OutputShape>> & get_curl_phi() const
254  calculate_curl_phi = calculate_dphiref = true; return curl_phi; }
255 
260  virtual_for_inffe
261  const std::vector<std::vector<OutputDivergence>> & get_div_phi() const
263  calculate_div_phi = calculate_dphiref = true; return div_phi; }
264 
269  const std::vector<std::vector<OutputShape>> & get_dphidx() const
271  calculate_dphi = calculate_dphiref = true; return dphidx; }
272 
277  const std::vector<std::vector<OutputShape>> & get_dphidy() const
279  calculate_dphi = calculate_dphiref = true; return dphidy; }
280 
285  const std::vector<std::vector<OutputShape>> & get_dphidz() const
287  calculate_dphi = calculate_dphiref = true; return dphidz; }
288 
293  const std::vector<std::vector<OutputShape>> & get_dphidxi() const
295  calculate_dphiref = true; return dphidxi; }
296 
301  const std::vector<std::vector<OutputShape>> & get_dphideta() const
303  calculate_dphiref = true; return dphideta; }
304 
309  const std::vector<std::vector<OutputShape>> & get_dphidzeta() const
311  calculate_dphiref = true; return dphidzeta; }
312 
313 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
314 
319  const std::vector<std::vector<OutputTensor>> & get_d2phi() const
321  calculate_d2phi = calculate_dphiref = true; return d2phi; }
322 
323  const std::vector<std::vector<OutputTensor>> & get_dual_d2phi() const
326 
331  const std::vector<std::vector<OutputShape>> & get_d2phidx2() const
333  calculate_d2phi = calculate_dphiref = true; return d2phidx2; }
334 
339  const std::vector<std::vector<OutputShape>> & get_d2phidxdy() const
341  calculate_d2phi = calculate_dphiref = true; return d2phidxdy; }
342 
347  const std::vector<std::vector<OutputShape>> & get_d2phidxdz() const
349  calculate_d2phi = calculate_dphiref = true; return d2phidxdz; }
350 
355  const std::vector<std::vector<OutputShape>> & get_d2phidy2() const
357  calculate_d2phi = calculate_dphiref = true; return d2phidy2; }
358 
363  const std::vector<std::vector<OutputShape>> & get_d2phidydz() const
365  calculate_d2phi = calculate_dphiref = true; return d2phidydz; }
366 
371  const std::vector<std::vector<OutputShape>> & get_d2phidz2() const
373  calculate_d2phi = calculate_dphiref = true; return d2phidz2; }
374 
379  const std::vector<std::vector<OutputShape>> & get_d2phidxi2() const
381  calculate_d2phi = calculate_dphiref = true; return d2phidxi2; }
382 
387  const std::vector<std::vector<OutputShape>> & get_d2phidxideta() const
390 
395  const std::vector<std::vector<OutputShape>> & get_d2phidxidzeta() const
398 
403  const std::vector<std::vector<OutputShape>> & get_d2phideta2() const
405  calculate_d2phi = calculate_dphiref = true; return d2phideta2; }
406 
411  const std::vector<std::vector<OutputShape>> & get_d2phidetadzeta() const
414 
419  const std::vector<std::vector<OutputShape>> & get_d2phidzeta2() const
421  calculate_d2phi = calculate_dphiref = true; return d2phidzeta2; }
422 
423 #endif //LIBMESH_ENABLE_SECOND_DERIVATIVES
424 
425 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
426 
437  const std::vector<OutputGradient> & get_dphase() const
438  { return dphase; }
439 
440 
453  virtual const std::vector<Real> & get_Sobolev_weight() const
454  { return weight; }
455 
461  virtual const std::vector<RealGradient> & get_Sobolev_dweight() const
462  { return dweight; }
463 
470  virtual const std::vector<Real> & get_Sobolev_weightxR_sq() const
471  { return weight; }
472 
480  virtual const std::vector<RealGradient> & get_Sobolev_dweightxR_sq() const
481  { return dweight; }
482 
493  virtual const std::vector<std::vector<OutputShape>> & get_phi_over_decayxR () const
494  { return get_phi();}
495 
501  virtual const std::vector<std::vector<OutputGradient>> & get_dphi_over_decayxR () const
502  { return get_dphi();}
503 
511  virtual const std::vector<std::vector<OutputGradient>> & get_dphi_over_decay () const
512  { return get_dphi();}
513 
514 #endif
515 
519  virtual void print_phi(std::ostream & os) const override;
520  virtual void print_dual_phi(std::ostream & os) const override;
521 
526  virtual void print_dphi(std::ostream & os) const override;
527  virtual void print_dual_dphi(std::ostream & os) const override;
528 
529 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
530 
535  virtual void print_d2phi(std::ostream & os) const override;
536  virtual void print_dual_d2phi(std::ostream & os) const override;
537 
538 #endif
539 
540 
541 protected:
542 
543 
544 
545 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
546 
552  virtual void init_base_shape_functions(const std::vector<Point> & qp,
553  const Elem * e) = 0;
554 
555 #endif
556 
561  virtual_for_inffe
562  void determine_calculations();
563 
568  bool calculating_nothing() const
569  {
570  return calculate_nothing &&
571  !this->calculate_phi && !this->calculate_dphi &&
572 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
573  !this->calculate_d2phi &&
574 #endif
575  !this->calculate_curl_phi && !this->calculate_div_phi &&
576  !this->calculate_map;
577  }
578 
589  virtual void compute_shape_functions(const Elem * elem, const std::vector<Point> & qp) override;
590 
596  void compute_dual_shape_coeffs(const std::vector<Real> & JxW, const std::vector<std::vector<OutputShape>> & phi);
597 
604 
609  std::unique_ptr<FETransformationBase<OutputType>> _fe_trans;
610 
614  std::vector<std::vector<OutputShape>> phi;
615  std::vector<std::vector<OutputShape>> dual_phi;
616 
620  std::vector<std::vector<OutputGradient>> dphi;
621  std::vector<std::vector<OutputGradient>> dual_dphi;
622 
627 
631  std::vector<std::vector<OutputShape>> curl_phi;
632 
636  std::vector<std::vector<OutputDivergence>> div_phi;
637 
641  std::vector<std::vector<OutputShape>> dphidxi;
642 
646  std::vector<std::vector<OutputShape>> dphideta;
647 
651  std::vector<std::vector<OutputShape>> dphidzeta;
652 
656  std::vector<std::vector<OutputShape>> dphidx;
657 
661  std::vector<std::vector<OutputShape>> dphidy;
662 
666  std::vector<std::vector<OutputShape>> dphidz;
667 
668 
669 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
670 
674  std::vector<std::vector<OutputTensor>> d2phi;
675  std::vector<std::vector<OutputTensor>> dual_d2phi;
676 
680  std::vector<std::vector<OutputShape>> d2phidxi2;
681 
685  std::vector<std::vector<OutputShape>> d2phidxideta;
686 
690  std::vector<std::vector<OutputShape>> d2phidxidzeta;
691 
695  std::vector<std::vector<OutputShape>> d2phideta2;
696 
700  std::vector<std::vector<OutputShape>> d2phidetadzeta;
701 
705  std::vector<std::vector<OutputShape>> d2phidzeta2;
706 
710  std::vector<std::vector<OutputShape>> d2phidx2;
711 
715  std::vector<std::vector<OutputShape>> d2phidxdy;
716 
720  std::vector<std::vector<OutputShape>> d2phidxdz;
721 
725  std::vector<std::vector<OutputShape>> d2phidy2;
726 
730  std::vector<std::vector<OutputShape>> d2phidydz;
731 
735  std::vector<std::vector<OutputShape>> d2phidz2;
736 
737 #endif
738 
739 
740 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
741 
742  //--------------------------------------------------------------
743  /* protected members for infinite elements, which are accessed
744  * from the outside through some inline functions
745  */
746 
747 
753  std::vector<OutputGradient> dphase;
754 
760  std::vector<RealGradient> dweight;
761 
767  std::vector<Real> weight;
768 
769 #endif
770 
771 private:
772 
773 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
774 
780  template <unsigned int friend_Dim, FEFamily friend_T_radial, InfMapType friend_T_map>
781  friend class InfFE;
782 
783 #endif
784 
785 
786 };
787 
788 // --------------------------------------------------------------------
789 // Generic templates. We specialize for OutputType = Real, so these are
790 // only used for OutputType = RealVectorValue
791 template <typename OutputType>
793 {
794  libmesh_error_msg(
795  "Computation of dual shape functions for vector finite element "
796  "families is not currently implemented");
797 }
798 
799 template <typename OutputType>
800 void FEGenericBase<OutputType>::compute_dual_shape_coeffs(const std::vector<Real> & /*JxW*/, const std::vector<std::vector<OutputShape>> & /*phi_vals*/)
801 {
802  libmesh_error_msg(
803  "Computation of dual shape functions for vector finite element "
804  "families is not currently implemented");
805 }
806 
807 // -----------------------------------------------------------
808 // Forward declaration of specialization
809 template <>
811 
812 template <>
813 void FEGenericBase<Real>::compute_dual_shape_coeffs(const std::vector<Real> & /*JxW*/, const std::vector<std::vector<OutputShape>> & /*phi_vals*/);
814 
815 
816 // Typedefs for convenience and backwards compatibility
819 
820 
821 
822 
823 // ------------------------------------------------------------
824 // FEGenericBase class inline members
825 template <typename OutputType>
826 inline
828  const FEType & fet) :
829  FEAbstract(d,fet),
830  _fe_trans( FETransformationBase<OutputType>::build(fet) ),
831  phi(),
832  dual_phi(),
833  dphi(),
834  dual_dphi(),
835  curl_phi(),
836  div_phi(),
837  dphidxi(),
838  dphideta(),
839  dphidzeta(),
840  dphidx(),
841  dphidy(),
842  dphidz()
843 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
844  ,d2phi(),
845  dual_d2phi(),
846  d2phidxi2(),
847  d2phidxideta(),
848  d2phidxidzeta(),
849  d2phideta2(),
850  d2phidetadzeta(),
851  d2phidzeta2(),
852  d2phidx2(),
853  d2phidxdy(),
854  d2phidxdz(),
855  d2phidy2(),
856  d2phidydz(),
857  d2phidz2()
858 #endif
859 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
860  ,dphase(),
861  dweight(),
862  weight()
863 #endif
864 {
865 }
866 
867 
868 
869 template <typename OutputType>
870 inline
872 {
873 }
874 
875 } // namespace libMesh
876 
877 #endif // LIBMESH_FE_BASE_H
FEGenericBase(const unsigned int dim, const FEType &fet)
Constructor.
Definition: fe_base.h:827
bool calculate_d2phi
Should we calculate shape function hessians?
Definition: fe_abstract.h:702
bool calculate_curl_phi
Should we calculate shape function curls?
Definition: fe_abstract.h:712
std::vector< std::vector< OutputTensor > > d2phi
Shape function second derivative values.
Definition: fe_base.h:674
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:196
std::vector< std::vector< OutputShape > > dphidxi
Shape function derivatives in the xi direction.
Definition: fe_base.h:641
std::vector< std::vector< OutputShape > > d2phidxdz
Shape function second derivatives in the x-z direction.
Definition: fe_base.h:720
const std::vector< std::vector< OutputShape > > & get_d2phidxdy() const
Definition: fe_base.h:339
virtual_for_inffe const std::vector< std::vector< OutputDivergence > > & get_div_phi() const
Definition: fe_base.h:261
std::vector< std::vector< OutputShape > > dphidzeta
Shape function derivatives in the zeta direction.
Definition: fe_base.h:651
bool calculations_started
Have calculations with this object already been started? Then all get_* functions should already have...
Definition: fe_abstract.h:666
bool calculating_nothing() const
Definition: fe_base.h:568
std::vector< std::vector< OutputShape > > d2phidydz
Shape function second derivatives in the y-z direction.
Definition: fe_base.h:730
bool calculate_phi
Should we calculate shape functions?
Definition: fe_abstract.h:691
A specific instantiation of the FEBase class.
Definition: fe.h:42
const std::vector< std::vector< OutputShape > > & get_dual_phi() const
Definition: fe_base.h:211
std::unique_ptr< FETransformationBase< OutputType > > _fe_trans
Object that handles computing shape function values, gradients, etc in the physical domain...
Definition: fe_base.h:609
TensorTools::DecrementRank< OutputNumber >::type OutputNumberDivergence
Definition: fe_base.h:126
virtual const std::vector< Real > & get_Sobolev_weightxR_sq() const
Definition: fe_base.h:470
virtual void print_dphi(std::ostream &os) const override
Prints the value of each shape function&#39;s derivative at each quadrature point.
Definition: fe_base.C:895
void compute_dual_shape_functions()
Compute dual_phi, dual_dphi, dual_d2phi It is only valid for this to be called after reinit has occur...
Definition: fe_base.h:792
virtual void request_dphi() const override
request dphi calculations
Definition: fe_base.h:238
const std::vector< std::vector< OutputShape > > & get_d2phideta2() const
Definition: fe_base.h:403
virtual void init_base_shape_functions(const std::vector< Point > &qp, const Elem *e)=0
Initialize the data fields for the base of an an infinite element.
We&#39;re using a class instead of a typedef to allow forward declarations and future flexibility...
const std::vector< std::vector< OutputShape > > & get_d2phidxideta() const
Definition: fe_base.h:387
std::vector< std::vector< OutputShape > > d2phidxideta
Shape function second derivatives in the xi-eta direction.
Definition: fe_base.h:685
std::vector< std::vector< OutputTensor > > dual_d2phi
Definition: fe_base.h:675
const std::vector< std::vector< OutputShape > > & get_d2phidydz() const
Definition: fe_base.h:363
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
TensorTools::IncrementRank< OutputNumberGradient >::type OutputNumberTensor
Definition: fe_base.h:125
virtual const std::vector< Real > & get_Sobolev_weight() const
Definition: fe_base.h:453
MeshBase & mesh
virtual void request_phi() const override
request phi calculations
Definition: fe_base.h:220
std::vector< Real > weight
Used for certain infinite element families: the additional radial weight in local coordinates...
Definition: fe_base.h:767
std::vector< std::vector< OutputGradient > > dual_dphi
Definition: fe_base.h:621
TensorTools::IncrementRank< OutputNumber >::type OutputNumberGradient
Definition: fe_base.h:124
The Node constraint storage format.
Definition: dof_map.h:156
This class defines a vector in LIBMESH_DIM dimensional Real or Complex space.
TensorTools::IncrementRank< OutputShape >::type OutputGradient
Definition: fe_base.h:120
std::vector< std::vector< OutputShape > > d2phidx2
Shape function second derivatives in the x direction.
Definition: fe_base.h:710
std::vector< std::vector< OutputShape > > curl_phi
Shape function curl values.
Definition: fe_base.h:631
The libMesh namespace provides an interface to certain functionality in the library.
const std::vector< std::vector< OutputShape > > & get_dphideta() const
Definition: fe_base.h:301
DenseMatrix< Real > dual_coeff
Coefficient matrix for the dual basis.
Definition: fe_base.h:626
std::vector< std::vector< OutputShape > > dphidy
Shape function derivatives in the y direction.
Definition: fe_base.h:661
FEGenericBase< RealGradient > FEVectorBase
Definition: fe_base.h:818
TensorTools::DecrementRank< OutputShape >::type OutputDivergence
Definition: fe_base.h:122
std::vector< std::vector< OutputShape > > d2phidy2
Shape function second derivatives in the y direction.
Definition: fe_base.h:725
virtual void request_dual_phi() const override
Definition: fe_base.h:223
bool calculate_div_phi
Should we calculate shape function divergences?
Definition: fe_abstract.h:717
This is the MeshBase class.
Definition: mesh_base.h:75
const std::vector< std::vector< OutputShape > > & get_dphidzeta() const
Definition: fe_base.h:309
const std::vector< OutputGradient > & get_dphase() const
Definition: fe_base.h:437
const DenseMatrix< Real > & get_dual_coeff() const
Definition: fe_base.h:244
std::vector< std::vector< OutputShape > > d2phidetadzeta
Shape function second derivatives in the eta-zeta direction.
Definition: fe_base.h:700
const std::vector< std::vector< OutputShape > > & get_dphidx() const
Definition: fe_base.h:269
static void coarsened_dof_values(const NumericVector< Number > &global_vector, const DofMap &dof_map, const Elem *coarse_elem, DenseVector< Number > &coarse_dofs, const unsigned int var, const bool use_old_dof_indices=false)
Creates a local projection on coarse_elem, based on the DoF values in global_vector for it&#39;s children...
Definition: fe_base.C:1012
std::vector< std::vector< OutputShape > > d2phidxidzeta
Shape function second derivatives in the xi-zeta direction.
Definition: fe_base.h:690
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:179
std::vector< std::vector< OutputShape > > d2phidxdy
Shape function second derivatives in the x-y direction.
Definition: fe_base.h:715
std::vector< std::vector< OutputShape > > dphidx
Shape function derivatives in the x direction.
Definition: fe_base.h:656
const std::vector< std::vector< OutputGradient > > & get_dphi() const
Definition: fe_base.h:230
std::vector< std::vector< OutputShape > > phi
Shape function values.
Definition: fe_base.h:614
void compute_dual_shape_coeffs(const std::vector< Real > &JxW, const std::vector< std::vector< OutputShape >> &phi)
Compute the dual basis coefficients dual_coeff we rely on the JxW (or weights) and the phi values...
Definition: fe_base.h:800
std::vector< std::vector< OutputShape > > d2phideta2
Shape function second derivatives in the eta direction.
Definition: fe_base.h:695
const unsigned int dim
The dimensionality of the object.
Definition: fe_abstract.h:660
std::vector< std::vector< OutputShape > > dual_phi
Definition: fe_base.h:615
dof_id_type weight(const MeshBase &mesh, const processor_id_type pid)
Definition: mesh_tools.C:437
const std::vector< std::vector< OutputGradient > > & get_dual_dphi() const
Definition: fe_base.h:234
TensorTools::IncrementRank< OutputGradient >::type OutputTensor
Definition: fe_base.h:121
virtual const std::vector< std::vector< OutputShape > > & get_phi_over_decayxR() const
Definition: fe_base.h:493
const std::vector< std::vector< OutputShape > > & get_d2phidx2() const
Definition: fe_base.h:331
virtual void print_phi(std::ostream &os) const override
Prints the value of each shape function at each quadrature point.
Definition: fe_base.C:876
virtual const std::vector< std::vector< OutputGradient > > & get_dphi_over_decay() const
Definition: fe_base.h:511
std::vector< OutputGradient > dphase
Used for certain infinite element families: the first derivatives of the phase term in global coordin...
Definition: fe_base.h:753
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
libmesh_assert(ctx)
virtual void print_dual_dphi(std::ostream &os) const override
Definition: fe_base.C:903
virtual ~FEGenericBase()
Destructor.
const std::vector< std::vector< OutputShape > > & get_dphidy() const
Definition: fe_base.h:277
static void compute_proj_constraints(DofConstraints &constraints, DofMap &dof_map, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for non-conforming adapted meshes) corresponding to var...
Definition: fe_base.C:1570
const std::vector< std::vector< OutputShape > > & get_d2phidy2() const
Definition: fe_base.h:355
This is the base class for point locators.
const std::vector< std::vector< OutputTensor > > & get_d2phi() const
Definition: fe_base.h:319
OutputType OutputShape
Convenient typedefs for gradients of output, hessians of output, and potentially-complex-valued versi...
Definition: fe_base.h:119
TensorTools::MakeNumber< OutputShape >::type OutputNumber
Definition: fe_base.h:123
std::vector< std::vector< OutputDivergence > > div_phi
Shape function divergence values.
Definition: fe_base.h:636
const std::vector< std::vector< OutputShape > > & get_d2phidxi2() const
Definition: fe_base.h:379
FEGenericBase< Real > FEBase
virtual_for_inffe void determine_calculations()
Determine which values are to be calculated, for both the FE itself and for the FEMap.
Definition: fe_base.C:913
const std::vector< std::vector< OutputShape > > & get_d2phidxdz() const
Definition: fe_base.h:347
bool calculate_dphiref
Should we calculate reference shape function gradients?
Definition: fe_abstract.h:722
std::vector< std::vector< OutputGradient > > dphi
Shape function derivative values.
Definition: fe_base.h:620
virtual void request_dual_dphi() const override
Definition: fe_base.h:241
static void compute_periodic_constraints(DofConstraints &constraints, DofMap &dof_map, const PeriodicBoundaries &boundaries, const MeshBase &mesh, const PointLocatorBase *point_locator, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for meshes with periodic boundary conditions) correspon...
Definition: fe_base.C:1878
const std::vector< std::vector< OutputShape > > & get_d2phidxidzeta() const
Definition: fe_base.h:395
const std::vector< std::vector< OutputShape > > & get_d2phidetadzeta() const
Definition: fe_base.h:411
bool calculate_dphi
Should we calculate shape function gradients?
Definition: fe_abstract.h:696
bool calculate_map
Are we calculating mapping functions?
Definition: fe_abstract.h:686
virtual const std::vector< RealGradient > & get_Sobolev_dweight() const
Definition: fe_base.h:461
virtual void compute_shape_functions(const Elem *elem, const std::vector< Point > &qp) override
After having updated the jacobian and the transformation from local to global coordinates in FEMap::c...
Definition: fe_base.C:762
std::vector< std::vector< OutputShape > > d2phidz2
Shape function second derivatives in the z direction.
Definition: fe_base.h:735
virtual void print_dual_phi(std::ostream &os) const override
Definition: fe_base.C:884
const std::vector< std::vector< OutputShape > > & get_d2phidzeta2() const
Definition: fe_base.h:419
This class handles the computation of the shape functions in the physical domain. ...
Definition: fe_base.h:54
static std::unique_ptr< FEGenericBase > build_InfFE(const unsigned int dim, const FEType &type)
Builds a specific infinite element type.
virtual const std::vector< RealGradient > & get_Sobolev_dweightxR_sq() const
Definition: fe_base.h:480
std::vector< std::vector< OutputShape > > d2phidxi2
Shape function second derivatives in the xi direction.
Definition: fe_base.h:680
const std::vector< std::vector< OutputShape > > & get_dphidxi() const
Definition: fe_base.h:293
This class forms the foundation from which generic finite elements may be derived.
Definition: fe_abstract.h:99
virtual void print_d2phi(std::ostream &os) const override
Prints the value of each shape function&#39;s second derivatives at each quadrature point.
Definition: fe_base.C:989
bool calculate_dual
Are we calculating dual basis?
Definition: fe_abstract.h:671
virtual void print_dual_d2phi(std::ostream &os) const override
Definition: fe_base.C:997
The constraint matrix storage format.
Definition: dof_map.h:108
bool calculate_nothing
Are we potentially deliberately calculating nothing?
Definition: fe_abstract.h:681
This class forms the foundation from which generic finite elements may be derived.
std::vector< std::vector< OutputShape > > d2phidzeta2
Shape function second derivatives in the zeta direction.
Definition: fe_base.h:705
virtual const std::vector< std::vector< OutputGradient > > & get_dphi_over_decayxR() const
Definition: fe_base.h:501
std::vector< std::vector< OutputShape > > dphidz
Shape function derivatives in the z direction.
Definition: fe_base.h:666
const std::vector< std::vector< OutputShape > > & get_d2phidz2() const
Definition: fe_base.h:371
virtual_for_inffe const std::vector< std::vector< OutputShape > > & get_curl_phi() const
Definition: fe_base.h:252
const std::vector< std::vector< OutputShape > > & get_phi() const
Definition: fe_base.h:207
std::vector< std::vector< OutputShape > > dphideta
Shape function derivatives in the eta direction.
Definition: fe_base.h:646
std::vector< RealGradient > dweight
Used for certain infinite element families: the global derivative of the additional radial weight ...
Definition: fe_base.h:760
const std::vector< std::vector< OutputTensor > > & get_dual_d2phi() const
Definition: fe_base.h:323
const std::vector< std::vector< OutputShape > > & get_dphidz() const
Definition: fe_base.h:285