libMesh
hp_coarsentest.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_HP_COARSENTEST_H
21 #define LIBMESH_HP_COARSENTEST_H
22 
23 // Local Includes
24 #include "libmesh/auto_ptr.h"
25 #include "libmesh/dense_matrix.h"
26 #include "libmesh/dense_vector.h"
27 #include "libmesh/hp_selector.h"
28 #include "libmesh/id_types.h"
29 #include "libmesh/libmesh_common.h"
30 
31 #include "libmesh/fe.h" // MipsPro requires fe.h and quadrature.h in order to
32 #include "libmesh/quadrature.h" // delete UniquePtrs<> upon destruction
33 
34 // C++ includes
35 #include <vector>
36 
37 #ifdef LIBMESH_ENABLE_AMR
38 
39 namespace libMesh
40 {
41 
42 // Forward Declarations
43 class Elem;
44 class Point;
45 class System;
46 template <typename T> class TensorValue;
47 template <typename T> class VectorValue;
50 typedef RealVectorValue RealGradient;
51 typedef RealTensorValue RealTensor;
52 
53 
68 class HPCoarsenTest : public HPSelector
69 {
70 public:
71 
76  {
77  libmesh_experimental();
78  }
79 
83  virtual ~HPCoarsenTest() {}
84 
85 
92  virtual void select_refinement (System & system) libmesh_override;
93 
99 
100 protected:
105  void add_projection(const System &, const Elem *, unsigned int var);
106 
110  const Elem * coarse;
111 
115  std::vector<dof_id_type> dof_indices;
116 
121 
125  const std::vector<std::vector<Real>> * phi, * phi_coarse;
126  const std::vector<std::vector<RealGradient>> * dphi, * dphi_coarse;
127  const std::vector<std::vector<RealTensor>> * d2phi, * d2phi_coarse;
128 
132  const std::vector<Real> * JxW;
133 
137  const std::vector<Point> * xyz_values;
138  std::vector<Point> coarse_qpoints;
139 
144 
156 };
157 
158 } // namespace libMesh
159 
160 #endif // #ifdef LIBMESH_ENABLE_AMR
161 
162 #endif // LIBMESH_HP_COARSENTEST_H
RealVectorValue RealGradient
const std::vector< std::vector< Real > > * phi
The shape functions and their derivatives.
virtual void select_refinement(System &system) libmesh_override
This pure virtual function must be redefined in derived classes to take a mesh flagged for h refineme...
void add_projection(const System &, const Elem *, unsigned int var)
The helper function which adds individual fine element data to the coarse element projection...
virtual ~HPCoarsenTest()
Destructor.
DenseVector< Number > Fe
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
const std::vector< std::vector< Real > > * phi_coarse
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
const std::vector< std::vector< RealGradient > > * dphi_coarse
RealTensorValue RealTensor
UniquePtr< QBase > qrule
The quadrature rule for the fine element.
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.
UniquePtr< FEBase > fe
The finite element objects for fine and coarse elements.
const std::vector< Point > * xyz_values
Quadrature locations.
DenseVector< Number > Uc
Coefficients for projected coarse and projected p-derefined solutions.
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
const std::vector< Real > * JxW
Mapping jacobians.
TensorValue< Real > RealTensorValue
Useful typedefs to allow transparent switching between Real and Complex data types.
HPCoarsenTest()
Constructor.
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
DenseVector< Number > Up
UniquePtr< FEBase > fe_coarse
This class uses the error estimate given by different types of derefinement (h coarsening or p reduct...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const Elem * coarse
The coarse element on which a solution projection is cached.
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
This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space.
DenseMatrix< Number > Ke
Linear system for projections.