1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2020 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_QUADRATURE_COMPOSITE_H
21 #define LIBMESH_QUADRATURE_COMPOSITE_H
22 
23 #include "libmesh/libmesh_config.h"
24 
25 #if defined(LIBMESH_HAVE_TRIANGLE) && defined(LIBMESH_HAVE_TETGEN)
26 
27 // Local includes
28 #include "libmesh/quadrature.h"
29 #include "libmesh/elem_cutter.h"
30 #include "libmesh/fe_base.h"
31 
32 // C++ includes
33 #include <memory>
34 
35 namespace libMesh
36 {
37 
38 /**
39  * This class implements generic composite quadrature rules.
40  * Composite quadrature rules are constructed from any of the
41  * supported rules by breaking an element into subelements and
42  * applying the base rule on each subelement.  This class uses the
43  * ElemCutter, which is only available if libmesh is configured with
44  * --disable-strict-lgpl.
45  *
46  * \author Benjamin Kirk
47  * \date 2013
48  * \brief A quadrature rule for subdivided elements.
49  */
50 template <class QSubCell>
51 class QComposite final : public QSubCell
52 {
53 public:
54 
55   /**
56    * Import base class data into namespace scope.
57    */
58   using QSubCell::_dim;
59   using QSubCell::_points;
60   using QSubCell::_weights;
61   using QSubCell::init;
62 
63   /**
64    * Constructor.  Declares the order of the quadrature rule.
65    */
66   QComposite (unsigned int dim,
67               Order order=INVALID_ORDER);
68 
69   /**
70    * This class contains a unique_ptr member, so it can't be default
71    * copy constructed or assigned.
72    */
73   QComposite (const QComposite &) = delete;
74   QComposite & operator= (const QComposite &) = delete;
75 
76   /**
77    * Copy/move ctor, copy/move assignment operator, and destructor are
78    * all explicitly defaulted for this simple class.
79    */
80   QComposite (QComposite &&) = default;
81   QComposite & operator= (QComposite &&) = default;
82   virtual ~QComposite() = default;
83 
84   /**
85    * \returns \p QCOMPOSITE.
86    */
87   virtual QuadratureType type() const override;
88 
89   /**
90    * Overrides the base class init() function, and uses the ElemCutter to
91    * subdivide the element into "inside" and "outside" subelements.
92    */
93   virtual void init (const Elem & elem,
94                      const std::vector<Real> & vertex_distance_func,
95                      unsigned int p_level=0) override;
96 
97 private:
98 
99   /**
100    * Helper function called from init() to collect all the points and
101    * weights of the subelement quadrature rules.
102    */
103   void add_subelem_values (const std::vector<Elem const *> & subelem);
104 
105   /**
106    * Subcell quadrature object.
107    */
108   QSubCell _q_subcell;
109 
110   /**
111    * ElemCutter object.
112    */
113   ElemCutter _elem_cutter;
114 
115   /**
116    * Lagrange FE to use for subcell mapping.
117    */
118   std::unique_ptr<FEBase> _lagrange_fe;
119 };
120 
121 } // namespace libMesh
122 
123 #endif // LIBMESH_HAVE_TRIANGLE && LIBMESH_HAVE_TETGEN
124 #endif // LIBMESH_QUADRATURE_COMPOSITE_H
125