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 
85 template <typename OutputType>
86 class FEGenericBase : public FEAbstract
87 {
88 protected:
89 
95  FEGenericBase (const unsigned int dim,
96  const FEType & fet);
97 
98 public:
99 
103  virtual ~FEGenericBase();
104 
113  static UniquePtr<FEGenericBase> build (const unsigned int dim,
114  const FEType & type);
115 
120  typedef OutputType OutputShape;
128 
129 
130 
131 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
132 
141  static UniquePtr<FEGenericBase> build_InfFE (const unsigned int dim,
142  const FEType & type);
143 
144 #endif
145 
146 #ifdef LIBMESH_ENABLE_AMR
147 
154  static void compute_proj_constraints (DofConstraints & constraints,
155  DofMap & dof_map,
156  const unsigned int variable_number,
157  const Elem * elem);
158 
166  static void coarsened_dof_values(const NumericVector<Number> & global_vector,
167  const DofMap & dof_map,
168  const Elem * coarse_elem,
169  DenseVector<Number> & coarse_dofs,
170  const unsigned int var,
171  const bool use_old_dof_indices = false);
172 
179  static void coarsened_dof_values(const NumericVector<Number> & global_vector,
180  const DofMap & dof_map,
181  const Elem * coarse_elem,
182  DenseVector<Number> & coarse_dofs,
183  const bool use_old_dof_indices = false);
184 
185 #endif // #ifdef LIBMESH_ENABLE_AMR
186 
187 #ifdef LIBMESH_ENABLE_PERIODIC
188 
194  static void compute_periodic_constraints (DofConstraints & constraints,
195  DofMap & dof_map,
196  const PeriodicBoundaries & boundaries,
197  const MeshBase & mesh,
198  const PointLocatorBase * point_locator,
199  const unsigned int variable_number,
200  const Elem * elem);
201 
202 #endif // LIBMESH_ENABLE_PERIODIC
203 
208  const std::vector<std::vector<OutputShape> > & get_phi() const
210  calculate_phi = true; return phi; }
211 
216  const std::vector<std::vector<OutputGradient> > & get_dphi() const
218  calculate_dphi = calculate_dphiref = true; return dphi; }
219 
224  const std::vector<std::vector<OutputShape> > & get_curl_phi() const
226  calculate_curl_phi = calculate_dphiref = true; return curl_phi; }
227 
232  const std::vector<std::vector<OutputDivergence> > & get_div_phi() const
234  calculate_div_phi = calculate_dphiref = true; return div_phi; }
235 
240  const std::vector<std::vector<OutputShape> > & get_dphidx() const
242  calculate_dphi = calculate_dphiref = true; return dphidx; }
243 
248  const std::vector<std::vector<OutputShape> > & get_dphidy() const
250  calculate_dphi = calculate_dphiref = true; return dphidy; }
251 
256  const std::vector<std::vector<OutputShape> > & get_dphidz() const
258  calculate_dphi = calculate_dphiref = true; return dphidz; }
259 
264  const std::vector<std::vector<OutputShape> > & get_dphidxi() const
266  calculate_dphiref = true; return dphidxi; }
267 
272  const std::vector<std::vector<OutputShape> > & get_dphideta() const
274  calculate_dphiref = true; return dphideta; }
275 
280  const std::vector<std::vector<OutputShape> > & get_dphidzeta() const
282  calculate_dphiref = true; return dphidzeta; }
283 
284 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
285 
290  const std::vector<std::vector<OutputTensor> > & get_d2phi() const
292  calculate_d2phi = calculate_dphiref = true; return d2phi; }
293 
298  const std::vector<std::vector<OutputShape> > & get_d2phidx2() const
300  calculate_d2phi = calculate_dphiref = true; return d2phidx2; }
301 
306  const std::vector<std::vector<OutputShape> > & get_d2phidxdy() const
308  calculate_d2phi = calculate_dphiref = true; return d2phidxdy; }
309 
314  const std::vector<std::vector<OutputShape> > & get_d2phidxdz() const
316  calculate_d2phi = calculate_dphiref = true; return d2phidxdz; }
317 
322  const std::vector<std::vector<OutputShape> > & get_d2phidy2() const
324  calculate_d2phi = calculate_dphiref = true; return d2phidy2; }
325 
330  const std::vector<std::vector<OutputShape> > & get_d2phidydz() const
332  calculate_d2phi = calculate_dphiref = true; return d2phidydz; }
333 
338  const std::vector<std::vector<OutputShape> > & get_d2phidz2() const
340  calculate_d2phi = calculate_dphiref = true; return d2phidz2; }
341 
346  const std::vector<std::vector<OutputShape> > & get_d2phidxi2() const
348  calculate_d2phi = calculate_dphiref = true; return d2phidxi2; }
349 
354  const std::vector<std::vector<OutputShape> > & get_d2phidxideta() const
357 
362  const std::vector<std::vector<OutputShape> > & get_d2phidxidzeta() const
365 
370  const std::vector<std::vector<OutputShape> > & get_d2phideta2() const
372  calculate_d2phi = calculate_dphiref = true; return d2phideta2; }
373 
378  const std::vector<std::vector<OutputShape> > & get_d2phidetadzeta() const
381 
386  const std::vector<std::vector<OutputShape> > & get_d2phidzeta2() const
388  calculate_d2phi = calculate_dphiref = true; return d2phidzeta2; }
389 
390 #endif //LIBMESH_ENABLE_SECOND_DERIVATIVES
391 
392 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
393 
404  const std::vector<OutputGradient> & get_dphase() const
405  { return dphase; }
406 
407 
420  const std::vector<Real> & get_Sobolev_weight() const
421  { return weight; }
422 
428  const std::vector<RealGradient> & get_Sobolev_dweight() const
429  { return dweight; }
430 
431 #endif
432 
433 
437  void print_phi(std::ostream & os) const;
438 
443  void print_dphi(std::ostream & os) const;
444 
445 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
446 
451  void print_d2phi(std::ostream & os) const;
452 
453 #endif
454 
455 
456 protected:
457 
458 
459 
460 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
461 
467  virtual void init_base_shape_functions(const std::vector<Point> & qp,
468  const Elem * e) = 0;
469 
470 #endif
471 
476  void determine_calculations();
477 
488  virtual void compute_shape_functions(const Elem * elem, const std::vector<Point> & qp);
489 
495 
499  std::vector<std::vector<OutputShape> > phi;
500 
504  std::vector<std::vector<OutputGradient> > dphi;
505 
509  std::vector<std::vector<OutputShape> > curl_phi;
510 
514  std::vector<std::vector<OutputDivergence> > div_phi;
515 
519  std::vector<std::vector<OutputShape> > dphidxi;
520 
524  std::vector<std::vector<OutputShape> > dphideta;
525 
529  std::vector<std::vector<OutputShape> > dphidzeta;
530 
534  std::vector<std::vector<OutputShape> > dphidx;
535 
539  std::vector<std::vector<OutputShape> > dphidy;
540 
544  std::vector<std::vector<OutputShape> > dphidz;
545 
546 
547 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
548 
552  std::vector<std::vector<OutputTensor> > d2phi;
553 
557  std::vector<std::vector<OutputShape> > d2phidxi2;
558 
562  std::vector<std::vector<OutputShape> > d2phidxideta;
563 
567  std::vector<std::vector<OutputShape> > d2phidxidzeta;
568 
572  std::vector<std::vector<OutputShape> > d2phideta2;
573 
577  std::vector<std::vector<OutputShape> > d2phidetadzeta;
578 
582  std::vector<std::vector<OutputShape> > d2phidzeta2;
583 
587  std::vector<std::vector<OutputShape> > d2phidx2;
588 
592  std::vector<std::vector<OutputShape> > d2phidxdy;
593 
597  std::vector<std::vector<OutputShape> > d2phidxdz;
598 
602  std::vector<std::vector<OutputShape> > d2phidy2;
603 
607  std::vector<std::vector<OutputShape> > d2phidydz;
608 
612  std::vector<std::vector<OutputShape> > d2phidz2;
613 
614 #endif
615 
616 
617 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
618 
619  //--------------------------------------------------------------
620  /* protected members for infinite elements, which are accessed
621  * from the outside through some inline functions
622  */
623 
624 
630  std::vector<OutputGradient> dphase;
631 
637  std::vector<RealGradient> dweight;
638 
644  std::vector<Real> weight;
645 
646 #endif
647 
648 private:
649 
650 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
651 
657  template <unsigned int friend_Dim, FEFamily friend_T_radial, InfMapType friend_T_map>
658  friend class InfFE;
659 
660 #endif
661 
662 
663 };
664 
665 
666 // Typedefs for convenience and backwards compatibility
669 
670 
671 
672 
673 // ------------------------------------------------------------
674 // FEGenericBase class inline members
675 template <typename OutputType>
676 inline
678  const FEType & fet) :
679  FEAbstract(d,fet),
680  _fe_trans( FETransformationBase<OutputType>::build(fet) ),
681  phi(),
682  dphi(),
683  curl_phi(),
684  div_phi(),
685  dphidxi(),
686  dphideta(),
687  dphidzeta(),
688  dphidx(),
689  dphidy(),
690  dphidz()
691 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
692  ,d2phi(),
693  d2phidxi2(),
694  d2phidxideta(),
695  d2phidxidzeta(),
696  d2phideta2(),
697  d2phidetadzeta(),
698  d2phidzeta2(),
699  d2phidx2(),
700  d2phidxdy(),
701  d2phidxdz(),
702  d2phidy2(),
703  d2phidydz(),
704  d2phidz2()
705 #endif
706 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
707  ,dphase(),
708  dweight(),
709  weight()
710 #endif
711 {
712 }
713 
714 
715 
716 template <typename OutputType>
717 inline
719 {
720 }
721 
722 } // namespace libMesh
723 
724 #endif // LIBMESH_FE_BASE_H
FEGenericBase(const unsigned int dim, const FEType &fet)
Definition: fe_base.h:677
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:552
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:519
const std::vector< std::vector< OutputShape > > & get_d2phidydz() const
Definition: fe_base.h:330
std::vector< std::vector< OutputShape > > d2phidxdz
Definition: fe_base.h:597
std::vector< std::vector< OutputShape > > dphidzeta
Definition: fe_base.h:529
std::vector< std::vector< OutputShape > > d2phidydz
Definition: fe_base.h:607
const std::vector< std::vector< OutputDivergence > > & get_div_phi() const
Definition: fe_base.h:232
Base class for all the infinite geometric element types.
Definition: fe.h:40
TensorTools::DecrementRank< OutputNumber >::type OutputNumberDivergence
Definition: fe_base.h:127
const std::vector< std::vector< OutputShape > > & get_dphidzeta() const
Definition: fe_base.h:280
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:256
const std::vector< std::vector< OutputShape > > & get_dphideta() const
Definition: fe_base.h:272
const std::vector< std::vector< OutputShape > > & get_d2phidxdz() const
Definition: fe_base.h:314
const std::vector< std::vector< OutputShape > > & get_d2phidxdy() const
Definition: fe_base.h:306
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:562
const std::vector< std::vector< OutputShape > > & get_d2phidxidzeta() const
Definition: fe_base.h:362
The base class for all geometric element types.
Definition: elem.h:86
TensorTools::IncrementRank< OutputNumberGradient >::type OutputNumberTensor
Definition: fe_base.h:126
MeshBase & mesh
std::vector< Real > weight
Definition: fe_base.h:644
TensorTools::IncrementRank< OutputNumber >::type OutputNumberGradient
Definition: fe_base.h:125
const std::vector< std::vector< OutputShape > > & get_d2phideta2() const
Definition: fe_base.h:370
TensorTools::IncrementRank< OutputShape >::type OutputGradient
Definition: fe_base.h:121
std::vector< std::vector< OutputShape > > d2phidx2
Definition: fe_base.h:587
std::vector< std::vector< OutputShape > > curl_phi
Definition: fe_base.h:509
std::vector< std::vector< OutputShape > > dphidy
Definition: fe_base.h:539
FEGenericBase< RealGradient > FEVectorBase
Definition: fe_base.h:668
TensorTools::DecrementRank< OutputShape >::type OutputDivergence
Definition: fe_base.h:123
std::vector< std::vector< OutputShape > > d2phidy2
Definition: fe_base.h:602
const std::vector< std::vector< OutputShape > > & get_d2phidzeta2() const
Definition: fe_base.h:386
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:577
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:567
Manages the degrees of freedom (DOFs) in a simulation.
Definition: dof_map.h:167
std::vector< std::vector< OutputShape > > d2phidxdy
Definition: fe_base.h:592
std::vector< std::vector< OutputShape > > dphidx
Definition: fe_base.h:534
std::vector< std::vector< OutputShape > > phi
Definition: fe_base.h:499
std::vector< std::vector< OutputShape > > d2phideta2
Definition: fe_base.h:572
const std::vector< Real > & get_Sobolev_weight() const
Definition: fe_base.h:420
const unsigned int dim
Definition: fe_abstract.h:517
const std::vector< std::vector< OutputShape > > & get_d2phidxi2() const
Definition: fe_base.h:346
TensorTools::IncrementRank< OutputGradient >::type OutputTensor
Definition: fe_base.h:122
const std::vector< std::vector< OutputShape > > & get_phi() const
Definition: fe_base.h:208
const std::vector< std::vector< OutputShape > > & get_d2phidz2() const
Definition: fe_base.h:338
std::vector< OutputGradient > dphase
Definition: fe_base.h:630
const std::vector< std::vector< OutputShape > > & get_dphidxi() const
Definition: fe_base.h:264
const std::vector< std::vector< OutputShape > > & get_curl_phi() const
Definition: fe_base.h:224
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:120
TensorTools::MakeNumber< OutputShape >::type OutputNumber
Definition: fe_base.h:124
std::vector< std::vector< OutputDivergence > > div_phi
Definition: fe_base.h:514
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:504
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:428
const std::vector< std::vector< OutputGradient > > & get_dphi() const
Definition: fe_base.h:216
UniquePtr< FETransformationBase< OutputType > > _fe_trans
Definition: fe_base.h:494
const std::vector< std::vector< OutputShape > > & get_d2phidx2() const
Definition: fe_base.h:298
const std::vector< std::vector< OutputShape > > & get_dphidx() const
Definition: fe_base.h:240
const std::vector< std::vector< OutputTensor > > & get_d2phi() const
Definition: fe_base.h:290
const std::vector< std::vector< OutputShape > > & get_d2phidxideta() const
Definition: fe_base.h:354
std::vector< std::vector< OutputShape > > d2phidz2
Definition: fe_base.h:612
std::vector< std::vector< OutputShape > > d2phidxi2
Definition: fe_base.h:557
const std::vector< std::vector< OutputShape > > & get_dphidy() const
Definition: fe_base.h:248
const std::vector< std::vector< OutputShape > > & get_d2phidy2() const
Definition: fe_base.h:322
std::vector< std::vector< OutputShape > > d2phidzeta2
Definition: fe_base.h:582
std::vector< std::vector< OutputShape > > dphidz
Definition: fe_base.h:544
const std::vector< OutputGradient > & get_dphase() const
Definition: fe_base.h:404
std::vector< std::vector< OutputShape > > dphideta
Definition: fe_base.h:524
std::vector< RealGradient > dweight
Definition: fe_base.h:637
const std::vector< std::vector< OutputShape > > & get_d2phidetadzeta() const
Definition: fe_base.h:378
static UniquePtr< FEGenericBase > build(const unsigned int dim, const FEType &type)