fe_base.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_FE_BASE_H
21 #define LIBMESH_FE_BASE_H
22 
23 // Local includes
24 #include "libmesh/libmesh_common.h"
25 #include "libmesh/auto_ptr.h"
26 #include "libmesh/compare_types.h"
27 #include "libmesh/enum_elem_type.h"
28 #include "libmesh/fe_abstract.h"
30 #include "libmesh/point.h"
32 #include "libmesh/tensor_tools.h"
33 #include "libmesh/type_n_tensor.h"
34 #include "libmesh/vector_value.h"
35 
36 // C++ includes
37 #include <cstddef>
38 #include <vector>
39 
40 
41 namespace libMesh
42 {
43 
44 
45 // forward declarations
46 template <typename T> class DenseMatrix;
47 template <typename T> class DenseVector;
48 class BoundaryInfo;
49 class DofConstraints;
50 class DofMap;
51 class Elem;
52 class MeshBase;
53 template <typename T> class NumericVector;
54 class QBase;
55 template <typename T> class FETransformationBase;
56 class FEType;
57 
58 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
59 class NodeConstraints;
60 #endif
61 
62 #ifdef LIBMESH_ENABLE_PERIODIC
63 class PeriodicBoundaries;
64 class PointLocatorBase;
65 #endif
66 
67 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
68 template <unsigned int Dim, FEFamily T_radial, InfMapType T_map>
69 class InfFE;
70 #endif
71 
97 template <typename OutputType>
98 class FEGenericBase : public FEAbstract
99 {
100 protected:
101 
107  FEGenericBase (const unsigned int dim,
108  const FEType & fet);
109 
110 public:
111 
115  virtual ~FEGenericBase();
116 
125  static UniquePtr<FEGenericBase> build (const unsigned int dim,
126  const FEType & type);
127 
132  typedef OutputType OutputShape;
140 
141 
142 
143 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
144 
153  static UniquePtr<FEGenericBase> build_InfFE (const unsigned int dim,
154  const FEType & type);
155 
156 #endif
157 
158 #ifdef LIBMESH_ENABLE_AMR
159 
166  static void compute_proj_constraints (DofConstraints & constraints,
167  DofMap & dof_map,
168  const unsigned int variable_number,
169  const Elem * elem);
170 
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 unsigned int var,
183  const bool use_old_dof_indices = false);
184 
191  static void coarsened_dof_values(const NumericVector<Number> & global_vector,
192  const DofMap & dof_map,
193  const Elem * coarse_elem,
194  DenseVector<Number> & coarse_dofs,
195  const bool use_old_dof_indices = false);
196 
197 #endif // #ifdef LIBMESH_ENABLE_AMR
198 
199 #ifdef LIBMESH_ENABLE_PERIODIC
200 
206  static void compute_periodic_constraints (DofConstraints & constraints,
207  DofMap & dof_map,
208  const PeriodicBoundaries & boundaries,
209  const MeshBase & mesh,
210  const PointLocatorBase * point_locator,
211  const unsigned int variable_number,
212  const Elem * elem);
213 
214 #endif // LIBMESH_ENABLE_PERIODIC
215 
220  const std::vector<std::vector<OutputShape> > & get_phi() const
222  calculate_phi = true; return phi; }
223 
228  const std::vector<std::vector<OutputGradient> > & get_dphi() const
230  calculate_dphi = calculate_dphiref = true; return dphi; }
231 
236  const std::vector<std::vector<OutputShape> > & get_curl_phi() const
238  calculate_curl_phi = calculate_dphiref = true; return curl_phi; }
239 
244  const std::vector<std::vector<OutputDivergence> > & get_div_phi() const
246  calculate_div_phi = calculate_dphiref = true; return div_phi; }
247 
252  const std::vector<std::vector<OutputShape> > & get_dphidx() const
254  calculate_dphi = calculate_dphiref = true; return dphidx; }
255 
260  const std::vector<std::vector<OutputShape> > & get_dphidy() const
262  calculate_dphi = calculate_dphiref = true; return dphidy; }
263 
268  const std::vector<std::vector<OutputShape> > & get_dphidz() const
270  calculate_dphi = calculate_dphiref = true; return dphidz; }
271 
276  const std::vector<std::vector<OutputShape> > & get_dphidxi() const
278  calculate_dphiref = true; return dphidxi; }
279 
284  const std::vector<std::vector<OutputShape> > & get_dphideta() const
286  calculate_dphiref = true; return dphideta; }
287 
292  const std::vector<std::vector<OutputShape> > & get_dphidzeta() const
294  calculate_dphiref = true; return dphidzeta; }
295 
296 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
297 
302  const std::vector<std::vector<OutputTensor> > & get_d2phi() const
304  calculate_d2phi = calculate_dphiref = true; return d2phi; }
305 
310  const std::vector<std::vector<OutputShape> > & get_d2phidx2() const
312  calculate_d2phi = calculate_dphiref = true; return d2phidx2; }
313 
318  const std::vector<std::vector<OutputShape> > & get_d2phidxdy() const
320  calculate_d2phi = calculate_dphiref = true; return d2phidxdy; }
321 
326  const std::vector<std::vector<OutputShape> > & get_d2phidxdz() const
328  calculate_d2phi = calculate_dphiref = true; return d2phidxdz; }
329 
334  const std::vector<std::vector<OutputShape> > & get_d2phidy2() const
336  calculate_d2phi = calculate_dphiref = true; return d2phidy2; }
337 
342  const std::vector<std::vector<OutputShape> > & get_d2phidydz() const
344  calculate_d2phi = calculate_dphiref = true; return d2phidydz; }
345 
350  const std::vector<std::vector<OutputShape> > & get_d2phidz2() const
352  calculate_d2phi = calculate_dphiref = true; return d2phidz2; }
353 
358  const std::vector<std::vector<OutputShape> > & get_d2phidxi2() const
360  calculate_d2phi = calculate_dphiref = true; return d2phidxi2; }
361 
366  const std::vector<std::vector<OutputShape> > & get_d2phidxideta() const
369 
374  const std::vector<std::vector<OutputShape> > & get_d2phidxidzeta() const
377 
382  const std::vector<std::vector<OutputShape> > & get_d2phideta2() const
384  calculate_d2phi = calculate_dphiref = true; return d2phideta2; }
385 
390  const std::vector<std::vector<OutputShape> > & get_d2phidetadzeta() const
393 
398  const std::vector<std::vector<OutputShape> > & get_d2phidzeta2() const
400  calculate_d2phi = calculate_dphiref = true; return d2phidzeta2; }
401 
402 #endif //LIBMESH_ENABLE_SECOND_DERIVATIVES
403 
404 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
405 
416  const std::vector<OutputGradient> & get_dphase() const
417  { return dphase; }
418 
419 
432  const std::vector<Real> & get_Sobolev_weight() const
433  { return weight; }
434 
440  const std::vector<RealGradient> & get_Sobolev_dweight() const
441  { return dweight; }
442 
443 #endif
444 
445 
449  void print_phi(std::ostream & os) const;
450 
455  void print_dphi(std::ostream & os) const;
456 
457 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
458 
463  void print_d2phi(std::ostream & os) const;
464 
465 #endif
466 
467 
468 protected:
469 
470 
471 
472 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
473 
479  virtual void init_base_shape_functions(const std::vector<Point> & qp,
480  const Elem * e) = 0;
481 
482 #endif
483 
488  void determine_calculations();
489 
500  virtual void compute_shape_functions(const Elem * elem, const std::vector<Point> & qp);
501 
507 
511  std::vector<std::vector<OutputShape> > phi;
512 
516  std::vector<std::vector<OutputGradient> > dphi;
517 
521  std::vector<std::vector<OutputShape> > curl_phi;
522 
526  std::vector<std::vector<OutputDivergence> > div_phi;
527 
531  std::vector<std::vector<OutputShape> > dphidxi;
532 
536  std::vector<std::vector<OutputShape> > dphideta;
537 
541  std::vector<std::vector<OutputShape> > dphidzeta;
542 
546  std::vector<std::vector<OutputShape> > dphidx;
547 
551  std::vector<std::vector<OutputShape> > dphidy;
552 
556  std::vector<std::vector<OutputShape> > dphidz;
557 
558 
559 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
560 
564  std::vector<std::vector<OutputTensor> > d2phi;
565 
569  std::vector<std::vector<OutputShape> > d2phidxi2;
570 
574  std::vector<std::vector<OutputShape> > d2phidxideta;
575 
579  std::vector<std::vector<OutputShape> > d2phidxidzeta;
580 
584  std::vector<std::vector<OutputShape> > d2phideta2;
585 
589  std::vector<std::vector<OutputShape> > d2phidetadzeta;
590 
594  std::vector<std::vector<OutputShape> > d2phidzeta2;
595 
599  std::vector<std::vector<OutputShape> > d2phidx2;
600 
604  std::vector<std::vector<OutputShape> > d2phidxdy;
605 
609  std::vector<std::vector<OutputShape> > d2phidxdz;
610 
614  std::vector<std::vector<OutputShape> > d2phidy2;
615 
619  std::vector<std::vector<OutputShape> > d2phidydz;
620 
624  std::vector<std::vector<OutputShape> > d2phidz2;
625 
626 #endif
627 
628 
629 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
630 
631  //--------------------------------------------------------------
632  /* protected members for infinite elements, which are accessed
633  * from the outside through some inline functions
634  */
635 
636 
642  std::vector<OutputGradient> dphase;
643 
649  std::vector<RealGradient> dweight;
650 
656  std::vector<Real> weight;
657 
658 #endif
659 
660 private:
661 
662 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
663 
669  template <unsigned int friend_Dim, FEFamily friend_T_radial, InfMapType friend_T_map>
670  friend class InfFE;
671 
672 #endif
673 
674 
675 };
676 
677 
678 // Typedefs for convenience and backwards compatibility
681 
682 
683 
684 
685 // ------------------------------------------------------------
686 // FEGenericBase class inline members
687 template <typename OutputType>
688 inline
690  const FEType & fet) :
691  FEAbstract(d,fet),
692  _fe_trans( FETransformationBase<OutputType>::build(fet) ),
693  phi(),
694  dphi(),
695  curl_phi(),
696  div_phi(),
697  dphidxi(),
698  dphideta(),
699  dphidzeta(),
700  dphidx(),
701  dphidy(),
702  dphidz()
703 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
704  ,d2phi(),
705  d2phidxi2(),
706  d2phidxideta(),
707  d2phidxidzeta(),
708  d2phideta2(),
709  d2phidetadzeta(),
710  d2phidzeta2(),
711  d2phidx2(),
712  d2phidxdy(),
713  d2phidxdz(),
714  d2phidy2(),
715  d2phidydz(),
716  d2phidz2()
717 #endif
718 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
719  ,dphase(),
720  dweight(),
721  weight()
722 #endif
723 {
724 }
725 
726 
727 
728 template <typename OutputType>
729 inline
731 {
732 }
733 
734 } // namespace libMesh
735 
736 #endif // LIBMESH_FE_BASE_H
FEGenericBase(const unsigned int dim, const FEType &fet)
Definition: fe_base.h:689
void print_phi(std::ostream &os) const
Definition: fe_base.C:720
void print_d2phi(std::ostream &os) const
Definition: fe_base.C:789
std::vector< std::vector< OutputTensor > > d2phi
Definition: fe_base.h:564
Manages the family, order, etc. parameters for a given FE.
Definition: fe_type.h:178
std::vector< std::vector< OutputShape > > dphidxi
Definition: fe_base.h:531
const std::vector< std::vector< OutputShape > > & get_d2phidydz() const
Definition: fe_base.h:342
std::vector< std::vector< OutputShape > > d2phidxdz
Definition: fe_base.h:609
std::vector< std::vector< OutputShape > > dphidzeta
Definition: fe_base.h:541
std::vector< std::vector< OutputShape > > d2phidydz
Definition: fe_base.h:619
const std::vector< std::vector< OutputDivergence > > & get_div_phi() const
Definition: fe_base.h:244
Base class for all the infinite geometric element types.
Definition: fe.h:40
TensorTools::DecrementRank< OutputNumber >::type OutputNumberDivergence
Definition: fe_base.h:139
const std::vector< std::vector< OutputShape > > & get_dphidzeta() const
Definition: fe_base.h:292
virtual void init_base_shape_functions(const std::vector< Point > &qp, const Elem *e)=0
const std::vector< std::vector< OutputShape > > & get_dphidz() const
Definition: fe_base.h:268
const std::vector< std::vector< OutputShape > > & get_dphideta() const
Definition: fe_base.h:284
const std::vector< std::vector< OutputShape > > & get_d2phidxdz() const
Definition: fe_base.h:326
const std::vector< std::vector< OutputShape > > & get_d2phidxdy() const
Definition: fe_base.h:318
Maps between boundary ids and PeriodicBoundaryBase objects.
virtual void compute_shape_functions(const Elem *elem, const std::vector< Point > &qp)
Definition: fe_base.C:682
std::vector< std::vector< OutputShape > > d2phidxideta
Definition: fe_base.h:574
const std::vector< std::vector< OutputShape > > & get_d2phidxidzeta() const
Definition: fe_base.h:374
The base class for all geometric element types.
Definition: elem.h:86
TensorTools::IncrementRank< OutputNumberGradient >::type OutputNumberTensor
Definition: fe_base.h:138
MeshBase & mesh
std::vector< Real > weight
Definition: fe_base.h:656
TensorTools::IncrementRank< OutputNumber >::type OutputNumberGradient
Definition: fe_base.h:137
const std::vector< std::vector< OutputShape > > & get_d2phideta2() const
Definition: fe_base.h:382
TensorTools::IncrementRank< OutputShape >::type OutputGradient
Definition: fe_base.h:133
std::vector< std::vector< OutputShape > > d2phidx2
Definition: fe_base.h:599
std::vector< std::vector< OutputShape > > curl_phi
Definition: fe_base.h:521
std::vector< std::vector< OutputShape > > dphidy
Definition: fe_base.h:551
FEGenericBase< RealGradient > FEVectorBase
Definition: fe_base.h:680
TensorTools::DecrementRank< OutputShape >::type OutputDivergence
Definition: fe_base.h:135
std::vector< std::vector< OutputShape > > d2phidy2
Definition: fe_base.h:614
const std::vector< std::vector< OutputShape > > & get_d2phidzeta2() const
Definition: fe_base.h:398
Base class for Mesh.
Definition: mesh_base.h:67
libmesh_assert(j)
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
std::vector< std::vector< OutputShape > > d2phidetadzeta
Definition: fe_base.h:589
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)
Definition: fe_base.C:804
std::vector< std::vector< OutputShape > > d2phidxidzeta
Definition: fe_base.h:579
Manages the degrees of freedom (DOFs) in a simulation.
Definition: dof_map.h:167
std::vector< std::vector< OutputShape > > d2phidxdy
Definition: fe_base.h:604
std::vector< std::vector< OutputShape > > dphidx
Definition: fe_base.h:546
std::vector< std::vector< OutputShape > > phi
Definition: fe_base.h:511
std::vector< std::vector< OutputShape > > d2phideta2
Definition: fe_base.h:584
const std::vector< Real > & get_Sobolev_weight() const
Definition: fe_base.h:432
const unsigned int dim
Definition: fe_abstract.h:516
const std::vector< std::vector< OutputShape > > & get_d2phidxi2() const
Definition: fe_base.h:358
TensorTools::IncrementRank< OutputGradient >::type OutputTensor
Definition: fe_base.h:134
const std::vector< std::vector< OutputShape > > & get_phi() const
Definition: fe_base.h:220
const std::vector< std::vector< OutputShape > > & get_d2phidz2() const
Definition: fe_base.h:350
std::vector< OutputGradient > dphase
Definition: fe_base.h:642
const std::vector< std::vector< OutputShape > > & get_dphidxi() const
Definition: fe_base.h:276
const std::vector< std::vector< OutputShape > > & get_curl_phi() const
Definition: fe_base.h:236
static void compute_proj_constraints(DofConstraints &constraints, DofMap &dof_map, const unsigned int variable_number, const Elem *elem)
Definition: fe_base.C:1374
OutputType OutputShape
Definition: fe_base.h:132
TensorTools::MakeNumber< OutputShape >::type OutputNumber
Definition: fe_base.h:136
std::vector< std::vector< OutputDivergence > > div_phi
Definition: fe_base.h:526
static UniquePtr< FEGenericBase > build_InfFE(const unsigned int dim, const FEType &type)
FEGenericBase< Real > FEBase
void determine_calculations()
Definition: fe_base.C:741
std::vector< std::vector< OutputGradient > > dphi
Definition: fe_base.h:516
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)
Definition: fe_base.C:1658
void print_dphi(std::ostream &os) const
Definition: fe_base.C:731
const std::vector< RealGradient > & get_Sobolev_dweight() const
Definition: fe_base.h:440
const std::vector< std::vector< OutputGradient > > & get_dphi() const
Definition: fe_base.h:228
UniquePtr< FETransformationBase< OutputType > > _fe_trans
Definition: fe_base.h:506
const std::vector< std::vector< OutputShape > > & get_d2phidx2() const
Definition: fe_base.h:310
const std::vector< std::vector< OutputShape > > & get_dphidx() const
Definition: fe_base.h:252
const std::vector< std::vector< OutputTensor > > & get_d2phi() const
Definition: fe_base.h:302
const std::vector< std::vector< OutputShape > > & get_d2phidxideta() const
Definition: fe_base.h:366
std::vector< std::vector< OutputShape > > d2phidz2
Definition: fe_base.h:624
std::vector< std::vector< OutputShape > > d2phidxi2
Definition: fe_base.h:569
const std::vector< std::vector< OutputShape > > & get_dphidy() const
Definition: fe_base.h:260
const std::vector< std::vector< OutputShape > > & get_d2phidy2() const
Definition: fe_base.h:334
std::vector< std::vector< OutputShape > > d2phidzeta2
Definition: fe_base.h:594
std::vector< std::vector< OutputShape > > dphidz
Definition: fe_base.h:556
const std::vector< OutputGradient > & get_dphase() const
Definition: fe_base.h:416
std::vector< std::vector< OutputShape > > dphideta
Definition: fe_base.h:536
std::vector< RealGradient > dweight
Definition: fe_base.h:649
const std::vector< std::vector< OutputShape > > & get_d2phidetadzeta() const
Definition: fe_base.h:390
static UniquePtr< FEGenericBase > build(const unsigned int dim, const FEType &type)