fe_clough.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 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 
25 
26 namespace libMesh
27 {
28 
29  // ------------------------------------------------------------
30  // Clough-specific implementations
31 
32  // Anonymous namespace for local helper functions
33  namespace {
34 
35  void clough_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, CLOUGH);
51 
52  switch (totalorder)
53  {
54  // Piecewise cubic shape functions with linear flux edges
55  case SECOND:
56  // Piecewise cubic shape functions
57  case THIRD:
58  {
59 
60  const unsigned int n_sf =
61  // FE<Dim,T>::n_shape_functions(elem_type, totalorder);
62  FEInterface::n_shape_functions(Dim, fe_type, elem_type);
63 
64  std::vector<Point> refspace_nodes;
65  FEBase::get_refspace_nodes(elem_type,refspace_nodes);
66  libmesh_assert_equal_to (refspace_nodes.size(), n_nodes);
67 
68  for (unsigned int n=0; n<n_nodes; n++)
69  {
70  libmesh_assert_equal_to (elem_soln.size(), n_sf);
71 
72  // Zero before summation
73  nodal_soln[n] = 0;
74 
75  // u_i = Sum (alpha_i phi_i)
76  for (unsigned int i=0; i<n_sf; i++)
77  nodal_soln[n] += elem_soln[i] *
78  // FE<Dim,T>::shape(elem, order, i, mapped_point);
79  FEInterface::shape(Dim, fe_type, elem, i, refspace_nodes[n]);
80  }
81 
82  return;
83  }
84 
85  default:
86  {
87  libmesh_error();
88  }
89  }
90  } // clough_nodal_soln()
91 
92 
93 
94 
95  unsigned int clough_n_dofs(const ElemType t, const Order o)
96  {
97  switch (o)
98  {
99  // Piecewise cubic shape functions with linear flux edges
100  case SECOND:
101  {
102  switch (t)
103  {
104  case TRI6:
105  return 9;
106 
107  default:
108  {
109 #ifdef DEBUG
110  libMesh::err << "ERROR: Bad ElemType = " << t
111  << " for " << o << "th order approximation!"
112  << std::endl;
113 #endif
114  libmesh_error();
115  }
116  }
117  }
118  // Piecewise cubic Clough-Tocher element
119  case THIRD:
120  {
121  switch (t)
122  {
123  case TRI6:
124  return 12;
125 
126  default:
127  {
128 #ifdef DEBUG
129  libMesh::err << "ERROR: Bad ElemType = " << t
130  << " for " << o << "th order approximation!"
131  << std::endl;
132 #endif
133  libmesh_error();
134  }
135  }
136  }
137 
138  default:
139  {
140  libmesh_error();
141  }
142  }
143 
144  libmesh_error();
145  return 0;
146  } // clough_n_dofs()
147 
148 
149 
150 
151 
152  unsigned int clough_n_dofs_at_node(const ElemType t,
153  const Order o,
154  const unsigned int n)
155  {
156  switch (o)
157  {
158  // Piecewise cubic shape functions with linear flux edges
159  case SECOND:
160  {
161  switch (t)
162  {
163  // The 2D Clough-Tocher defined on a 6-noded triangle
164  case TRI6:
165  {
166  switch (n)
167  {
168  case 0:
169  case 1:
170  case 2:
171  return 3;
172 
173  case 3:
174  case 4:
175  case 5:
176  return 0;
177 
178  default:
179  libmesh_error();
180  }
181  }
182 
183  default:
184  {
185 #ifdef DEBUG
186  libMesh::err << "ERROR: Bad ElemType = " << t
187  << " for " << o << "th order approximation!"
188  << std::endl;
189 #endif
190  libmesh_error();
191  }
192 
193  }
194  }
195  // The third-order clough shape functions
196  case THIRD:
197  {
198  switch (t)
199  {
200  // The 2D Clough-Tocher defined on a 6-noded triangle
201  case TRI6:
202  {
203  switch (n)
204  {
205  case 0:
206  case 1:
207  case 2:
208  return 3;
209 
210  case 3:
211  case 4:
212  case 5:
213  return 1;
214 
215  default:
216  libmesh_error();
217  }
218  }
219 
220  default:
221  {
222 #ifdef DEBUG
223  libMesh::err << "ERROR: Bad ElemType = " << t
224  << " for " << o << "th order approximation!"
225  << std::endl;
226 #endif
227  libmesh_error();
228  }
229 
230  }
231  }
232  default:
233  {
234  libmesh_error();
235  }
236  }
237 
238  libmesh_error();
239 
240  return 0;
241  } // clough_n_dofs_at_node()
242 
243 
244 
245 
246  unsigned int clough_n_dofs_per_elem(const ElemType t, const Order o)
247  {
248  switch (o)
249  {
250  // Piecewise cubic shape functions with linear flux edges
251  case SECOND:
252  // The third-order Clough-Tocher shape functions
253  case THIRD:
254  {
255  switch (t)
256  {
257  // The 2D clough defined on a 6-noded triangle
258  case TRI6:
259  return 0;
260 
261  default:
262  {
263 #ifdef DEBUG
264  libMesh::err << "ERROR: Bad ElemType = " << t
265  << " for " << o << "th order approximation!"
266  << std::endl;
267 #endif
268  libmesh_error();
269  }
270 
271  }
272  }
273  // Otherwise no DOFS per element
274  default:
275  libmesh_error();
276  return 0;
277  }
278  } // clough_n_dofs_per_elem()
279 
280  } // anonymous
281 
282 
283 
284 
285  // Do full-specialization of nodal_soln() function for every
286  // dimension, instead of explicit instantiation at the end of this
287  // file.
288  // This could be macro-ified so that it fits on one line...
289  template <>
290  void FE<0,CLOUGH>::nodal_soln(const Elem* elem,
291  const Order order,
292  const std::vector<Number>& elem_soln,
293  std::vector<Number>& nodal_soln)
294  { clough_nodal_soln(elem, order, elem_soln, nodal_soln, /*Dim=*/0); }
295 
296  template <>
297  void FE<1,CLOUGH>::nodal_soln(const Elem* elem,
298  const Order order,
299  const std::vector<Number>& elem_soln,
300  std::vector<Number>& nodal_soln)
301  { clough_nodal_soln(elem, order, elem_soln, nodal_soln, /*Dim=*/1); }
302 
303  template <>
304  void FE<2,CLOUGH>::nodal_soln(const Elem* elem,
305  const Order order,
306  const std::vector<Number>& elem_soln,
307  std::vector<Number>& nodal_soln)
308  { clough_nodal_soln(elem, order, elem_soln, nodal_soln, /*Dim=*/2); }
309 
310  template <>
311  void FE<3,CLOUGH>::nodal_soln(const Elem* elem,
312  const Order order,
313  const std::vector<Number>& elem_soln,
314  std::vector<Number>& nodal_soln)
315  { clough_nodal_soln(elem, order, elem_soln, nodal_soln, /*Dim=*/3); }
316 
317 
318  // Full specialization of n_dofs() function for every dimension
319  template <> unsigned int FE<0,CLOUGH>::n_dofs(const ElemType t, const Order o) { return clough_n_dofs(t, o); }
320  template <> unsigned int FE<1,CLOUGH>::n_dofs(const ElemType t, const Order o) { return clough_n_dofs(t, o); }
321  template <> unsigned int FE<2,CLOUGH>::n_dofs(const ElemType t, const Order o) { return clough_n_dofs(t, o); }
322  template <> unsigned int FE<3,CLOUGH>::n_dofs(const ElemType t, const Order o) { return clough_n_dofs(t, o); }
323 
324 
325  // Full specialization of n_dofs_at_node() function for every dimension.
326  template <> unsigned int FE<0,CLOUGH>::n_dofs_at_node(const ElemType t, const Order o, const unsigned int n) { return clough_n_dofs_at_node(t, o, n); }
327  template <> unsigned int FE<1,CLOUGH>::n_dofs_at_node(const ElemType t, const Order o, const unsigned int n) { return clough_n_dofs_at_node(t, o, n); }
328  template <> unsigned int FE<2,CLOUGH>::n_dofs_at_node(const ElemType t, const Order o, const unsigned int n) { return clough_n_dofs_at_node(t, o, n); }
329  template <> unsigned int FE<3,CLOUGH>::n_dofs_at_node(const ElemType t, const Order o, const unsigned int n) { return clough_n_dofs_at_node(t, o, n); }
330 
331  // Full specialization of n_dofs_per_elem() function for every dimension.
332  template <> unsigned int FE<0,CLOUGH>::n_dofs_per_elem(const ElemType t, const Order o) { return clough_n_dofs_per_elem(t, o); }
333  template <> unsigned int FE<1,CLOUGH>::n_dofs_per_elem(const ElemType t, const Order o) { return clough_n_dofs_per_elem(t, o); }
334  template <> unsigned int FE<2,CLOUGH>::n_dofs_per_elem(const ElemType t, const Order o) { return clough_n_dofs_per_elem(t, o); }
335  template <> unsigned int FE<3,CLOUGH>::n_dofs_per_elem(const ElemType t, const Order o) { return clough_n_dofs_per_elem(t, o); }
336 
337  // Clough FEMs are C^1 continuous
338  template <> FEContinuity FE<0,CLOUGH>::get_continuity() const { return C_ONE; }
339  template <> FEContinuity FE<1,CLOUGH>::get_continuity() const { return C_ONE; }
340  template <> FEContinuity FE<2,CLOUGH>::get_continuity() const { return C_ONE; }
341  template <> FEContinuity FE<3,CLOUGH>::get_continuity() const { return C_ONE; }
342 
343  // Clough FEMs are not (currently) hierarchic
344  template <> bool FE<0,CLOUGH>::is_hierarchic() const { return false; } // FIXME - this will be changed
345  template <> bool FE<1,CLOUGH>::is_hierarchic() const { return false; } // FIXME - this will be changed
346  template <> bool FE<2,CLOUGH>::is_hierarchic() const { return false; } // FIXME - this will be changed
347  template <> bool FE<3,CLOUGH>::is_hierarchic() const { return false; } // FIXME - this will be changed
348 
349 #ifdef LIBMESH_ENABLE_AMR
350  // compute_constraints() specializations are only needed for 2 and 3D
351  template <>
353  DofMap &dof_map,
354  const unsigned int variable_number,
355  const Elem* elem)
356  { compute_proj_constraints(constraints, dof_map, variable_number, elem); }
357 
358  template <>
360  DofMap &dof_map,
361  const unsigned int variable_number,
362  const Elem* elem)
363  { compute_proj_constraints(constraints, dof_map, variable_number, elem); }
364 #endif // #ifdef LIBMESH_ENABLE_AMR
365 
366  // Clough FEM shapes need reinit
367  template <> bool FE<0,CLOUGH>::shapes_need_reinit() const { return true; }
368  template <> bool FE<1,CLOUGH>::shapes_need_reinit() const { return true; }
369  template <> bool FE<2,CLOUGH>::shapes_need_reinit() const { return true; }
370  template <> bool FE<3,CLOUGH>::shapes_need_reinit() const { return true; }
371 
372 } // namespace libMesh

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:05 UTC

Hosted By:
SourceForge.net Logo