fe_abstract.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_ABSTRACT_H
21 #define LIBMESH_FE_ABSTRACT_H
22 
23 // Local includes
25 #include "libmesh/point.h"
26 #include "libmesh/vector_value.h"
27 #include "libmesh/enum_elem_type.h"
28 #include "libmesh/fe_type.h"
29 #include "libmesh/auto_ptr.h"
30 #include "libmesh/fe_map.h"
31 
32 // C++ includes
33 #include <cstddef>
34 #include <vector>
35 
36 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
37 #include "libmesh/tensor_value.h"
38 #endif
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 
55 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
56 class NodeConstraints;
57 #endif
58 
59 #ifdef LIBMESH_ENABLE_PERIODIC
60 class PeriodicBoundaries;
61 class PointLocatorBase;
62 #endif
63 
64 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
65 template <unsigned int Dim, FEFamily T_radial, InfMapType T_map>
66 class InfFE;
67 #endif
68 
69 
70 
95 // ------------------------------------------------------------
96 // FEAbstract class definition
97 class FEAbstract : public ReferenceCountedObject<FEAbstract>
98 {
99 protected:
100 
106  FEAbstract (const unsigned int dim,
107  const FEType& fet);
108 
109 public:
110 
114  virtual ~FEAbstract();
115 
121  static AutoPtr<FEAbstract> build (const unsigned int dim,
122  const FEType& type);
123 
138  virtual void reinit (const Elem* elem,
139  const std::vector<Point>* const pts = NULL,
140  const std::vector<Real>* const weights = NULL) = 0;
141 
150  virtual void reinit (const Elem* elem,
151  const unsigned int side,
152  const Real tolerance = TOLERANCE,
153  const std::vector<Point>* const pts = NULL,
154  const std::vector<Real>* const weights = NULL) = 0;
155 
164  virtual void edge_reinit (const Elem* elem,
165  const unsigned int edge,
166  const Real tolerance = TOLERANCE,
167  const std::vector<Point>* pts = NULL,
168  const std::vector<Real>* weights = NULL) = 0;
169 
174  virtual void side_map (const Elem* elem,
175  const Elem* side,
176  const unsigned int s,
177  const std::vector<Point>& reference_side_points,
178  std::vector<Point>& reference_points) = 0;
179 
187  static bool on_reference_element(const Point& p,
188  const ElemType t,
189  const Real eps = TOLERANCE);
194  static void get_refspace_nodes(const ElemType t,
195  std::vector<Point>& nodes);
196 
197 
198 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
199 
203  static void compute_node_constraints (NodeConstraints &constraints,
204  const Elem* elem);
205 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
206 
207 #ifdef LIBMESH_ENABLE_PERIODIC
208 
209 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
210 
214  static void compute_periodic_node_constraints (NodeConstraints &constraints,
215  const PeriodicBoundaries &boundaries,
216  const MeshBase& mesh,
217  const PointLocatorBase* point_locator,
218  const Elem* elem);
219 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
220 
221 #endif // LIBMESH_ENABLE_PERIODIC
222 
227  const std::vector<Point>& get_xyz() const
228  { return this->_fe_map->get_xyz(); }
229 
234  const std::vector<Real>& get_JxW() const
235  { return this->_fe_map->get_JxW(); }
236 
241  const std::vector<RealGradient>& get_dxyzdxi() const
242  { return this->_fe_map->get_dxyzdxi(); }
243 
248  const std::vector<RealGradient>& get_dxyzdeta() const
249  { return this->_fe_map->get_dxyzdeta(); }
250 
255  const std::vector<RealGradient>& get_dxyzdzeta() const
256  { return _fe_map->get_dxyzdzeta(); }
257 
261  const std::vector<RealGradient>& get_d2xyzdxi2() const
262  { return this->_fe_map->get_d2xyzdxi2(); }
263 
267  const std::vector<RealGradient>& get_d2xyzdeta2() const
268  { return this->_fe_map->get_d2xyzdeta2(); }
269 
270 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
271 
275  const std::vector<RealGradient>& get_d2xyzdzeta2() const
276  { return this->_fe_map->get_d2xyzdzeta2(); }
277 
278 #endif
279 
283  const std::vector<RealGradient>& get_d2xyzdxideta() const
284  { return this->_fe_map->get_d2xyzdxideta(); }
285 
286 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
287 
291  const std::vector<RealGradient>& get_d2xyzdxidzeta() const
292  { return this->_fe_map->get_d2xyzdxidzeta(); }
293 
297  const std::vector<RealGradient>& get_d2xyzdetadzeta() const
298  { return this->_fe_map->get_d2xyzdetadzeta(); }
299 
300 #endif
301 
306  const std::vector<Real>& get_dxidx() const
307  { return this->_fe_map->get_dxidx(); }
308 
313  const std::vector<Real>& get_dxidy() const
314  { return this->_fe_map->get_dxidy(); }
315 
320  const std::vector<Real>& get_dxidz() const
321  { return this->_fe_map->get_dxidz(); }
322 
327  const std::vector<Real>& get_detadx() const
328  { return this->_fe_map->get_detadx(); }
329 
334  const std::vector<Real>& get_detady() const
335  { return this->_fe_map->get_detady(); }
336 
341  const std::vector<Real>& get_detadz() const
342  { return this->_fe_map->get_detadz(); }
343 
348  const std::vector<Real>& get_dzetadx() const
349  { return this->_fe_map->get_dzetadx(); }
350 
355  const std::vector<Real>& get_dzetady() const
356  { return this->_fe_map->get_dzetady(); }
357 
362  const std::vector<Real>& get_dzetadz() const
363  { return this->_fe_map->get_dzetadz(); }
364 
368  const std::vector<std::vector<Point> >& get_tangents() const
369  { return this->_fe_map->get_tangents(); }
370 
374  const std::vector<Point>& get_normals() const
375  { return this->_fe_map->get_normals(); }
376 
380  const std::vector<Real>& get_curvatures() const
381  { return this->_fe_map->get_curvatures();}
382 
387  virtual void attach_quadrature_rule (QBase* q) = 0;
388 
394  virtual unsigned int n_shape_functions () const = 0;
395 
400  virtual unsigned int n_quadrature_points () const = 0;
401 
407  ElemType get_type() const { return elem_type; }
408 
413  unsigned int get_p_level() const { return _p_level; }
414 
418  FEType get_fe_type() const { return fe_type; }
419 
423  Order get_order() const { return static_cast<Order>(fe_type.order + _p_level); }
424 
428  virtual FEContinuity get_continuity() const = 0;
429 
434  virtual bool is_hierarchic() const = 0;
435 
439  FEFamily get_family() const { return fe_type.family; }
440 
444  const FEMap& get_fe_map() const { return *_fe_map.get(); }
445 
449  void print_JxW(std::ostream& os) const;
450 
456  virtual void print_phi(std::ostream& os) const =0;
457 
463  virtual void print_dphi(std::ostream& os) const =0;
464 
465 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
466 
472  virtual void print_d2phi(std::ostream& os) const =0;
473 
474 #endif
475 
480  void print_xyz(std::ostream& os) const;
481 
485  void print_info(std::ostream& os) const;
486 
490  friend std::ostream& operator << (std::ostream& os, const FEAbstract& fe);
491 
492 
493 protected:
494 
507  virtual void compute_shape_functions(const Elem*, const std::vector<Point>& ) =0;
508 
510 
511 
515  const unsigned int dim;
516 
521  mutable bool calculations_started;
522 
526  mutable bool calculate_phi;
527 
531  mutable bool calculate_dphi;
532 
536  mutable bool calculate_d2phi;
537 
541  mutable bool calculate_curl_phi;
542 
546  mutable bool calculate_div_phi;
547 
551  mutable bool calculate_dphiref;
552 
553 
559 
565 
570  unsigned int _p_level;
571 
576 
582 
589  virtual bool shapes_need_reinit() const = 0;
590 
591 };
592 
593 
594 
595 
596 // ------------------------------------------------------------
597 // FEAbstract class inline members
598 inline
599 FEAbstract::FEAbstract(const unsigned int d,
600  const FEType& fet) :
601  _fe_map( FEMap::build(fet) ),
602  dim(d),
603  calculations_started(false),
604  calculate_phi(false),
605  calculate_dphi(false),
606  calculate_d2phi(false),
607  calculate_curl_phi(false),
608  calculate_div_phi(false),
609  calculate_dphiref(false),
610  fe_type(fet),
611  elem_type(INVALID_ELEM),
612  _p_level(0),
613  qrule(NULL),
614  shapes_on_quadrature(false)
615 {
616 }
617 
618 
619 inline
621 {
622 }
623 
624 } // namespace libMesh
625 
626 #endif // LIBMESH_FE_ABSTRACT_H

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

Hosted By:
SourceForge.net Logo