libMesh
hp_coarsentest.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2025 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_HP_COARSENTEST_H
21 #define LIBMESH_HP_COARSENTEST_H
22 
23 // Local Includes
24 #include "libmesh/dense_matrix.h"
25 #include "libmesh/dense_vector.h"
26 #include "libmesh/hp_selector.h"
27 #include "libmesh/id_types.h"
28 #include "libmesh/libmesh_common.h"
29 #include "libmesh/fe.h" // MipsPro requires fe.h and quadrature.h
30 #include "libmesh/quadrature.h" // Required for inline deletion std::unique_ptrs<> in destructor
31 
32 // C++ includes
33 #include <vector>
34 #include <memory>
35 
36 #ifdef LIBMESH_ENABLE_AMR
37 
38 namespace libMesh
39 {
40 
41 // Forward Declarations
42 class Elem;
43 class Point;
44 class System;
45 template <typename T> class TensorValue;
46 template <typename T> class VectorValue;
51 
52 
67 class HPCoarsenTest : public HPSelector
68 {
69 public:
70 
75  {
76  libmesh_experimental();
77  }
78 
84  HPCoarsenTest (const HPCoarsenTest &) = delete;
85  HPCoarsenTest & operator= (const HPCoarsenTest &) = delete;
86 
90  HPCoarsenTest (HPCoarsenTest &&) = default;
91  HPCoarsenTest & operator= (HPCoarsenTest &&) = default;
92  virtual ~HPCoarsenTest() = default;
93 
100  virtual void select_refinement (System & system) override;
101 
107 
115  void extra_quadrature_order (const int extraorder)
116  { _extra_order = extraorder; }
117 
118 protected:
123  void add_projection(const System &, const Elem *, unsigned int var);
124 
129 
133  std::vector<dof_id_type> dof_indices;
134 
138  std::unique_ptr<FEBase> fe, fe_coarse;
139 
143  const std::vector<std::vector<Real>> * phi, * phi_coarse;
144  const std::vector<std::vector<RealGradient>> * dphi, * dphi_coarse;
145  const std::vector<std::vector<RealTensor>> * d2phi, * d2phi_coarse;
146 
150  const std::vector<Real> * JxW;
151 
155  const std::vector<Point> * xyz_values;
156  std::vector<Point> coarse_qpoints;
157 
161  std::unique_ptr<QBase> qrule;
162 
174 
179 };
180 
181 } // namespace libMesh
182 
183 #endif // #ifdef LIBMESH_ENABLE_AMR
184 
185 #endif // LIBMESH_HP_COARSENTEST_H
RealVectorValue RealGradient
const std::vector< std::vector< Real > > * phi
The shape functions and their derivatives.
void add_projection(const System &, const Elem *, unsigned int var)
The helper function which adds individual fine element data to the coarse element projection...
DenseVector< Number > Fe
HPCoarsenTest & operator=(const HPCoarsenTest &)=delete
virtual void select_refinement(System &system) override
This pure virtual function must be redefined in derived classes to take a mesh flagged for h refineme...
VectorValue< Real > RealVectorValue
Useful typedefs to allow transparent switching between Real and Complex data types.
Subclasses of this abstract base class choose between h refining and p elevation. ...
Definition: hp_selector.h:48
std::unique_ptr< FEBase > fe_coarse
const std::vector< std::vector< Real > > * phi_coarse
int _extra_order
Extra order to use for quadrature rule.
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
const std::vector< std::vector< RealGradient > > * dphi_coarse
RealTensorValue RealTensor
This class defines a vector in LIBMESH_DIM dimensional Real or Complex space.
The libMesh namespace provides an interface to certain functionality in the library.
const std::vector< Point > * xyz_values
Quadrature locations.
DenseVector< Number > Uc
Coefficients for projected coarse and projected p-derefined solutions.
const std::vector< Real > * JxW
Mapping jacobians.
TensorValue< Real > RealTensorValue
Useful typedefs to allow transparent switching between Real and Complex data types.
virtual ~HPCoarsenTest()=default
HPCoarsenTest()
Constructor.
Manages consistently variables, degrees of freedom, and coefficient vectors.
Definition: system.h:96
DenseVector< Number > Up
This class uses the error estimate given by different types of derefinement (h coarsening or p reduct...
std::unique_ptr< FEBase > fe
The finite element objects for fine and coarse elements.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const std::vector< std::vector< RealTensor > > * d2phi_coarse
const std::vector< std::vector< RealTensor > > * d2phi
const std::vector< std::vector< RealGradient > > * dphi
std::vector< dof_id_type > dof_indices
Global DOF indices for fine elements.
Real p_weight
Because the coarsening test seems to always choose p refinement, we&#39;re providing an option to make h ...
std::vector< Point > coarse_qpoints
Elem * coarse
The coarse element on which a solution projection is cached.
std::unique_ptr< QBase > qrule
The quadrature rule for the fine element.
This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space.
void extra_quadrature_order(const int extraorder)
Increases or decreases the order of the quadrature rule used for numerical integration.
DenseMatrix< Number > Ke
Linear system for projections.