fe_lagrange_shape_1D.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2018 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 
114 
115 
116 template <>
118  const Order order,
119  const unsigned int i,
120  const Point & p)
121 {
122  libmesh_assert(elem);
123 
124  return FE<1,LAGRANGE>::shape(elem->type(), static_cast<Order>(order + elem->p_level()), i, p);
125 }
126 
127 
128 
129 template <>
131  const Order order,
132  const unsigned int i,
133  const unsigned int libmesh_dbg_var(j),
134  const Point & p)
135 {
136  // only d()/dxi in 1D!
137 
138  libmesh_assert_equal_to (j, 0);
139 
140  const Real xi = p(0);
141 
142 
143  switch (order)
144  {
145  // Lagrange linear shape function derivatives
146  case FIRST:
147  {
148  libmesh_assert_less (i, 2);
149 
150  switch (i)
151  {
152  case 0:
153  return -.5;
154 
155  case 1:
156  return .5;
157 
158  default:
159  libmesh_error_msg("Invalid shape function index i = " << i);
160  }
161  }
162 
163 
164  // Lagrange quadratic shape function derivatives
165  case SECOND:
166  {
167  libmesh_assert_less (i, 3);
168 
169  switch (i)
170  {
171  case 0:
172  return xi-.5;
173 
174  case 1:
175  return xi+.5;
176 
177  case 2:
178  return -2.*xi;
179 
180  default:
181  libmesh_error_msg("Invalid shape function index i = " << i);
182  }
183  }
184 
185 
186  // Lagrange cubic shape function derivatives
187  case THIRD:
188  {
189  libmesh_assert_less (i, 4);
190 
191  switch (i)
192  {
193  case 0:
194  return -9./16.*(3.*xi*xi-2.*xi-1./9.);
195 
196  case 1:
197  return -9./16.*(-3.*xi*xi-2.*xi+1./9.);
198 
199  case 2:
200  return 27./16.*(3.*xi*xi-2./3.*xi-1.);
201 
202  case 3:
203  return 27./16.*(-3.*xi*xi-2./3.*xi+1.);
204 
205  default:
206  libmesh_error_msg("Invalid shape function index i = " << i);
207  }
208  }
209 
210 
211  default:
212  libmesh_error_msg("ERROR: Unsupported polynomial order = " << order);
213  }
214 }
215 
216 
217 
218 template <>
220  const Order order,
221  const unsigned int i,
222  const unsigned int j,
223  const Point & p)
224 {
225  libmesh_assert(elem);
226 
227  return FE<1,LAGRANGE>::shape_deriv(elem->type(),
228  static_cast<Order>(order + elem->p_level()), i, j, p);
229 }
230 
231 
232 
233 
234 template <>
236  const Order order,
237  const unsigned int i,
238  const unsigned int libmesh_dbg_var(j),
239  const Point & p)
240 {
241  // Don't need to switch on j. 1D shape functions
242  // depend on xi only!
243 
244  const Real xi = p(0);
245  libmesh_assert_equal_to (j, 0);
246 
247  switch (order)
248  {
249  // linear Lagrange shape functions
250  case FIRST:
251  {
252  // All second derivatives of linears are zero....
253  return 0.;
254  }
255 
256  // quadratic Lagrange shape functions
257  case SECOND:
258  {
259  switch (i)
260  {
261  case 0:
262  return 1.;
263 
264  case 1:
265  return 1.;
266 
267  case 2:
268  return -2.;
269 
270  default:
271  libmesh_error_msg("Invalid shape function index i = " << i);
272  }
273  } // end case SECOND
274 
275  case THIRD:
276  {
277  switch (i)
278  {
279  case 0:
280  return -9./16.*(6.*xi-2);
281 
282  case 1:
283  return -9./16.*(-6*xi-2.);
284 
285  case 2:
286  return 27./16.*(6*xi-2./3.);
287 
288  case 3:
289  return 27./16.*(-6*xi-2./3.);
290 
291  default:
292  libmesh_error_msg("Invalid shape function index i = " << i);
293  }
294  } // end case THIRD
295 
296 
297  default:
298  libmesh_error_msg("ERROR: Unsupported polynomial order = " << order);
299  } // end switch (order)
300 }
301 
302 
303 
304 template <>
306  const Order order,
307  const unsigned int i,
308  const unsigned int j,
309  const Point & p)
310 {
311  libmesh_assert(elem);
312 
314  static_cast<Order>(order + elem->p_level()), i, j, p);
315 }
316 
317 } // namespace libMesh
static OutputShape shape(const ElemType t, const Order o, const unsigned int i, const Point &p)
The base class for all geometric element types.
Definition: elem.h:100
static OutputShape shape_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
unsigned int p_level() const
Definition: elem.h:2555
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual ElemType type() const =0
A geometric point in (x,y,z) 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)