fe_base.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 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
25 #include "libmesh/fe_abstract.h"
26 #include "libmesh/point.h"
27 #include "libmesh/vector_value.h"
28 #include "libmesh/enum_elem_type.h"
29 #include "libmesh/fe_type.h"
30 #include "libmesh/auto_ptr.h"
32 #include "libmesh/tensor_tools.h"
33 #include "libmesh/type_n_tensor.h"
34 
35 // C++ includes
36 #include <cstddef>
37 #include <vector>
38 
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 
56 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
57 class NodeConstraints;
58 #endif
59 
60 #ifdef LIBMESH_ENABLE_PERIODIC
61 class PeriodicBoundaries;
62 class PointLocatorBase;
63 #endif
64 
65 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
66 template <unsigned int Dim, FEFamily T_radial, InfMapType T_map>
67 class InfFE;
68 #endif
69 
70  // TypesEqual takes two types as parameters.
71  // If they are the exact same type, then TypesEqual::value is the boolean true,
72  // otherwise TypesEqual::value is the boolean false.
73  // FIXME: Need to put this in a common place
74  template <typename T1, typename T2>
75  struct TypesEqual {
76  static const bool value = false;
77  };
78 
79  template <typename T>
80  struct TypesEqual<T,T> {
81  static const bool value = true;
82  };
83 
111 // ------------------------------------------------------------
112 // FEBase class definition
113 template <typename OutputType>
114 class FEGenericBase : public FEAbstract
115 {
116 protected:
117 
123  FEGenericBase (const unsigned int dim,
124  const FEType& fet);
125 
126 public:
127 
131  virtual ~FEGenericBase();
132 
141  static AutoPtr<FEGenericBase> build (const unsigned int dim,
142  const FEType& type);
143 
148  typedef OutputType OutputShape;
156 
157 
158 
159 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
160 
169  static AutoPtr<FEGenericBase> build_InfFE (const unsigned int dim,
170  const FEType& type);
171 
172 #endif
173 
174 #ifdef LIBMESH_ENABLE_AMR
175 
182  static void compute_proj_constraints (DofConstraints &constraints,
183  DofMap &dof_map,
184  const unsigned int variable_number,
185  const Elem* elem);
186 
192  static void coarsened_dof_values(const NumericVector<Number> &global_vector,
193  const DofMap &dof_map,
194  const Elem *coarse_elem,
195  DenseVector<Number> &coarse_dofs,
196  const unsigned int var,
197  const bool use_old_dof_indices = false);
198 
199 #endif // #ifdef LIBMESH_ENABLE_AMR
200 
201 #ifdef LIBMESH_ENABLE_PERIODIC
202 
208  static void compute_periodic_constraints (DofConstraints &constraints,
209  DofMap &dof_map,
210  const PeriodicBoundaries &boundaries,
211  const MeshBase& mesh,
212  const PointLocatorBase* point_locator,
213  const unsigned int variable_number,
214  const Elem* elem);
215 
216 #endif // LIBMESH_ENABLE_PERIODIC
217 
222  const std::vector<std::vector<OutputShape> >& get_phi() const
224  calculate_phi = true; return phi; }
225 
230  const std::vector<std::vector<OutputGradient> >& get_dphi() const
232  calculate_dphi = calculate_dphiref = true; return dphi; }
233 
238  const std::vector<std::vector<OutputShape> >& get_curl_phi() const
240  calculate_curl_phi = calculate_dphiref = true; return curl_phi; }
241 
246  const std::vector<std::vector<OutputDivergence> >& get_div_phi() const
248  calculate_div_phi = calculate_dphiref = true; return div_phi; }
249 
254  const std::vector<std::vector<OutputShape> >& get_dphidx() const
256  calculate_dphi = calculate_dphiref = true; return dphidx; }
257 
262  const std::vector<std::vector<OutputShape> >& get_dphidy() const
264  calculate_dphi = calculate_dphiref = true; return dphidy; }
265 
270  const std::vector<std::vector<OutputShape> >& get_dphidz() const
272  calculate_dphi = calculate_dphiref = true; return dphidz; }
273 
278  const std::vector<std::vector<OutputShape> >& get_dphidxi() const
280  calculate_dphiref = true; return dphidxi; }
281 
286  const std::vector<std::vector<OutputShape> >& get_dphideta() const
288  calculate_dphiref = true; return dphideta; }
289 
294  const std::vector<std::vector<OutputShape> >& get_dphidzeta() const
296  calculate_dphiref = true; return dphidzeta; }
297 
298 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
299 
304  const std::vector<std::vector<OutputTensor> >& get_d2phi() const
306  calculate_d2phi = calculate_dphiref = true; return d2phi; }
307 
312  const std::vector<std::vector<OutputShape> >& get_d2phidx2() const
314  calculate_d2phi = calculate_dphiref = true; return d2phidx2; }
315 
320  const std::vector<std::vector<OutputShape> >& get_d2phidxdy() const
322  calculate_d2phi = calculate_dphiref = true; return d2phidxdy; }
323 
328  const std::vector<std::vector<OutputShape> >& get_d2phidxdz() const
330  calculate_d2phi = calculate_dphiref = true; return d2phidxdz; }
331 
336  const std::vector<std::vector<OutputShape> >& get_d2phidy2() const
338  calculate_d2phi = calculate_dphiref = true; return d2phidy2; }
339 
344  const std::vector<std::vector<OutputShape> >& get_d2phidydz() const
346  calculate_d2phi = calculate_dphiref = true; return d2phidydz; }
347 
352  const std::vector<std::vector<OutputShape> >& get_d2phidz2() const
354  calculate_d2phi = calculate_dphiref = true; return d2phidz2; }
355 
360  const std::vector<std::vector<OutputShape> >& get_d2phidxi2() const
362  calculate_d2phi = calculate_dphiref = true; return d2phidxi2; }
363 
368  const std::vector<std::vector<OutputShape> >& get_d2phidxideta() const
371 
376  const std::vector<std::vector<OutputShape> >& get_d2phidxidzeta() const
379 
384  const std::vector<std::vector<OutputShape> >& get_d2phideta2() const
386  calculate_d2phi = calculate_dphiref = true; return d2phideta2; }
387 
392  const std::vector<std::vector<OutputShape> >& get_d2phidetadzeta() const
395 
400  const std::vector<std::vector<OutputShape> >& get_d2phidzeta2() const
402  calculate_d2phi = calculate_dphiref = true; return d2phidzeta2; }
403 
404 #endif //LIBMESH_ENABLE_SECOND_DERIVATIVES
405 
406 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
407 
418  const std::vector<OutputGradient>& get_dphase() const
419  { return dphase; }
420 
421 
434  const std::vector<Real>& get_Sobolev_weight() const
435  { return weight; }
436 
442  const std::vector<RealGradient>& get_Sobolev_dweight() const
443  { return dweight; }
444 
445 #endif
446 
447 
451  void print_phi(std::ostream& os) const;
452 
457  void print_dphi(std::ostream& os) const;
458 
459 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
460 
465  void print_d2phi(std::ostream& os) const;
466 
467 #endif
468 
469 
470 protected:
471 
472 
473 
474 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
475 
481  virtual void init_base_shape_functions(const std::vector<Point>& qp,
482  const Elem* e) = 0;
483 
484 #endif
485 
496  virtual void compute_shape_functions(const Elem* elem, const std::vector<Point>& qp);
497 
503 
507  std::vector<std::vector<OutputShape> > phi;
508 
512  std::vector<std::vector<OutputGradient> > dphi;
513 
517  std::vector<std::vector<OutputShape> > curl_phi;
518 
522  std::vector<std::vector<OutputDivergence> > div_phi;
523 
527  std::vector<std::vector<OutputShape> > dphidxi;
528 
532  std::vector<std::vector<OutputShape> > dphideta;
533 
537  std::vector<std::vector<OutputShape> > dphidzeta;
538 
542  std::vector<std::vector<OutputShape> > dphidx;
543 
547  std::vector<std::vector<OutputShape> > dphidy;
548 
552  std::vector<std::vector<OutputShape> > dphidz;
553 
554 
555 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
556 
560  std::vector<std::vector<OutputTensor> > d2phi;
561 
565  std::vector<std::vector<OutputShape> > d2phidxi2;
566 
570  std::vector<std::vector<OutputShape> > d2phidxideta;
571 
575  std::vector<std::vector<OutputShape> > d2phidxidzeta;
576 
580  std::vector<std::vector<OutputShape> > d2phideta2;
581 
585  std::vector<std::vector<OutputShape> > d2phidetadzeta;
586 
590  std::vector<std::vector<OutputShape> > d2phidzeta2;
591 
595  std::vector<std::vector<OutputShape> > d2phidx2;
596 
600  std::vector<std::vector<OutputShape> > d2phidxdy;
601 
605  std::vector<std::vector<OutputShape> > d2phidxdz;
606 
610  std::vector<std::vector<OutputShape> > d2phidy2;
611 
615  std::vector<std::vector<OutputShape> > d2phidydz;
616 
620  std::vector<std::vector<OutputShape> > d2phidz2;
621 
622 #endif
623 
624 
625 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
626 
627  //--------------------------------------------------------------
628  /* protected members for infinite elements, which are accessed
629  * from the outside through some inline functions
630  */
631 
632 
638  std::vector<OutputGradient> dphase;
639 
645  std::vector<RealGradient> dweight;
646 
652  std::vector<Real> weight;
653 
654 #endif
655 
656 private:
657 
658 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
659 
665  template <unsigned int friend_Dim, FEFamily friend_T_radial, InfMapType friend_T_map>
666  friend class InfFE;
667 
668 #endif
669 
670 
671 };
672 
673 
674 // Typedefs for convenience and backwards compatibility
677 
678 
679 
680 
681 // ------------------------------------------------------------
682 // FEGenericBase class inline members
683 template <typename OutputType>
684 inline
686  const FEType& fet) :
687  FEAbstract(d,fet),
688  _fe_trans( FETransformationBase<OutputType>::build(fet) ),
689  phi(),
690  dphi(),
691  curl_phi(),
692  div_phi(),
693  dphidxi(),
694  dphideta(),
695  dphidzeta(),
696  dphidx(),
697  dphidy(),
698  dphidz()
699 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
700  ,d2phi(),
701  d2phidxi2(),
702  d2phidxideta(),
703  d2phidxidzeta(),
704  d2phideta2(),
705  d2phidetadzeta(),
706  d2phidzeta2(),
707  d2phidx2(),
708  d2phidxdy(),
709  d2phidxdz(),
710  d2phidy2(),
711  d2phidydz(),
712  d2phidz2()
713 #endif
714 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
715  ,dphase(),
716  dweight(),
717  weight()
718 #endif
719 {
720 }
721 
722 
723 
724 template <typename OutputType>
725 inline
727 {
728 }
729 
730 } // namespace libMesh
731 
732 #endif // LIBMESH_FE_BASE_H

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:04 UTC

Hosted By:
SourceForge.net Logo