v0.9.0
h1_hdiv_hcurl_l2.h
Go to the documentation of this file.
1 /** \file h1_hdiv_hcurl_l2.h
2 \brief Functions to approximate hierarchical spaces
3 
4 \FIXME: Name Shape Functions is used, in that context is more appropriate
5 to use base functions. Need to be changed.
6 
7 */
8 
9 /*
10  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
11  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
13  * License for more details.
14  *
15  * You should have received a copy of the GNU Lesser General Public
16  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
17  */
18 
19 #ifndef __H1_HDIV_HCURL_L2_H__
20 #define __H1_HDIV_HCURL_L2_H__
21 
22 #ifndef DEPRECATED
23  #define DEPRECATED
24 #endif
25 
26 #ifdef __cplusplus
27 extern "C" {
28 #endif
29 
30 // L2
31 
32 // Number of dofs for L2 space
33 
34 /**
35  * @brief Number of base functions on tetrahedron for L2 space
36  */
37 #define NBVOLUMETET_L2(P) ((P + 1) * (P + 2) * (P + 3) / 6)
38 
39 /**
40  * @brief Number of base functions on triangle for L2 space
41  */
42 #define NBFACETRI_L2(P) ((P + 1) * (P + 2) / 2)
43 
44 /**
45  * @brief Number of base functions on edge fro L2 space
46  *
47  */
48 #define NBEDGE_L2(P) (P + 1)
49 
50 // H1
51 
52 /**
53  * @brief Numer of base function on edge for H1 space
54  */
55 #define NBEDGE_H1(P) (((P) > 0) ? (P - 1) : 0)
56 
57 /**
58  * @brief Number of base function on triangle for H1 space
59  */
60 #define NBFACETRI_H1(P) (((P) > 1) ? ((P - 2) * (P - 1) / 2) : 0)
61 
62 /**
63  * @brief Number of base functions on quad for H1 space
64  */
65 #define NBFACEQUAD_H1(P) (((P) > 2) ? ((P - 3) * (P - 2) / 2) : 0)
66 
67 /**
68  * @brief Number of base functions on tetrahedron for H1 space
69  */
70 #define NBVOLUMETET_H1(P) (((P) > 2) ? ((P - 3) * (P - 2) * (P - 1) / 6) : 0)
71 
72 /**
73  * @brief Number of base functions on prism for H1 space
74  */
75 #define NBVOLUMEPRISM_H1(P) ((P > 4) ? ((P - 2) * (P - 3) * (P - 4) / 6) : 0)
76 
77 // H curl
78 
79 #define NBEDGE_AINSWORTH_HCURL(P) (((P) > 0) ? (P + 1) : 0)
80 #define NBFACETRI_AINSWORTH_EDGE_HCURL(P) (((P) > 1) ? P - 1 : 0)
81 #define NBFACETRI_AINSWORTH_FACE_HCURL(P) (((P) > 2) ? (P - 1) * (P - 2) : 0)
82 #define NBFACETRI_AINSWORTH_HCURL(P) ((P > 1) ? ((P)-1) * (P + 1) : 0)
83 #define NBVOLUMETET_AINSWORTH_FACE_HCURL(P) \
84  (((P) > 2) ? (2 * (P - 1) * (P - 2)) : 0)
85 #define NBVOLUMETET_AINSWORTH_TET_HCURL(P) \
86  (((P) > 3) ? ((P - 3) * (P - 2) * (P - 1) / 2) : 0)
87 #define NBVOLUMETET_AINSWORTH_HCURL(P) \
88  (((P) > 2) ? (P - 2) * (P - 1) * (P + 1) / 2 : 0)
89 
90 #define NBEDGE_DEMKOWICZ_HCURL(P) (((P) > 0) ? (P) : 0)
91 #define NBFACETRI_DEMKOWICZ_HCURL(P) (((P) > 1) ? (P) * ((P)-1) : 0)
92 #define NBVOLUMETET_DEMKOWICZ_HCURL(P) \
93  (((P) > 2) ? ((P) * ((P)-1) * ((P)-2) / 2) : 0)
94 
95 // H div
96 
97 #define NBEDGE_HDIV(P) (0)
98 #define NBFACETRI_AINSWORTH_EDGE_HDIV(P) (((P) > 0) ? (P) : 0)
99 #define NBFACETRI_AINSWORTH_FACE_HDIV(P) (((P) > 2) ? (P - 1) * (P - 2) / 2 : 0)
100 #define NBFACETRI_AINSWORTH_HDIV(P) (((P) > 0) ? (P + 1) * (P + 2) / 2 : 0)
101 #define NBVOLUMETET_AINSWORTH_EDGE_HDIV(P) (((P) > 1) ? (P - 1) : 0)
102 #define NBVOLUMETET_AINSWORTH_FACE_HDIV(P) (((P) > 2) ? (P - 1) * (P - 2) : 0)
103 #define NBVOLUMETET_AINSWORTH_VOLUME_HDIV(P) \
104  (((P) > 3) ? (P - 3) * (P - 2) * (P - 1) / 2 : 0)
105 #define NBVOLUMETET_AINSWORTH_HDIV(P) \
106  (((P) > 1) ? (P - 1) * (P + 1) * (P + 2) / 2 : 0)
107 #define NBFACETRI_DEMKOWICZ_HDIV(P) ((P > 0) ? (P) * (P + 1) / 2 : 0)
108 #define NBVOLUMETET_DEMKOWICZ_HDIV(P) \
109  (((P) > 1) ? (P) * (P - 1) * (P + 1) / 2 : 0)
110 
111 // Bubbles for H div space
112 
113 /**
114  * @brief Get base functions on triangle for L2 space
115  *
116  * @param p polynomial order
117  * @param N barycentric coordinates (shape functions) at integration points
118  * @param diffN direvatives of barycentric coordinates, i.e. direvatives of shape functions
119  * @param L2N values of L2 base at integration points
120  * @param diff_L2N dirvatives of base functions at integration points
121  * @param GDIM number of integration points
122  * @param base_polynomials polynomial base used to construct L2 base on element
123  * @return PetscErrorCode
124  */
125 PetscErrorCode L2_Ainsworth_ShapeFunctions_MBTRI(
126  int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM,
127  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s,
128  double *L, double *diffL,
129  const int dim));
130 /**
131  * @brief Get base functions on tetrahedron for L2 space
132  *
133  * @param p polynomial order
134  * @param N barycentric coordinates (shape functions) at integration points
135  * @param diffN direvatives of barycentric coordinates, i.e. direvatives of shape functions
136  * @param L2N values of L2 base at integration points
137  * @param diff_L2N dirvatives of base functions at integration points
138  * @param GDIM number of integration points
139  * @param base_polynomials polynomial base used to construct L2 base on element
140  * @return PetscErrorCode
141  */
142 PetscErrorCode L2_Ainsworth_ShapeFunctions_MBTET(
143  int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM,
144  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s,
145  double *L, double *diffL,
146  const int dim));
147 
148 /**
149  * \deprecated Use L2_Ainsworth_ShapeFunctions_MBTRI
150  */
151 DEPRECATED PetscErrorCode L2_ShapeFunctions_MBTRI(
152  int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM,
153  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s,
154  double *L, double *diffL,
155  const int dim));
156 /**
157  * \deprecated Use L2_Ainsworth_ShapeFunctions_MBTET
158  */
159 DEPRECATED PetscErrorCode L2_ShapeFunctions_MBTET(
160  int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM,
161  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s,
162  double *L, double *diffL,
163  const int dim));
164 
165 PetscErrorCode L2_VolumeShapeDiffMBTETinvJ(int base_p, int p,
166  double *volume_diffN, double *invJac,
167  double *volume_diffNinvJac,
168  int GDIM);
169 
170 /**
171  * \brief H1_EdgeShapeFunctions_MBTRI
172  *
173  * \param sense of edges, it is array of integers dim 3 (3-edges of triangle)
174  * \param p of edges
175  */
176 PetscErrorCode H1_EdgeShapeFunctions_MBTRI(
177  int *sense, int *p, double *N, double *diffN, double *edgeN[3],
178  double *diff_edgeN[3], int GDIM,
179  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s,
180  double *L, double *diffL,
181  const int dim));
182 PetscErrorCode H1_FaceShapeFunctions_MBTRI(
183  const int *face_nodes, int p, double *N, double *diffN, double *faceN,
184  double *diff_faceN, int GDIM,
185  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s,
186  double *L, double *diffL,
187  const int dim));
188 PetscErrorCode H1_EdgeShapeFunctions_MBTET(
189  int *sense, int *p, double *N, double *diffN, double *edgeN[],
190  double *diff_edgeN[], int GDIM,
191  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s,
192  double *L, double *diffL,
193  const int dim));
194 PetscErrorCode H1_FaceShapeFunctions_MBTET(
195  int *faces_nodes, int *p, double *N, double *diffN, double *faceN[],
196  double *diff_faceN[], int GDIM,
197  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s,
198  double *L, double *diffL,
199  const int dim));
200 PetscErrorCode H1_VolumeShapeFunctions_MBTET(
201  int p, double *N, double *diffN, double *volumeN, double *diff_volumeN,
202  int GDIM,
203  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s,
204  double *L, double *diffL,
205  const int dim));
206 PetscErrorCode H1_EdgeShapeDiffMBTETinvJ(int *base_p, int *p,
207  double *edge_diffN[], double *invJac,
208  double *edge_diffNinvJac[], int GDIM);
209 PetscErrorCode H1_FaceShapeDiffMBTETinvJ(int *base_p, int *p,
210  double *face_diffN[], double *invJac,
211  double *face_diffNinvJac[], int GDIM);
212 PetscErrorCode H1_VolumeShapeDiffMBTETinvJ(int base_p, int p,
213  double *volume_diffN, double *invJac,
214  double *volume_diffNinvJac,
215  int GDIM);
216 
217 PetscErrorCode H1_EdgeGradientOfDeformation_hierarchical(int p, double *diffN,
218  double *dofs,
219  double *F);
220 PetscErrorCode H1_FaceGradientOfDeformation_hierarchical(int p, double *diffN,
221  double *dofs,
222  double *F);
223 PetscErrorCode H1_VolumeGradientOfDeformation_hierarchical(int p, double *diffN,
224  double *dofs, double *F);
225 
226 /**
227  * \deprecated use H1_EdgeGradientOfDeformation_hierarchical
228  */
230  int p, double *diffN, double *dofs, double *F);
231 
232 /**
233  * \deprecated use H1_FaceGradientOfDeformation_hierarchical
234  */
236  int p, double *diffN, double *dofs, double *F);
237 
238 /**
239  * \deprecated use H1_VolumeGradientOfDeformation_hierarchical
240  */
242  int p, double *diffN, double *dofs, double *F);
243 
244 PetscErrorCode H1_QuadShapeFunctions_MBPRISM(
245  int *faces_nodes, int *p, double *N, double *diffN, double *faceN[],
246  double *diff_faceN[], int GDIM,
247  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s,
248  double *L, double *diffL,
249  const int dim));
250 PetscErrorCode H1_VolumeShapeFunctions_MBPRISM(
251  int p, double *N, double *diffN, double *volumeN, double *diff_volumeN,
252  int GDIM,
253  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s,
254  double *L, double *diffL,
255  const int dim));
256 PetscErrorCode H1_QuadShapeFunctions_MBQUAD(
257  int *faces_nodes, int p, double *N, double *diffN, double *faceN,
258  double *diff_faceN, int GDIM,
259  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s,
260  double *L, double *diffL,
261  const int dim));
262 PetscErrorCode H1_EdgeShapeFunctions_MBQUAD(
263  int *sense, int *p, double *N, double *diffN, double *edgeN[4],
264  double *diff_edgeN[4], int GDIM,
265  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s,
266  double *L, double *diffL, const int dim));
267 
268 // Hdiv and Hcurl are implemented and declared in other files
269 
270 
271 
272 #ifdef __cplusplus
273 }
274 #endif
275 
276 #endif // __H1_HDIV_HCURL_L2_H__
PetscErrorCode H1_FaceShapeFunctions_MBTET(int *faces_nodes, int *p, double *N, double *diffN, double *faceN[], double *diff_faceN[], int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: h1.c:317
#define DEPRECATED
PetscErrorCode H1_VolumeShapeFunctions_MBTET(int p, double *N, double *diffN, double *volumeN, double *diff_volumeN, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: h1.c:419
PetscErrorCode H1_EdgeShapeFunctions_MBQUAD(int *sense, int *p, double *N, double *diffN, double *edgeN[4], double *diff_edgeN[4], int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: h1.c:937
PetscErrorCode H1_EdgeShapeFunctions_MBTET(int *sense, int *p, double *N, double *diffN, double *edgeN[], double *diff_edgeN[], int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: h1.c:245
PetscErrorCode H1_EdgeShapeFunctions_MBTRI(int *sense, int *p, double *N, double *diffN, double *edgeN[3], double *diff_edgeN[3], int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
H1_EdgeShapeFunctions_MBTRI.
Definition: h1.c:27
PetscErrorCode H1_QuadShapeFunctions_MBQUAD(int *faces_nodes, int p, double *N, double *diffN, double *faceN, double *diff_faceN, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: h1.c:827
PetscErrorCode H1_VolumeGradientOfDeformation_hierarchical(int p, double *diffN, double *dofs, double *F)
Definition: h1.c:575
PetscErrorCode H1_VolumeShapeDiffMBTETinvJ(int base_p, int p, double *volume_diffN, double *invJac, double *volume_diffNinvJac, int GDIM)
Definition: h1.c:536
PetscErrorCode L2_VolumeShapeDiffMBTETinvJ(int base_p, int p, double *volume_diffN, double *invJac, double *volume_diffNinvJac, int GDIM)
Definition: l2.c:158
PetscErrorCode H1_FaceShapeFunctions_MBTRI(const int *face_nodes, int p, double *N, double *diffN, double *faceN, double *diff_faceN, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: h1.c:162
DEPRECATED PetscErrorCode L2_ShapeFunctions_MBTET(int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: l2.c:185
PetscErrorCode H1_EdgeGradientOfDeformation_hierarchical(int p, double *diffN, double *dofs, double *F)
Definition: h1.c:553
PetscErrorCode H1_FaceGradientOfDeformation_hierarchical(int p, double *diffN, double *dofs, double *F)
Definition: h1.c:564
DEPRECATED PetscErrorCode H1_VolumeGradientOfDeformation_hierachical(int p, double *diffN, double *dofs, double *F)
Definition: h1.c:1132
PetscErrorCode H1_FaceShapeDiffMBTETinvJ(int *base_p, int *p, double *face_diffN[], double *invJac, double *face_diffNinvJac[], int GDIM)
Definition: h1.c:518
DEPRECATED PetscErrorCode H1_FaceGradientOfDeformation_hierachical(int p, double *diffN, double *dofs, double *F)
Definition: h1.c:1127
PetscErrorCode H1_QuadShapeFunctions_MBPRISM(int *faces_nodes, int *p, double *N, double *diffN, double *faceN[], double *diff_faceN[], int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: h1.c:586
PetscErrorCode H1_VolumeShapeFunctions_MBPRISM(int p, double *N, double *diffN, double *volumeN, double *diff_volumeN, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: h1.c:707
DEPRECATED PetscErrorCode H1_EdgeGradientOfDeformation_hierachical(int p, double *diffN, double *dofs, double *F)
Definition: h1.c:1122
PetscErrorCode H1_EdgeShapeDiffMBTETinvJ(int *base_p, int *p, double *edge_diffN[], double *invJac, double *edge_diffNinvJac[], int GDIM)
Definition: h1.c:500
const int N
Definition: speed_test.cpp:3
DEPRECATED PetscErrorCode L2_ShapeFunctions_MBTRI(int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: l2.c:176
PetscErrorCode L2_Ainsworth_ShapeFunctions_MBTET(int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Get base functions on tetrahedron for L2 space.
Definition: l2.c:84
PetscErrorCode L2_Ainsworth_ShapeFunctions_MBTRI(int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Get base functions on triangle for L2 space.
Definition: l2.c:29