libMesh
fe_l2_lagrange_shape_1D.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 // C++ includes
20 
21 // Local includes
22 #include "libmesh/fe.h"
23 #include "libmesh/elem.h"
24 
25 namespace libMesh
26 {
27 
28 
29 
30 
31 template <>
33  const Order order,
34  const unsigned int i,
35  const Point & p)
36 {
37  const Real xi = p(0);
38 
39 
40  switch (order)
41  {
42  // Lagrange linears
43  case FIRST:
44  {
45  libmesh_assert_less (i, 2);
46 
47  switch (i)
48  {
49  case 0:
50  return .5*(1. - xi);
51 
52  case 1:
53  return .5*(1. + xi);
54 
55  default:
56  libmesh_error_msg("Invalid shape function index i = " << i);
57  }
58  }
59 
60 
61 
62  // Lagrange quadratics
63  case SECOND:
64  {
65  libmesh_assert_less (i, 3);
66 
67  switch (i)
68  {
69  case 0:
70  return .5*xi*(xi - 1.);
71 
72  case 1:
73  return .5*xi*(xi + 1);
74 
75  case 2:
76  return (1. - xi*xi);
77 
78  default:
79  libmesh_error_msg("Invalid shape function index i = " << i);
80  }
81  }
82 
83 
84 
85  // Lagrange cubics
86  case THIRD:
87  {
88  libmesh_assert_less (i, 4);
89 
90  switch (i)
91  {
92  case 0:
93  return 9./16.*(1./9.-xi*xi)*(xi-1.);
94 
95  case 1:
96  return -9./16.*(1./9.-xi*xi)*(xi+1.);
97 
98  case 2:
99  return 27./16.*(1.-xi*xi)*(1./3.-xi);
100 
101  case 3:
102  return 27./16.*(1.-xi*xi)*(1./3.+xi);
103 
104  default:
105  libmesh_error_msg("Invalid shape function index i = " << i);
106  }
107  }
108 
109  default:
110  libmesh_error_msg("ERROR: Unsupported polynomial order = " << order);
111  }
112 
113  libmesh_error_msg("We'll never get here!");
114  return 0.;
115 }
116 
117 
118 
119 template <>
121  const Order order,
122  const unsigned int i,
123  const Point & p)
124 {
125  libmesh_assert(elem);
126 
127  return FE<1,L2_LAGRANGE>::shape(elem->type(), static_cast<Order>(order + elem->p_level()), i, p);
128 }
129 
130 
131 
132 template <>
134  const Order order,
135  const unsigned int i,
136  const unsigned int libmesh_dbg_var(j),
137  const Point & p)
138 {
139  // only d()/dxi in 1D!
140 
141  libmesh_assert_equal_to (j, 0);
142 
143  const Real xi = p(0);
144 
145 
146  switch (order)
147  {
148  // Lagrange linear shape function derivatives
149  case FIRST:
150  {
151  libmesh_assert_less (i, 2);
152 
153  switch (i)
154  {
155  case 0:
156  return -.5;
157 
158  case 1:
159  return .5;
160 
161  default:
162  libmesh_error_msg("Invalid shape function index i = " << i);
163  }
164  }
165 
166 
167  // Lagrange quadratic shape function derivatives
168  case SECOND:
169  {
170  libmesh_assert_less (i, 3);
171 
172  switch (i)
173  {
174  case 0:
175  return xi-.5;
176 
177  case 1:
178  return xi+.5;
179 
180  case 2:
181  return -2.*xi;
182 
183  default:
184  libmesh_error_msg("Invalid shape function index i = " << i);
185  }
186  }
187 
188 
189  // Lagrange cubic shape function derivatives
190  case THIRD:
191  {
192  libmesh_assert_less (i, 4);
193 
194  switch (i)
195  {
196  case 0:
197  return -9./16.*(3.*xi*xi-2.*xi-1./9.);
198 
199  case 1:
200  return -9./16.*(-3.*xi*xi-2.*xi+1./9.);
201 
202  case 2:
203  return 27./16.*(3.*xi*xi-2./3.*xi-1.);
204 
205  case 3:
206  return 27./16.*(-3.*xi*xi-2./3.*xi+1.);
207 
208  default:
209  libmesh_error_msg("Invalid shape function index i = " << i);
210  }
211  }
212 
213 
214  default:
215  libmesh_error_msg("ERROR: Unsupported polynomial order = " << order);
216  }
217 
218  libmesh_error_msg("We'll never get here!");
219  return 0.;
220 }
221 
222 
223 
224 template <>
226  const Order order,
227  const unsigned int i,
228  const unsigned int j,
229  const Point & p)
230 {
231  libmesh_assert(elem);
232 
233  return FE<1,L2_LAGRANGE>::shape_deriv(elem->type(),
234  static_cast<Order>(order + elem->p_level()), i, j, p);
235 }
236 
237 
238 
239 
240 template <>
242  const Order order,
243  const unsigned int i,
244  const unsigned int libmesh_dbg_var(j),
245  const Point & p)
246 {
247  // Don't need to switch on j. 1D shape functions
248  // depend on xi only!
249 
250  const Real xi = p(0);
251  libmesh_assert_equal_to (j, 0);
252 
253  switch (order)
254  {
255  // linear Lagrange shape functions
256  case FIRST:
257  {
258  // All second derivatives of linears are zero....
259  return 0.;
260  }
261 
262  // quadratic Lagrange shape functions
263  case SECOND:
264  {
265  switch (i)
266  {
267  case 0:
268  return 1.;
269 
270  case 1:
271  return 1.;
272 
273  case 2:
274  return -2.;
275 
276  default:
277  libmesh_error_msg("Invalid shape function index i = " << i);
278  }
279  } // end case SECOND
280 
281  case THIRD:
282  {
283  switch (i)
284  {
285  case 0:
286  return -9./16.*(6.*xi-2);
287 
288  case 1:
289  return -9./16.*(-6*xi-2.);
290 
291  case 2:
292  return 27./16.*(6*xi-2./3.);
293 
294  case 3:
295  return 27./16.*(-6*xi-2./3.);
296 
297  default:
298  libmesh_error_msg("Invalid shape function index i = " << i);
299  }
300  } // end case THIRD
301 
302 
303  default:
304  libmesh_error_msg("ERROR: Unsupported polynomial order = " << order);
305  } // end switch (order)
306 
307  libmesh_error_msg("We'll never get here!");
308  return 0.;
309 }
310 
311 
312 
313 template <>
315  const Order order,
316  const unsigned int i,
317  const unsigned int j,
318  const Point & p)
319 {
320  libmesh_assert(elem);
321 
323  static_cast<Order>(order + elem->p_level()), i, j, p);
324 }
325 
326 } // namespace libMesh
static OutputShape shape(const ElemType t, const Order o, const unsigned int i, const Point &p)
virtual ElemType type() const =0
unsigned int p_level() const
Definition: elem.h:2422
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 OutputShape shape_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
The libMesh namespace provides an interface to certain functionality in the library.
libmesh_assert(j)
PetscErrorCode Vec Mat libmesh_dbg_var(j)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Order
defines an enum for polynomial orders.
Definition: enum_order.h:32
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
static OutputShape shape_second_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)