libMesh
fe_monomial_shape_1D.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2024 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 
26 namespace libMesh
27 {
28 
29 
31 
32 
33 template <>
35  const Order libmesh_dbg_var(order),
36  const unsigned int i,
37  const Point & p)
38 {
39  const Real xi = p(0);
40 
41  libmesh_assert_less_equal (i, static_cast<unsigned int>(order));
42 
43  // monomials. since they are hierarchic we only need one case block.
44  switch (i)
45  {
46  case 0:
47  return 1.;
48 
49  case 1:
50  return xi;
51 
52  case 2:
53  return xi*xi;
54 
55  case 3:
56  return xi*xi*xi;
57 
58  case 4:
59  return xi*xi*xi*xi;
60 
61  default:
62  Real val = 1.;
63  for (unsigned int index = 0; index != i; ++index)
64  val *= xi;
65  return val;
66  }
67 }
68 
69 
70 
71 template <>
73  const Order order,
74  const unsigned int i,
75  const Point & p,
76  const bool add_p_level)
77 {
78  libmesh_assert(elem);
79 
80  return FE<1,MONOMIAL>::shape(elem->type(), static_cast<Order>(order + add_p_level * elem->p_level()), i, p);
81 }
82 
83 
84 template <>
86  const Elem * elem,
87  const unsigned int i,
88  const Point & p,
89  const bool add_p_level)
90 {
91  libmesh_assert(elem);
92  return FE<1,MONOMIAL>::shape(elem->type(), static_cast<Order>(fet.order + add_p_level * elem->p_level()), i, p);
93 }
94 
95 
96 template <>
98  const Order libmesh_dbg_var(order),
99  const unsigned int i,
100  const unsigned int libmesh_dbg_var(j),
101  const Point & p)
102 {
103  // only d()/dxi in 1D!
104 
105  libmesh_assert_equal_to (j, 0);
106 
107  const Real xi = p(0);
108 
109  libmesh_assert_less_equal (i, static_cast<unsigned int>(order));
110 
111  // monomials. since they are hierarchic we only need one case block.
112  switch (i)
113  {
114  case 0:
115  return 0.;
116 
117  case 1:
118  return 1.;
119 
120  case 2:
121  return 2.*xi;
122 
123  case 3:
124  return 3.*xi*xi;
125 
126  case 4:
127  return 4.*xi*xi*xi;
128 
129  default:
130  Real val = i;
131  for (unsigned int index = 1; index != i; ++index)
132  val *= xi;
133  return val;
134  }
135 }
136 
137 
138 
139 template <>
141  const Order order,
142  const unsigned int i,
143  const unsigned int j,
144  const Point & p,
145  const bool add_p_level)
146 {
147  libmesh_assert(elem);
148 
149  return FE<1,MONOMIAL>::shape_deriv(elem->type(),
150  static_cast<Order>(order + add_p_level * elem->p_level()), i, j, p);
151 }
152 
153 
154 template <>
156  const Elem * elem,
157  const unsigned int i,
158  const unsigned int j,
159  const Point & p,
160  const bool add_p_level)
161 {
162  libmesh_assert(elem);
163  return FE<1,MONOMIAL>::shape_deriv(elem->type(), static_cast<Order>(fet.order + add_p_level * elem->p_level()), i, j, p);
164 
165 
166 }
167 
168 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
169 
170 template <>
172  const Order libmesh_dbg_var(order),
173  const unsigned int i,
174  const unsigned int libmesh_dbg_var(j),
175  const Point & p)
176 {
177  // only d()/dxi in 1D!
178 
179  libmesh_assert_equal_to (j, 0);
180 
181  const Real xi = p(0);
182 
183  libmesh_assert_less_equal (i, static_cast<unsigned int>(order));
184 
185  switch (i)
186  {
187  case 0:
188  case 1:
189  return 0.;
190 
191  case 2:
192  return 2.;
193 
194  case 3:
195  return 6.*xi;
196 
197  case 4:
198  return 12.*xi*xi;
199 
200  default:
201  Real val = 2.;
202  for (unsigned int index = 2; index != i; ++index)
203  val *= (index+1) * xi;
204  return val;
205  }
206 }
207 
208 
209 
210 template <>
212  const Order order,
213  const unsigned int i,
214  const unsigned int j,
215  const Point & p,
216  const bool add_p_level)
217 {
218  libmesh_assert(elem);
219 
221  static_cast<Order>(order + add_p_level * elem->p_level()), i, j, p);
222 }
223 
224 template <>
226  const Elem * elem,
227  const unsigned int i,
228  const unsigned int j,
229  const Point & p,
230  const bool add_p_level)
231 {
232  libmesh_assert(elem);
233  return FE<1,MONOMIAL>::shape_second_deriv(elem->type(), static_cast<Order>(fet.order + add_p_level * elem->p_level()), i, j, p);
234 }
235 
236 #endif
237 
238 } // namespace libMesh
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:182
ElemType
Defines an enum for geometric element types.
Order
defines an enum for polynomial orders.
Definition: enum_order.h:40
static OutputShape shape(const ElemType t, const Order o, const unsigned int i, const Point &p)
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
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:2945
OrderWrapper order
The approximation order of the element.
Definition: fe_type.h:201
The libMesh namespace provides an interface to certain functionality in the library.
LIBMESH_DEFAULT_VECTORIZED_FE(template<>Real FE< 0, BERNSTEIN)
libmesh_assert(ctx)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual ElemType type() const =0
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39
static OutputShape shape_second_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)