libMesh
inf_fe.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_INF_FE_H
21 #define LIBMESH_INF_FE_H
22 
23 #include "libmesh/libmesh_config.h"
24 
25 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
26 
27 // Local includes
28 #include "libmesh/fe_base.h"
29 
30 // C++ includes
31 #include <cstddef>
32 
33 namespace libMesh
34 {
35 
36 
37 // forward declarations
38 class Elem;
39 class FEComputeData;
40 
41 
42 
74 template <unsigned int Dim, FEFamily T_radial, InfMapType T_map>
75 class InfFE : public FEBase
76 {
77 
78  /*
79  * Protect the nested class
80  */
81 protected:
82 
93  class Radial
94  {
95  private:
96 
100  Radial() {}
101 
102  public:
103 
108  static Real decay (const Real v);
109 
114  static Real decay_deriv (const Real) { return -.5; }
115 
120  static Real D (const Real v) { return (1.-v)*(1.-v)/4.; }
121 
125  static Real D_deriv (const Real v) { return (v-1.)/2.; }
126 
131  static Order mapping_order() { return FIRST; }
132 
147  static unsigned int n_dofs (const Order o_radial)
148  { return static_cast<unsigned int>(o_radial)+1; }
149 
159  static unsigned int n_dofs_at_node (const Order o_radial,
160  const unsigned int n_onion);
161 
170  static unsigned int n_dofs_per_elem (const Order o_radial)
171  { return static_cast<unsigned int>(o_radial)+1; }
172 
173  };
174 
175 
176 
185  class Base
186  {
187  private:
188 
192  Base() {}
193 
194  public:
195 
201  static Elem * build_elem (const Elem * inf_elem);
202 
208  static ElemType get_elem_type (const ElemType type);
209 
215  static unsigned int n_base_mapping_sf (const ElemType base_elem_type,
216  const Order base_mapping_order);
217 
218  };
219 
220 
221 
222 
223 
224 
225 
226 public:
227 
228  // InfFE continued
229 
243  explicit
244  InfFE(const FEType & fet);
245  ~InfFE() {}
246 
247  // The static public members for access from FEInterface etc
248 
261  static Real shape(const FEType & fet,
262  const ElemType t,
263  const unsigned int i,
264  const Point & p);
265 
278  static Real shape(const FEType & fet,
279  const Elem * elem,
280  const unsigned int i,
281  const Point & p);
282 
293  static void compute_data(const FEType & fe_t,
294  const Elem * inf_elem,
295  FEComputeData & data);
296 
301  static unsigned int n_shape_functions (const FEType & fet,
302  const ElemType t)
303  { return n_dofs(fet, t); }
304 
311  static unsigned int n_dofs(const FEType & fet,
312  const ElemType inf_elem_type);
313 
318  static unsigned int n_dofs_at_node(const FEType & fet,
319  const ElemType inf_elem_type,
320  const unsigned int n);
321 
326  static unsigned int n_dofs_per_elem(const FEType & fet,
327  const ElemType inf_elem_type);
328 
332  virtual FEContinuity get_continuity() const libmesh_override
333  { return C_ZERO; } // FIXME - is this true??
334 
339  virtual bool is_hierarchic() const libmesh_override
340  { return false; } // FIXME - Inf FEs don't handle p elevation yet
341 
349  static void nodal_soln(const FEType & fet,
350  const Elem * elem,
351  const std::vector<Number> & elem_soln,
352  std::vector<Number> & nodal_soln);
353 
358  static Point map (const Elem * inf_elem,
359  const Point & reference_point);
360 
376  static Point inverse_map (const Elem * elem,
377  const Point & p,
378  const Real tolerance = TOLERANCE,
379  const bool secure = true);
380 
381 
388  static void inverse_map (const Elem * elem,
389  const std::vector<Point> & physical_points,
390  std::vector<Point> & reference_points,
391  const Real tolerance = TOLERANCE,
392  const bool secure = true);
393 
394 
395  // The workhorses of InfFE. These are often used during matrix assembly.
396 
403  virtual void reinit (const Elem * elem,
404  const std::vector<Point> * const pts = libmesh_nullptr,
405  const std::vector<Real> * const weights = libmesh_nullptr) libmesh_override;
406 
412  virtual void reinit (const Elem * elem,
413  const unsigned int side,
414  const Real tolerance = TOLERANCE,
415  const std::vector<Point> * const pts = libmesh_nullptr,
416  const std::vector<Real> * const weights = libmesh_nullptr) libmesh_override;
417 
423  virtual void edge_reinit (const Elem * elem,
424  const unsigned int edge,
425  const Real tolerance = TOLERANCE,
426  const std::vector<Point> * const pts = libmesh_nullptr,
427  const std::vector<Real> * const weights = libmesh_nullptr) libmesh_override;
428 
433  virtual void side_map (const Elem * /* elem */,
434  const Elem * /* side */,
435  const unsigned int /* s */,
436  const std::vector<Point> & /* reference_side_points */,
437  std::vector<Point> & /* reference_points */) libmesh_override
438  {
439  libmesh_not_implemented();
440  }
441 
453  virtual void attach_quadrature_rule (QBase * q) libmesh_override;
454 
459  virtual unsigned int n_shape_functions () const libmesh_override
460  { return _n_total_approx_sf; }
461 
467  virtual unsigned int n_quadrature_points () const libmesh_override
469 
470 
471 protected:
472 
473  // static members used by the workhorses
474 
491  static Real eval(Real v,
492  Order o_radial,
493  unsigned int i);
494 
500  static Real eval_deriv(Real v,
501  Order o_radial,
502  unsigned int i);
503 
504 
505 
506  // Non-static members used by the workhorses
507 
512  void update_base_elem (const Elem * inf_elem);
513 
517  virtual void init_base_shape_functions(const std::vector<Point> &,
518  const Elem *) libmesh_override
519  { libmesh_not_implemented(); }
520 
526  void init_radial_shape_functions(const Elem * inf_elem,
527  const std::vector<Point> * radial_pts = libmesh_nullptr);
528 
535  void init_shape_functions(const std::vector<Point> & radial_qp,
536  const std::vector<Point> & base_qp,
537  const Elem * inf_elem);
538 
543  void init_face_shape_functions (const std::vector<Point> & qp,
544  const Elem * side);
545 
554  void combine_base_radial(const Elem * inf_elem);
555 
566  virtual void compute_shape_functions(const Elem *, const std::vector<Point> &) libmesh_override;
567 
568 
569 
570  // Miscellaneous static members
571 
578  static void compute_node_indices (const ElemType inf_elem_type,
579  const unsigned int outer_node_index,
580  unsigned int & base_node,
581  unsigned int & radial_node);
582 
591  static void compute_node_indices_fast (const ElemType inf_elem_type,
592  const unsigned int outer_node_index,
593  unsigned int & base_node,
594  unsigned int & radial_node);
595 
602  static void compute_shape_indices (const FEType & fet,
603  const ElemType inf_elem_type,
604  const unsigned int i,
605  unsigned int & base_shape,
606  unsigned int & radial_shape);
607 
611  std::vector<Real> dist;
612 
619  std::vector<Real> dweightdv;
620 
628  std::vector<Real> som;
633  std::vector<Real> dsomdv;
634 
639  std::vector<std::vector<Real>> mode;
640 
645  std::vector<std::vector<Real>> dmodedv;
646 
650  std::vector<std::vector<Real>> radial_map;
651 
652 
656  std::vector<std::vector<Real>> dradialdv_map;
657 
663  std::vector<Real> dphasedxi;
664 
670  std::vector<Real> dphasedeta;
671 
677  std::vector<Real> dphasedzeta;
678 
679 
680 
681 
682  // numbering scheme maps
683 
692  std::vector<unsigned int> _radial_node_index;
693 
702  std::vector<unsigned int> _base_node_index;
703 
712  std::vector<unsigned int> _radial_shape_index;
713 
722  std::vector<unsigned int> _base_shape_index;
723 
724 
725 
726 
727  // some more protected members
728 
733  unsigned int _n_total_approx_sf;
734 
739  unsigned int _n_total_qp;
740 
745  std::vector<Real> _total_qrule_weights;
746 
752 
758 
764 
772 
782 
783 
784 private:
785 
789  virtual bool shapes_need_reinit() const libmesh_override;
790 
799 
800 
801 #ifdef DEBUG
802 
807  static bool _warned_for_shape;
808 
809 #endif
810 
816  template <unsigned int friend_Dim, FEFamily friend_T_radial, InfMapType friend_T_map>
817  friend class InfFE;
818 
819 };
820 
821 
822 
823 // InfFE::Radial class inline members
824 template <unsigned int Dim, FEFamily T_radial, InfMapType T_map>
825 inline
827 {
828  switch (Dim)
829  //TODO:[DD] What decay do i have in 2D and 1D?
830  {
831  case 3:
832  return (1.-v)/2.;
833 
834  case 2:
835  return 0.;
836 
837  case 1:
838  return 0.;
839 
840  default:
841  libmesh_error_msg("Invalid Dim = " << Dim);
842  }
843 }
844 
845 } // namespace libMesh
846 
847 
848 #endif // ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
849 
850 
851 #endif // LIBMESH_INF_FE_H
void init_face_shape_functions(const std::vector< Point > &qp, const Elem *side)
Not implemented yet.
static void compute_node_indices_fast(const ElemType inf_elem_type, const unsigned int outer_node_index, unsigned int &base_node, unsigned int &radial_node)
Does the same as compute_node_indices(), but stores the maps for the current element type...
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
virtual void init_base_shape_functions(const std::vector< Point > &, const Elem *) libmesh_override
Do not use this derived member in InfFE<Dim,T_radial,T_map>.
Definition: inf_fe.h:517
virtual bool is_hierarchic() const libmesh_override
Definition: inf_fe.h:339
This nested class contains most of the static methods related to the base part of an infinite element...
Definition: inf_fe.h:185
UniquePtr< QBase > radial_qrule
The quadrature rule for the base element associated with the current infinite element.
Definition: inf_fe.h:757
static Real decay_deriv(const Real)
Definition: inf_fe.h:114
std::vector< Real > dsomdv
the first local derivative of the radial decay in local coordinates.
Definition: inf_fe.h:633
static bool _warned_for_shape
Definition: inf_fe.h:807
std::vector< std::vector< Real > > dradialdv_map
the first local derivative of the radial mapping shapes
Definition: inf_fe.h:656
std::vector< std::vector< Real > > radial_map
the radial mapping shapes in local coordinates
Definition: inf_fe.h:650
static Real D(const Real v)
Definition: inf_fe.h:120
A specific instantiation of the FEBase class.
Definition: fe.h:40
virtual unsigned int n_quadrature_points() const libmesh_override
Definition: inf_fe.h:467
std::vector< Real > dphasedxi
the first local derivative (for 3D, the first in the base) of the phase term in local coordinates...
Definition: inf_fe.h:663
void combine_base_radial(const Elem *inf_elem)
Combines the shape functions, which were formed in init_shape_functions(Elem *), with geometric data...
Definition: inf_fe.C:715
virtual bool shapes_need_reinit() const libmesh_override
Definition: inf_fe.C:962
static Real eval_deriv(Real v, Order o_radial, unsigned int i)
std::vector< Real > dweightdv
the additional radial weight in local coordinates, over all quadrature points.
Definition: inf_fe.h:619
class FEComputeData hides arbitrary data to be passed to and from children of FEBase through the FEIn...
std::vector< unsigned int > _base_node_index
The internal structure of the InfFE – tensor product of base element times radial nodes – has to be...
Definition: inf_fe.h:702
ElemType
Defines an enum for geometric element types.
void init_shape_functions(const std::vector< Point > &radial_qp, const std::vector< Point > &base_qp, const Elem *inf_elem)
Initialize all the data fields like weight, mode, phi, dphidxi, dphideta, dphidzeta, etc.
Definition: inf_fe.C:386
unsigned short int side
Definition: xdr_io.C:49
Base()
Never use an object of this type.
Definition: inf_fe.h:192
virtual unsigned int n_shape_functions() const libmesh_override
Definition: inf_fe.h:459
unsigned int _n_total_qp
The total number of quadrature points for the current configuration.
Definition: inf_fe.h:739
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
void init_radial_shape_functions(const Elem *inf_elem, const std::vector< Point > *radial_pts=libmesh_nullptr)
Some of the member data only depend on the radial part of the infinite element.
Definition: inf_fe.C:302
const class libmesh_nullptr_t libmesh_nullptr
static Real shape(const FEType &fet, const ElemType t, const unsigned int i, const Point &p)
static unsigned int n_dofs_at_node(const Order o_radial, const unsigned int n_onion)
static const Real TOLERANCE
static unsigned int n_shape_functions(const FEType &fet, const ElemType t)
Definition: inf_fe.h:301
The libMesh namespace provides an interface to certain functionality in the library.
static void compute_node_indices(const ElemType inf_elem_type, const unsigned int outer_node_index, unsigned int &base_node, unsigned int &radial_node)
Computes the indices in the base base_node and in radial direction radial_node (either 0 or 1) associ...
std::vector< std::vector< Real > > mode
the radial approximation shapes in local coordinates Needed when setting up the overall shape functio...
Definition: inf_fe.h:639
UniquePtr< QBase > base_qrule
The quadrature rule for the base element associated with the current infinite element.
Definition: inf_fe.h:751
UniquePtr< Elem > base_elem
The base element associated with the current infinite element.
Definition: inf_fe.h:763
libmesh_assert(j)
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
Radial()
Never use an object of this type.
Definition: inf_fe.h:100
friend class InfFE
Make all InfFE<Dim,T_radial,T_map> classes friends of each other, so that the protected eval() may be...
Definition: inf_fe.h:817
static unsigned int n_dofs_per_elem(const Order o_radial)
Definition: inf_fe.h:170
static Real decay(const Real v)
Definition: inf_fe.h:826
virtual void reinit(const Elem *elem, const std::vector< Point > *const pts=libmesh_nullptr, const std::vector< Real > *const weights=libmesh_nullptr) libmesh_override
This is at the core of this class.
Definition: inf_fe.C:107
std::vector< Real > dphasedzeta
the third local derivative (for 3D, the derivative in radial direction) of the phase term in local co...
Definition: inf_fe.h:677
static Real eval(Real v, Order o_radial, unsigned int i)
FEType current_fe_type
This FEType stores the characteristics for which the data structures phi, phi_map etc are currently i...
Definition: inf_fe.h:781
virtual void compute_shape_functions(const Elem *, const std::vector< Point > &) libmesh_override
After having updated the jacobian and the transformation from local to global coordinates in FEAbstra...
Definition: inf_fe.C:860
std::vector< Real > dphasedeta
the second local derivative (for 3D, the second in the base) of the phase term in local coordinates...
Definition: inf_fe.h:670
static Real D_deriv(const Real v)
Definition: inf_fe.h:125
static unsigned int n_dofs(const Order o_radial)
Definition: inf_fe.h:147
FEGenericBase< Real > FEBase
FEContinuity
defines an enum for finite element types to libmesh_assert a certain level (or type? Hcurl?) of continuity.
unsigned int _n_total_approx_sf
The number of total approximation shape functions for the current configuration.
Definition: inf_fe.h:733
std::vector< unsigned int > _radial_node_index
The internal structure of the InfFE – tensor product of base element times radial nodes – has to be...
Definition: inf_fe.h:692
static Point inverse_map(const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: inf_fe_map.C:89
std::vector< Real > som
the radial decay in local coordinates.
Definition: inf_fe.h:628
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static void compute_shape_indices(const FEType &fet, const ElemType inf_elem_type, const unsigned int i, unsigned int &base_shape, unsigned int &radial_shape)
Computes the indices of shape functions in the base base_shape and in radial direction radial_shape (...
static Point map(const Elem *inf_elem, const Point &reference_point)
Definition: inf_fe_map.C:40
std::vector< unsigned int > _radial_shape_index
The internal structure of the InfFE – tensor product of base element shapes times radial shapes – h...
Definition: inf_fe.h:712
Order
defines an enum for polynomial orders.
Definition: enum_order.h:32
Infinite elements are in some sense directional, compared to conventional finite elements.
Definition: inf_fe.h:93
static bool _warned_for_nodal_soln
static members that are used to issue warning messages only once.
Definition: inf_fe.h:806
static void compute_data(const FEType &fe_t, const Elem *inf_elem, FEComputeData &data)
Generalized version of shape(), takes an Elem *.
UniquePtr< FEBase > base_fe
Have a FE<Dim-1,T_base> handy for base approximation.
Definition: inf_fe.h:771
static void nodal_soln(const FEType &fet, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln)
Usually, this method would build the nodal soln from the element soln.
IterBase * data
Ideally this private member data should have protected access.
std::vector< unsigned int > _base_shape_index
The internal structure of the InfFE – tensor product of base element shapes times radial shapes – h...
Definition: inf_fe.h:722
virtual void side_map(const Elem *, const Elem *, const unsigned int, const std::vector< Point > &, std::vector< Point > &) libmesh_override
Computes the reference space quadrature points on the side of an element based on the side quadrature...
Definition: inf_fe.h:433
std::vector< Real > dist
the radial distance of the base nodes from the origin
Definition: inf_fe.h:611
virtual void edge_reinit(const Elem *elem, const unsigned int edge, const Real tolerance=TOLERANCE, const std::vector< Point > *const pts=libmesh_nullptr, const std::vector< Real > *const weights=libmesh_nullptr) libmesh_override
Not implemented yet.
static Order mapping_order()
Definition: inf_fe.h:131
virtual FEContinuity get_continuity() const libmesh_override
Definition: inf_fe.h:332
static ElemType _compute_node_indices_fast_current_elem_type
When compute_node_indices_fast() is used, this static variable remembers the element type for which t...
Definition: inf_fe.h:798
std::vector< std::vector< Real > > dmodedv
the first local derivative of the radial approximation shapes.
Definition: inf_fe.h:645
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
The QBase class provides the basic functionality from which various quadrature rules can be derived...
Definition: quadrature.h:53
virtual void attach_quadrature_rule(QBase *q) libmesh_override
The use of quadrature rules with the InfFE class is somewhat different from the approach of the FE cl...
Definition: inf_fe.C:67
std::vector< Real > _total_qrule_weights
this vector contains the combined integration weights, so that FEAbstract::compute_map() can still be...
Definition: inf_fe.h:745
void update_base_elem(const Elem *inf_elem)
Updates the protected member base_elem to the appropriate base element for the given inf_elem...
Definition: inf_fe.C:96