libMesh
fe_l2_hierarchic.C
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 // Local includes
21 #include "libmesh/elem.h"
22 #include "libmesh/fe.h"
23 #include "libmesh/fe_interface.h"
24 #include "libmesh/string_to_enum.h"
25 
26 namespace libMesh
27 {
28 
29 // ------------------------------------------------------------
30 // Hierarchic-specific implementations
31 
32 // Anonymous namespace for local helper functions
33 namespace {
34 
35 void l2_hierarchic_nodal_soln(const Elem * elem,
36  const Order order,
37  const std::vector<Number> & elem_soln,
38  std::vector<Number> & nodal_soln,
39  unsigned Dim)
40 {
41  const unsigned int n_nodes = elem->n_nodes();
42 
43  const ElemType elem_type = elem->type();
44 
45  nodal_soln.resize(n_nodes);
46 
47  const Order totalorder = static_cast<Order>(order + elem->p_level());
48 
49  // FEType object to be passed to various FEInterface functions below.
50  FEType fe_type(totalorder, L2_HIERARCHIC);
51 
52  switch (totalorder)
53  {
54  // Constant shape functions
55  case CONSTANT:
56  {
57  libmesh_assert_equal_to (elem_soln.size(), 1);
58 
59  const Number val = elem_soln[0];
60 
61  for (unsigned int n=0; n<n_nodes; n++)
62  nodal_soln[n] = val;
63 
64  return;
65  }
66 
67 
68  // For other orders do interpolation at the nodes
69  // explicitly.
70  default:
71  {
72 
73  const unsigned int n_sf =
74  // FE<Dim,T>::n_shape_functions(elem_type, totalorder);
75  FEInterface::n_shape_functions(Dim, fe_type, elem_type);
76 
77  std::vector<Point> refspace_nodes;
78  FEBase::get_refspace_nodes(elem_type,refspace_nodes);
79  libmesh_assert_equal_to (refspace_nodes.size(), n_nodes);
80 
81  for (unsigned int n=0; n<n_nodes; n++)
82  {
83  libmesh_assert_equal_to (elem_soln.size(), n_sf);
84 
85  // Zero before summation
86  nodal_soln[n] = 0;
87 
88  // u_i = Sum (alpha_i phi_i)
89  for (unsigned int i=0; i<n_sf; i++)
90  nodal_soln[n] += elem_soln[i] *
91  // FE<Dim,T>::shape(elem, order, i, mapped_point);
92  FEInterface::shape(Dim, fe_type, elem, i, refspace_nodes[n]);
93  }
94 
95  return;
96  }
97  }
98 } // l2_hierarchic_nodal_soln()
99 
100 
101 
102 
103 unsigned int l2_hierarchic_n_dofs(const ElemType t, const Order o)
104 {
105  libmesh_assert_greater (o, 0);
106  switch (t)
107  {
108  case NODEELEM:
109  return 1;
110  case EDGE2:
111  case EDGE3:
112  return (o+1);
113  case QUAD4:
114  case QUADSHELL4:
115  case QUAD8:
116  case QUADSHELL8:
117  case QUAD9:
118  return ((o+1)*(o+1));
119  case HEX8:
120  case HEX20:
121  case HEX27:
122  return ((o+1)*(o+1)*(o+1));
123  case TRI3:
124  libmesh_assert_less (o, 2);
125  libmesh_fallthrough();
126  case TRI6:
127  return ((o+1)*(o+2)/2);
128  case INVALID_ELEM:
129  return 0;
130  default:
131  libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for L2_HIERARCHIC FE family!");
132  }
133 
134  libmesh_error_msg("We'll never get here!");
135  return 0;
136 } // l2_hierarchic_n_dofs()
137 
138 
139 } // anonymous namespace
140 
141 
142 
143 
144  // Do full-specialization of nodal_soln() function for every
145  // dimension, instead of explicit instantiation at the end of this
146  // file.
147  // This could be macro-ified so that it fits on one line...
148 template <>
150  const Order order,
151  const std::vector<Number> & elem_soln,
152  std::vector<Number> & nodal_soln)
153 { l2_hierarchic_nodal_soln(elem, order, elem_soln, nodal_soln, /*Dim=*/0); }
154 
155 template <>
157  const Order order,
158  const std::vector<Number> & elem_soln,
159  std::vector<Number> & nodal_soln)
160 { l2_hierarchic_nodal_soln(elem, order, elem_soln, nodal_soln, /*Dim=*/1); }
161 
162 template <>
164  const Order order,
165  const std::vector<Number> & elem_soln,
166  std::vector<Number> & nodal_soln)
167 { l2_hierarchic_nodal_soln(elem, order, elem_soln, nodal_soln, /*Dim=*/2); }
168 
169 template <>
171  const Order order,
172  const std::vector<Number> & elem_soln,
173  std::vector<Number> & nodal_soln)
174 { l2_hierarchic_nodal_soln(elem, order, elem_soln, nodal_soln, /*Dim=*/3); }
175 
176 // Full specialization of n_dofs() function for every dimension
177 template <> unsigned int FE<0,L2_HIERARCHIC>::n_dofs(const ElemType t, const Order o) { return l2_hierarchic_n_dofs(t, o); }
178 template <> unsigned int FE<1,L2_HIERARCHIC>::n_dofs(const ElemType t, const Order o) { return l2_hierarchic_n_dofs(t, o); }
179 template <> unsigned int FE<2,L2_HIERARCHIC>::n_dofs(const ElemType t, const Order o) { return l2_hierarchic_n_dofs(t, o); }
180 template <> unsigned int FE<3,L2_HIERARCHIC>::n_dofs(const ElemType t, const Order o) { return l2_hierarchic_n_dofs(t, o); }
181 
182 // Full specialization of n_dofs_at_node() function for every dimension.
183 // Discontinuous L2 elements only have interior nodes
184 template <> unsigned int FE<0,L2_HIERARCHIC>::n_dofs_at_node(const ElemType, const Order, const unsigned int) { return 0; }
185 template <> unsigned int FE<1,L2_HIERARCHIC>::n_dofs_at_node(const ElemType, const Order, const unsigned int) { return 0; }
186 template <> unsigned int FE<2,L2_HIERARCHIC>::n_dofs_at_node(const ElemType, const Order, const unsigned int) { return 0; }
187 template <> unsigned int FE<3,L2_HIERARCHIC>::n_dofs_at_node(const ElemType, const Order, const unsigned int) { return 0; }
188 
189 // Full specialization of n_dofs_per_elem() function for every dimension.
190 template <> unsigned int FE<0,L2_HIERARCHIC>::n_dofs_per_elem(const ElemType t, const Order o) { return l2_hierarchic_n_dofs(t, o); }
191 template <> unsigned int FE<1,L2_HIERARCHIC>::n_dofs_per_elem(const ElemType t, const Order o) { return l2_hierarchic_n_dofs(t, o); }
192 template <> unsigned int FE<2,L2_HIERARCHIC>::n_dofs_per_elem(const ElemType t, const Order o) { return l2_hierarchic_n_dofs(t, o); }
193 template <> unsigned int FE<3,L2_HIERARCHIC>::n_dofs_per_elem(const ElemType t, const Order o) { return l2_hierarchic_n_dofs(t, o); }
194 
195 // L2 Hierarchic FEMs are C^0 continuous
200 
201 // L2 Hierarchic FEMs are hierarchic (duh!)
202 template <> bool FE<0,L2_HIERARCHIC>::is_hierarchic() const { return true; }
203 template <> bool FE<1,L2_HIERARCHIC>::is_hierarchic() const { return true; }
204 template <> bool FE<2,L2_HIERARCHIC>::is_hierarchic() const { return true; }
205 template <> bool FE<3,L2_HIERARCHIC>::is_hierarchic() const { return true; }
206 
207 #ifdef LIBMESH_ENABLE_AMR
208 // compute_constraints() is a NOOP for DISCONTINUOUS FE's
209 template <>
211  DofMap &,
212  const unsigned int,
213  const Elem *)
214 { }
215 
216 template <>
218  DofMap &,
219  const unsigned int,
220  const Elem *)
221 { }
222 #endif // #ifdef LIBMESH_ENABLE_AMR
223 
224 // L2-Hierarchic FEM shapes need reinit
225 template <> bool FE<0,L2_HIERARCHIC>::shapes_need_reinit() const { return true; }
226 template <> bool FE<1,L2_HIERARCHIC>::shapes_need_reinit() const { return true; }
227 template <> bool FE<2,L2_HIERARCHIC>::shapes_need_reinit() const { return true; }
228 template <> bool FE<3,L2_HIERARCHIC>::shapes_need_reinit() const { return true; }
229 
230 } // namespace libMesh
static unsigned int n_dofs(const ElemType t, const Order o)
virtual bool shapes_need_reinit() const libmesh_override
ElemType
Defines an enum for geometric element types.
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
static unsigned int n_dofs_at_node(const ElemType t, const Order o, const unsigned int n)
The libMesh namespace provides an interface to certain functionality in the library.
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:167
const dof_id_type n_nodes
Definition: tecplot_io.C:67
virtual FEContinuity get_continuity() const libmesh_override
static unsigned int n_shape_functions(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:385
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
Definition: fe_interface.C:641
static void get_refspace_nodes(const ElemType t, std::vector< Point > &nodes)
Definition: fe_abstract.C:259
virtual bool is_hierarchic() const libmesh_override
static unsigned int n_dofs_per_elem(const ElemType t, const Order o)
FEContinuity
defines an enum for finite element types to libmesh_assert a certain level (or type? Hcurl?) of continuity.
std::string enum_to_string(const T e)
static void compute_constraints(DofConstraints &constraints, DofMap &dof_map, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for non-conforming adapted meshes) corresponding to var...
static void nodal_soln(const Elem *elem, const Order o, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln)
Build the nodal soln from the element soln.
Order
defines an enum for polynomial orders.
Definition: enum_order.h:32
The constraint matrix storage format.
Definition: dof_map.h:96