v0.13.0
fem_tools.h
Go to the documentation of this file.
1 /** \file fem_tools.h
2  * \brief Loose implementation of some useful functions
3  *
4  * FIXME: Implementation here is very unstructured, need cleaning and pruning
5  *
6  */
7 
8 /*
9  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
10  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
11  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
12  * License for more details.
13  *
14  * You should have received a copy of the GNU Lesser General Public
15  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
16  */
17 
18 #ifndef __FEM_H__
19 #define __FEM_H__
20 
21 #include <petscsys.h>
22 
23 #ifdef __cplusplus
24 extern "C" {
25 #endif
26 #include <cblas.h>
27 #include <lapack_wrap.h>
28 #ifdef __cplusplus
29 }
30 #endif
31 
32 #define LAMBDA(E, NU) (E * NU / ((1. + NU) * (1. - 2. * NU)))
33 #define MU(E, NU) (0.5 * E / (1. + NU))
34 #define DELTA(NU_P, NU_PZ, E_P, E_Z) \
35  (((1 + NU_P) * (1 - NU_P - 2 * NU_PZ * (NU_PZ * E_Z / E_P))) / \
36  (E_P * E_P * E_Z))
37 
38 #define N_MBTET0(x, y, z) (1. - x - y - z) ///< tetrahedral shape function
39 #define N_MBTET1(x, y, z) (x) ///< tetrahedral shape function
40 #define N_MBTET2(x, y, z) (y) ///< tetrahedral shape function
41 #define N_MBTET3(x, y, z) (z) ///< tetrahedral shape function
42 #define diffN_MBTET0x (-1.) ///< derivative of tetrahedral shape function
43 #define diffN_MBTET0y (-1.) ///< derivative of tetrahedral shape function
44 #define diffN_MBTET0z (-1.) ///< derivative of tetrahedral shape function
45 #define diffN_MBTET1x (1.) ///< derivative of tetrahedral shape function
46 #define diffN_MBTET1y (0.) ///< derivative of tetrahedral shape function
47 #define diffN_MBTET1z (0.) ///< derivative of tetrahedral shape function
48 #define diffN_MBTET2x (0.) ///< derivative of tetrahedral shape function
49 #define diffN_MBTET2y (1.) ///< derivative of tetrahedral shape function
50 #define diffN_MBTET2z (0.) ///< derivative of tetrahedral shape function
51 #define diffN_MBTET3x (0.) ///< derivative of tetrahedral shape function
52 #define diffN_MBTET3y (0.) ///< derivative of tetrahedral shape function
53 #define diffN_MBTET3z (1.) ///< derivative of tetrahedral shape function
54 
55 // MBTRI
56 #define N_MBTRI0(x, y) (1. - x - y) ///< triangle shape function
57 #define N_MBTRI1(x, y) (x) ///< triangle shape function
58 #define N_MBTRI2(x, y) (y) ///< triangle shape function
59 #define diffN_MBTRI0x (-1.) ///< derivative of triangle shape function
60 #define diffN_MBTRI0y (-1.) ///< derivative of triangle shape function
61 #define diffN_MBTRI1x (1.) ///< derivative of triangle shape function
62 #define diffN_MBTRI1y (0.) ///< derivative of triangle shape function
63 #define diffN_MBTRI2x (0.) ///< derivative of triangle shape function
64 #define diffN_MBTRI2y (1.) ///< derivative of triangle shape function
65 
66 // MBQUAD
67 #define N_MBQUAD0(x, y) ((1. - x) * (1. - y)) ///< quad shape function
68 #define N_MBQUAD1(x, y) ((x) * (1. - y)) ///< quad shape function
69 #define N_MBQUAD2(x, y) ((x) * (y)) ///< quad shape function
70 #define N_MBQUAD3(x, y) ((1. - x) * (y)) ///< quad shape function
71 #define diffN_MBQUAD0x(y) (-(1. - y))
72 #define diffN_MBQUAD0y(x) (-(1. - x))
73 #define diffN_MBQUAD1x(y) ((1. - y))
74 #define diffN_MBQUAD1y(x) (-x)
75 #define diffN_MBQUAD2x(y) (y)
76 #define diffN_MBQUAD2y(x) (x)
77 #define diffN_MBQUAD3x(y) (-y)
78 #define diffN_MBQUAD3y(x) ((1. - x))
79 
80 // MBHEX
81 #define N_MBHEX0(x, y, z) (N_MBQUAD0(x, y) * (1 - z))
82 #define N_MBHEX1(x, y, z) (N_MBQUAD1(x, y) * (1 - z))
83 #define N_MBHEX2(x, y, z) (N_MBQUAD2(x, y) * (1 - z))
84 #define N_MBHEX3(x, y, z) (N_MBQUAD3(x, y) * (1 - z))
85 #define N_MBHEX4(x, y, z) (N_MBQUAD0(x, y) * (z))
86 #define N_MBHEX5(x, y, z) (N_MBQUAD1(x, y) * (z))
87 #define N_MBHEX6(x, y, z) (N_MBQUAD2(x, y) * (z))
88 #define N_MBHEX7(x, y, z) (N_MBQUAD3(x, y) * (z))
89 #define diffN_MBHEX0x(y, z) (diffN_MBQUAD0x(y) * (1 - z))
90 #define diffN_MBHEX1x(y, z) (diffN_MBQUAD1x(y) * (1 - z))
91 #define diffN_MBHEX2x(y, z) (diffN_MBQUAD2x(y) * (1 - z))
92 #define diffN_MBHEX3x(y, z) (diffN_MBQUAD3x(y) * (1 - z))
93 #define diffN_MBHEX4x(y, z) (diffN_MBQUAD0x(y) * (z))
94 #define diffN_MBHEX5x(y, z) (diffN_MBQUAD1x(y) * (z))
95 #define diffN_MBHEX6x(y, z) (diffN_MBQUAD2x(y) * (z))
96 #define diffN_MBHEX7x(y, z) (diffN_MBQUAD3x(y) * (z))
97 #define diffN_MBHEX0y(x, z) (diffN_MBQUAD0y(x) * (1 - z))
98 #define diffN_MBHEX1y(x, z) (diffN_MBQUAD1y(x) * (1 - z))
99 #define diffN_MBHEX2y(x, z) (diffN_MBQUAD2y(x) * (1 - z))
100 #define diffN_MBHEX3y(x, z) (diffN_MBQUAD3y(x) * (1 - z))
101 #define diffN_MBHEX4y(x, z) (diffN_MBQUAD0y(x) * (z))
102 #define diffN_MBHEX5y(x, z) (diffN_MBQUAD1y(x) * (z))
103 #define diffN_MBHEX6y(x, z) (diffN_MBQUAD2y(x) * (z))
104 #define diffN_MBHEX7y(x, z) (diffN_MBQUAD3y(x) * (z))
105 #define diffN_MBHEX0z(x, y) (-N_MBQUAD0(x, y))
106 #define diffN_MBHEX1z(x, y) (-N_MBQUAD1(x, y))
107 #define diffN_MBHEX2z(x, y) (-N_MBQUAD2(x, y))
108 #define diffN_MBHEX3z(x, y) (-N_MBQUAD3(x, y))
109 #define diffN_MBHEX4z(x, y) (N_MBQUAD0(x, y))
110 #define diffN_MBHEX5z(x, y) (N_MBQUAD1(x, y))
111 #define diffN_MBHEX6z(x, y) (N_MBQUAD2(x, y))
112 #define diffN_MBHEX7z(x, y) (N_MBQUAD3(x, y))
113 
114 // MBEDGE
115 #define N_MBEDGE0(x) (1. - (x)) ///< edge shape function
116 #define N_MBEDGE1(x) (x) ///< edge shape function
117 #define diffN_MBEDGE0 (-1.) ///< derivative of edge shape function
118 #define diffN_MBEDGE1 (1.) ///< derivative of edge shape function
119 
120 // MBTRIQ
121 #define N_MBTRIQ0(x, y) ((1. - x - y) * (2 * (1. - x - y) - 1.))
122 #define N_MBTRIQ1(x, y) (x * (2. * x - 1.))
123 #define N_MBTRIQ2(x, y) (y * (2. * y - 1.))
124 #define N_MBTRIQ3(x, y) (4. * (1. - x - y) * x)
125 #define N_MBTRIQ4(x, y) (4. * x * y)
126 #define N_MBTRIQ5(x, y) (4. * (1. - x - y) * y)
127 #define diffN_MBTRIQ0x(x, y) (x + y - 3. * (1. - x - y))
128 #define diffN_MBTRIQ0y(x, y) (x + y - 3. * (1. - x - y))
129 #define diffN_MBTRIQ1x(x, y) (-1. + 4. * x)
130 #define diffN_MBTRIQ1y(x, y) (0.)
131 #define diffN_MBTRIQ2x(x, y) (0.)
132 #define diffN_MBTRIQ2y(x, y) (-1. + 4. * y)
133 #define diffN_MBTRIQ3x(x, y) (4. * ((1. - x - y) - x))
134 #define diffN_MBTRIQ3y(x, y) (-4. * x)
135 #define diffN_MBTRIQ4x(x, y) (4. * y)
136 #define diffN_MBTRIQ4y(x, y) (4. * x)
137 #define diffN_MBTRIQ5x(x, y) (-4. * y)
138 #define diffN_MBTRIQ5y(x, y) (4. * ((1. - x - y) - y))
139 
140 #ifdef __cplusplus
141 extern "C" {
142 #endif
143 
144 /// print matric M
145 void print_mat(double *M, int m, int n);
146 /// print upper part of the symmetric matrix
147 void print_mat_sym_upper(double *M, int m, int n);
148 /// priint complex matrix
150 
151 /// \brief calculate shape functions on triangle
152 /// \param N shape function array
153 /// \param X array of Gauss X coordinates
154 /// \param Y array of Gauss Y coordinates
155 /// \param G_DIM number of Gauss points
156 PetscErrorCode ShapeMBTRI(double *N, const double *X, const double *Y,
157  const int G_DIM);
158 /// calculate derivatives of shape functions
159 PetscErrorCode ShapeDiffMBTRI(double *diffN);
160 
161 /// calculate face normal
162 /// \param diffN derivatives of shape functions
163 /// \param coords is position of the nodes
164 /// \param normal vector
165 PetscErrorCode ShapeFaceNormalMBTRI(double *diffN, const double *coords,
166  double *normal);
167 PetscErrorCode ShapeFaceBaseMBTRI(double *diffN, const double *coords,
168  double *normal, double *s1, double *s2);
169 
170 /// calculate derivative of normal in respect to nodal positions
171 PetscErrorCode ShapeFaceDiffNormalMBTRI(double *diffN, const double *coords,
172  double *diff_normal);
173 /// calculate jacobioan
174 void ShapeJacMBTRI(double *diffN, const double *coords, double *Jac);
175 /// calculate derivatives of shape functions in space
176 void ShapeDiffMBTRIinvJ(double *diffN, double *invJac, double *diffNinvJac);
177 /// calculate shape functions
178 PetscErrorCode ShapeMBTET(double *N, const double *G_X, const double *G_Y,
179  const double *G_Z, int DIM);
180 /// calculate derivatives of shape functions
181 PetscErrorCode ShapeDiffMBTET(double *diffN);
182 /// determined of jacobian
183 double ShapeDetJacVolume(double *jac);
184 /// calculate jacobian
185 PetscErrorCode ShapeJacMBTET(double *diffN, const double *coords, double *jac);
186 // calculate inverse of jacobian
187 PetscErrorCode ShapeInvJacVolume(double *jac);
188 /// calculate TET volume
189 double ShapeVolumeMBTET(double *diffN, const double *coords);
190 /// calculate shape functions derivatives in space
191 PetscErrorCode ShapeDiffMBTETinvJ(double *diffN, double *invJac,
192  double *diffNinvJac);
193 
194 /// calculate spin matrix from vector
195 // \param spinOmega is a spin matrix
196 // \param vecOmega is a spin vector
197 PetscErrorCode Spin(double *spinOmega, double *vecOmega);
198 
199 /// Compose complex matrix (3x3) from two real matrices
200 PetscErrorCode make_complex_matrix(double *reA, double *imA,
202 /// Complex normal
203 PetscErrorCode
204 Normal_hierarchical(int order_approx, int *order_edge_approx, int order,
205  int *order_edge, double *diffN, double *diffN_face,
206  double *diffN_edge[], double *dofs, double *dofs_edge[],
207  double *dofs_face, double *idofs, double *idofs_edge[],
208  double *idofs_face, __CLPK_doublecomplex *xnormal,
209  __CLPK_doublecomplex *s1, __CLPK_doublecomplex *s2, int gg);
210 PetscErrorCode Base_scale(__CLPK_doublecomplex *xnormal,
212 
213 /**
214  * \brief calculate local coordinates for given global coordinates
215  *
216  * new version for multiple points need to be implemented
217  */
218 PetscErrorCode ShapeMBTET_inverse(double *N, double *diffN,
219  const double *elem_coords,
220  const double *glob_coords,
221  double *loc_coords);
222 
223 /**
224  * \brief calculate local coordinates of triangle element for given global
225 coordinates in 2D (Assume e.g. z=0) \f[ \left[\begin{array}{cc}
226 \frac{\partial N_{1}}{\partial\xi}x_{N_{1}}+\frac{\partial
227 N_{2}}{\partial\xi}x_{N_{2}}+\frac{\partial N_{3}}{\partial\xi}x_{N_{3}} &
228 \frac{\partial N_{1}}{\partial\eta}x_{N_{1}}+\frac{\partial
229 N_{2}}{\partial\eta}x_{N_{2}}+\frac{\partial N_{3}}{\partial\eta}x_{N_{3}}\\
230 \frac{\partial N_{1}}{\partial\xi}y_{N_{1}}+\frac{\partial
231 N_{2}}{\partial\xi}y_{N_{2}}+\frac{\partial N_{3}}{\partial\xi}y_{N_{3}} &
232 \frac{\partial N_{1}}{\partial\eta}y_{N_{1}}+\frac{\partial
233 N_{2}}{\partial\eta}y_{N_{2}}+\frac{\partial N_{3}}{\partial\eta}y_{N_{3}}
234 \end{array}\right]\left\{ \begin{array}{c}
235 \xi\\
236 \eta
237 \end{array}\right\} =\left\{ \begin{array}{c}
238 x_{gp}-\left(N_{1}x_{N_{1}}+N_{2}x_{N_{2}}+N_{3}x_{N_{3}}\right)\\
239 y_{gp}-\left(N_{1}y_{N_{1}}+N_{2}y_{N_{2}}+N_{3}y_{N_{3}}\right)
240 \end{array}\right\}
241  \f]
242 
243  /// \param N shape function array
244  /// \param diffN array of shape function derivative w.r.t to local coordinates
245  /// \param elem_coords global coordinates of element
246  /// \param glob_coords global coordinates of required point
247  /// \param loc_coords local coordinates of required point
248  */
249 PetscErrorCode ShapeMBTRI_inverse(double *N, double *diffN,
250  const double *elem_coords,
251  const double *glob_coords,
252  double *loc_coords);
253 
254 /// calculate gradient of deformation
255 PetscErrorCode GradientOfDeformation(double *diffN, double *dofs, double *F);
256 
257 // 2 Node edge
258 PetscErrorCode ShapeMBEDGE(double *N, const double *G_X, int DIM);
259 PetscErrorCode ShapeDiffMBEDGE(double *diffN);
260 
261 // 10 Node Tet
262 PetscErrorCode ShapeMBTRIQ(double *N, const double *X, const double *Y,
263  const int G_DIM);
264 PetscErrorCode ShapeDiffMBTRIQ(double *diffN, const double *X, const double *Y,
265  const int G_DIM);
266 PetscErrorCode ShapeMBTETQ(double *N, const double x, const double y,
267  const double z);
268 PetscErrorCode ShapeDiffMBTETQ(double *diffN, const double x, const double y,
269  const double z);
270 PetscErrorCode ShapeMBTETQ_GAUSS(double *N, const double *X, const double *Y,
271  const double *Z, const int G_DIM);
272 PetscErrorCode ShapeDiffMBTETQ_GAUSS(double *diffN, const double *X,
273  const double *Y, const double *Z,
274  const int G_DIM);
275 PetscErrorCode ShapeJacMBTETQ(const double *diffN, const double *coords,
276  double *Jac);
277 PetscErrorCode
278 ShapeMBTETQ_detJac_at_Gauss_Points(double *detJac_at_Gauss_Points,
279  const double *diffN, const double *coords,
280  int G_DIM);
281 double ShapeVolumeMBTETQ(const double *diffN, const double *coords, int G_DIM,
282  double *G_TET_W);
283 PetscErrorCode ShapeMBTETQ_inverse(double *N, double *diffN,
284  const double *elem_coords,
285  const double *glob_coords,
286  double *loc_coords, const double eps);
287 
288 // complex part
290  __CLPK_doublecomplex *diffNinvJac,
291  CBLAS_TRANSPOSE Trans);
292 PetscErrorCode ShapeFaceNormalMBTRI_complex(double *diffN,
293  __CLPK_doublecomplex *xcoords,
294  __CLPK_doublecomplex *xnormal);
295 PetscErrorCode MakeComplexTensor(double *reA, double *imA,
297 PetscErrorCode InvertComplexGradient(__CLPK_doublecomplex *xF);
299  __CLPK_doublecomplex *det_xF);
300 
301 // integration
302 /// Compute weights and integration points for edge using Grundmann_Moeller rule
303 PetscErrorCode Grundmann_Moeller_integration_points_1D_EDGE(int rule,
304  double *G_TRI_X,
305  double *G_TRI_W);
306 /// Compute weights and integration points for 2D Triangle using
307 /// Grundmann_Moeller rule
308 PetscErrorCode Grundmann_Moeller_integration_points_2D_TRI(int rule,
309  double *G_TRI_X,
310  double *G_TRI_Y,
311  double *G_TRI_W);
312 /// Compute weights and integration points for 3D Tet using Grundmann_Moeller
313 /// rule
314 PetscErrorCode Grundmann_Moeller_integration_points_3D_TET(int rule,
315  double *G_TET_X,
316  double *G_TET_Y,
317  double *G_TET_Z,
318  double *G_TET_W);
319 
320 // http://people.sc.fsu.edu/~jburkardt/datasets/quadrature_rules_tet/quadrature_rules_tet.html
321 // TRI
322 static const double G_TRI_X1[] = {3.3333333333333331e-01};
323 static const double G_TRI_Y1[] = {3.3333333333333331e-01};
324 static const double G_TRI_W1[] = {1.};
325 static const double G_TRI_X3[] = {0.5, 0., 0.5};
326 static const double G_TRI_Y3[] = {0., 0.5, 0.5};
327 static const double G_TRI_W3[] = {
328  3.3333333333333331e-01, 3.3333333333333331e-01, 3.3333333333333331e-01};
329 static const double G_TRI_X4[] = {
330  7.503111022260811058e-02, 1.785587282636164064e-01,
331  2.800199154990741235e-01, 6.663902460147014262e-01};
332 static const double G_TRI_Y4[] = {
333  2.800199154990741235e-01, 6.663902460147014262e-01,
334  7.503111022260811058e-02, 1.785587282636164064e-01};
335 static const double G_TRI_W4[] = {
336  1.8195861825602258066e-01, 3.1804138174397683647e-01,
337  1.8195861825602258066e-01, 3.1804138174397683647e-01};
338 static const double G_TRI_X7[] = {
339  0.333333333333333, 0.736712498968435, 0.736712498968435, 0.237932366472434,
340  0.237932366472434, 0.025355134551932, 0.025355134551932};
341 static const double G_TRI_Y7[] = {
342  0.333333333333333, 0.237932366472434, 0.025355134551932, 0.736712498968435,
343  0.025355134551932, 0.736712498968435, 0.237932366472434};
344 static const double G_TRI_W7[] = {
345  0.375000000000000, 0.104166666666667, 0.104166666666667, 0.104166666666667,
346  0.104166666666667, 0.104166666666667, 0.104166666666667};
347 static const double G_TRI_X13[] = {
348  0.333333333333333, 0.479308067841923, 0.260345966079038, 0.260345966079038,
349  0.869739794195568, 0.065130102902216, 0.065130102902216, 0.638444188569809,
350  0.638444188569809, 0.312865496004875, 0.312865496004875, 0.048690315425316,
351  0.048690315425316};
352 static const double G_TRI_Y13[] = {
353  0.333333333333333, 0.260345966079038, 0.479308067841923, 0.260345966079038,
354  0.065130102902216, 0.869739794195568, 0.065130102902216, 0.312865496004875,
355  0.048690315425316, 0.638444188569809, 0.048690315425316, 0.638444188569809,
356  0.312865496004875};
357 static const double G_TRI_W13[] = {
358  -0.149570044467670, 0.175615257433204, 0.175615257433204, 0.175615257433204,
359  0.053347235608839, 0.053347235608839, 0.053347235608839, 0.077113760890257,
360  0.077113760890257, 0.077113760890257, 0.077113760890257, 0.077113760890257,
361  0.077113760890257};
362 
363 static const double G_TRI_X19[] = {
364  0.333333333333333, 0.797426985353087, 0.101286507323456, 0.101286507323456,
365  0.059715871789770, 0.470142064105115, 0.470142064105115, 0.535795346449899,
366  0.232102326775050, 0.232102326775050, 0.941038278231121, 0.029480860884440,
367  0.029480860884440, 0.738416812340510, 0.738416812340510, 0.232102326775050,
368  0.232102326775050, 0.029480860884440, 0.029480860884440};
369 static const double G_TRI_Y19[] = {
370  0.333333333333333, 0.101286507323456, 0.797426985353087, 0.101286507323456,
371  0.470142064105115, 0.059715871789770, 0.470142064105115, 0.232102326775050,
372  0.535795346449899, 0.232102326775050, 0.029480860884440, 0.941038278231121,
373  0.029480860884440, 0.232102326775050, 0.029480860884440, 0.738416812340510,
374  0.029480860884440, 0.738416812340510, 0.232102326775050};
375 static const double G_TRI_W19[] = {
376  9.71357962827961025E-002, 3.13347002271398278E-002,
377  3.13347002271398278E-002, 3.13347002271398278E-002,
378  7.78275410047754301E-002, 7.78275410047754301E-002,
379  7.78275410047754301E-002, 7.96477389272090969E-002,
380  7.96477389272090969E-002, 7.96477389272090969E-002,
381  2.55776756586981006E-002, 2.55776756586981006E-002,
382  2.55776756586981006E-002, 4.32835393772893970E-002,
383  4.32835393772893970E-002, 4.32835393772893970E-002,
384  4.32835393772893970E-002, 4.32835393772893970E-002,
385  4.32835393772893970E-002};
386 
387 static const double G_TRI_X28[] = {
388  0.333333333333333, 0.948021718143423, 0.025989140928288,
389  0.025989140928288, 0.811424994704155, 0.094287502647923,
390  0.094287502647923, 0.010726449965571, 0.494636775017215,
391  0.494636775017215, 0.585313234770972, 0.207343382614514,
392  0.207343382614514, 0.122184388599019, 0.438907805700491,
393  0.438907805700491, 0.677937654882590, 0.677937654882590,
394  0.044841677589131, 0.044841677589131, 0.277220667528279,
395  0.277220667528279, 0.858870281282636, 0.858870281282636,
396  0.0000000000000000, 0.0000000000000000, 0.141129718717364,
397  0.141129718717364};
398 static const double G_TRI_Y28[] = {
399  0.333333333333333, 0.025989140928288, 0.948021718143423, 0.025989140928288,
400  0.094287502647923, 0.811424994704155, 0.094287502647923, 0.494636775017215,
401  0.010726449965571, 0.494636775017215, 0.207343382614514, 0.585313234770972,
402  0.207343382614514, 0.438907805700491, 0.122184388599019, 0.438907805700491,
403  0.044841677589131, 0.277220667528279, 0.677937654882590, 0.277220667528279,
404  0.677937654882590, 0.044841677589131, 0.000000000000000, 0.141129718717364,
405  0.858870281282636, 0.141129718717364, 0.858870281282636, 0.000000000000000};
406 static const double G_TRI_W28[] = {
407  0.08797730116222190, 0.008744311553736190, 0.008744311553736190,
408  0.008744311553736190, 0.03808157199393533, 0.03808157199393533,
409  0.03808157199393533, 0.01885544805613125, 0.01885544805613125,
410  0.01885544805613125, 0.07215969754474100, 0.07215969754474100,
411  0.07215969754474100, 0.06932913870553720, 0.06932913870553720,
412  0.06932913870553720, 0.04105631542928860, 0.04105631542928860,
413  0.04105631542928860, 0.04105631542928860, 0.04105631542928860,
414  0.04105631542928860, 0.007362383783300573, 0.007362383783300573,
415  0.007362383783300573, 0.007362383783300573, 0.007362383783300573,
416  0.007362383783300573};
417 
418 static const double G_TRI_X37[] = {
419  0.333333333333333, 0.950275662924106, 0.024862168537947, 0.024862168537947,
420  0.171614914923835, 0.414192542538082, 0.414192542538082, 0.539412243677190,
421  0.230293878161405, 0.230293878161405, 0.772160036676533, 0.113919981661734,
422  0.113919981661734, 0.009085399949835, 0.495457300025082, 0.495457300025082,
423  0.062277290305887, 0.468861354847056, 0.468861354847056, 0.022076289653624,
424  0.022076289653624, 0.851306504174348, 0.851306504174348, 0.126617206172027,
425  0.126617206172027, 0.018620522802521, 0.018620522802521, 0.689441970728591,
426  0.689441970728591, 0.291937506468888, 0.291937506468888, 0.096506481292159,
427  0.096506481292159, 0.635867859433873, 0.635867859433873, 0.267625659273968,
428  0.267625659273968};
429 static const double G_TRI_Y37[] = {
430  0.333333333333333, 0.024862168537947, 0.950275662924106, 0.024862168537947,
431  0.414192542538082, 0.171614914923835, 0.414192542538082, 0.230293878161405,
432  0.539412243677190, 0.230293878161405, 0.113919981661734, 0.772160036676533,
433  0.113919981661734, 0.495457300025082, 0.009085399949835, 0.495457300025082,
434  0.468861354847056, 0.062277290305887, 0.468861354847056, 0.851306504174348,
435  0.126617206172027, 0.022076289653624, 0.126617206172027, 0.022076289653624,
436  0.851306504174348, 0.689441970728591, 0.291937506468888, 0.018620522802521,
437  0.291937506468888, 0.018620522802521, 0.689441970728591, 0.635867859433873,
438  0.267625659273968, 0.096506481292159, 0.267625659273968, 0.096506481292159,
439  0.635867859433873};
440 static const double G_TRI_W37[] = {
441  0.051739766065744, 0.008007799555565, 0.008007799555565, 0.008007799555565,
442  0.046868898981822, 0.046868898981822, 0.046868898981822, 0.046590940183976,
443  0.046590940183976, 0.046590940183976, 0.031016943313796, 0.031016943313796,
444  0.031016943313796, 0.010791612736631, 0.010791612736631, 0.010791612736631,
445  0.032195534242432, 0.032195534242432, 0.032195534242432, 0.015445834210702,
446  0.015445834210702, 0.015445834210702, 0.015445834210702, 0.015445834210702,
447  0.015445834210702, 0.017822989923179, 0.017822989923179, 0.017822989923179,
448  0.017822989923179, 0.017822989923179, 0.017822989923179, 0.037038683681385,
449  0.037038683681385, 0.037038683681385, 0.037038683681385, 0.037038683681385,
450  0.037038683681385};
451 static const double G_TRI_X286[] = {0.04347826086956522,
452  0.1304347826086956,
453  0.2173913043478261,
454  0.3043478260869565,
455  0.391304347826087,
456  0.4782608695652174,
457  0.5652173913043478,
458  0.6521739130434783,
459  0.7391304347826086,
460  0.8260869565217391,
461  0.9130434782608695,
462  0.04347826086956522,
463  0.1304347826086956,
464  0.2173913043478261,
465  0.3043478260869565,
466  0.391304347826087,
467  0.4782608695652174,
468  0.5652173913043478,
469  0.6521739130434783,
470  0.7391304347826086,
471  0.8260869565217391,
472  0.04347826086956522,
473  0.1304347826086956,
474  0.2173913043478261,
475  0.3043478260869565,
476  0.391304347826087,
477  0.4782608695652174,
478  0.5652173913043478,
479  0.6521739130434783,
480  0.7391304347826086,
481  0.04347826086956522,
482  0.1304347826086956,
483  0.2173913043478261,
484  0.3043478260869565,
485  0.391304347826087,
486  0.4782608695652174,
487  0.5652173913043478,
488  0.6521739130434783,
489  0.04347826086956522,
490  0.1304347826086956,
491  0.2173913043478261,
492  0.3043478260869565,
493  0.391304347826087,
494  0.4782608695652174,
495  0.5652173913043478,
496  0.04347826086956522,
497  0.1304347826086956,
498  0.2173913043478261,
499  0.3043478260869565,
500  0.391304347826087,
501  0.4782608695652174,
502  0.04347826086956522,
503  0.1304347826086956,
504  0.2173913043478261,
505  0.3043478260869565,
506  0.391304347826087,
507  0.04347826086956522,
508  0.1304347826086956,
509  0.2173913043478261,
510  0.3043478260869565,
511  0.04347826086956522,
512  0.1304347826086956,
513  0.2173913043478261,
514  0.04347826086956522,
515  0.1304347826086956,
516  0.04347826086956522,
517  0.04761904761904762,
518  0.1428571428571428,
519  0.2380952380952381,
520  0.3333333333333333,
521  0.4285714285714285,
522  0.5238095238095238,
523  0.6190476190476191,
524  0.7142857142857143,
525  0.8095238095238095,
526  0.9047619047619048,
527  0.04761904761904762,
528  0.1428571428571428,
529  0.2380952380952381,
530  0.3333333333333333,
531  0.4285714285714285,
532  0.5238095238095238,
533  0.6190476190476191,
534  0.7142857142857143,
535  0.8095238095238095,
536  0.04761904761904762,
537  0.1428571428571428,
538  0.2380952380952381,
539  0.3333333333333333,
540  0.4285714285714285,
541  0.5238095238095238,
542  0.6190476190476191,
543  0.7142857142857143,
544  0.04761904761904762,
545  0.1428571428571428,
546  0.2380952380952381,
547  0.3333333333333333,
548  0.4285714285714285,
549  0.5238095238095238,
550  0.6190476190476191,
551  0.04761904761904762,
552  0.1428571428571428,
553  0.2380952380952381,
554  0.3333333333333333,
555  0.4285714285714285,
556  0.5238095238095238,
557  0.04761904761904762,
558  0.1428571428571428,
559  0.2380952380952381,
560  0.3333333333333333,
561  0.4285714285714285,
562  0.04761904761904762,
563  0.1428571428571428,
564  0.2380952380952381,
565  0.3333333333333333,
566  0.04761904761904762,
567  0.1428571428571428,
568  0.2380952380952381,
569  0.04761904761904762,
570  0.1428571428571428,
571  0.04761904761904762,
572  0.05263157894736842,
573  0.1578947368421053,
574  0.2631578947368421,
575  0.3684210526315789,
576  0.4736842105263158,
577  0.5789473684210527,
578  0.6842105263157895,
579  0.7894736842105263,
580  0.8947368421052632,
581  0.05263157894736842,
582  0.1578947368421053,
583  0.2631578947368421,
584  0.3684210526315789,
585  0.4736842105263158,
586  0.5789473684210527,
587  0.6842105263157895,
588  0.7894736842105263,
589  0.05263157894736842,
590  0.1578947368421053,
591  0.2631578947368421,
592  0.3684210526315789,
593  0.4736842105263158,
594  0.5789473684210527,
595  0.6842105263157895,
596  0.05263157894736842,
597  0.1578947368421053,
598  0.2631578947368421,
599  0.3684210526315789,
600  0.4736842105263158,
601  0.5789473684210527,
602  0.05263157894736842,
603  0.1578947368421053,
604  0.2631578947368421,
605  0.3684210526315789,
606  0.4736842105263158,
607  0.05263157894736842,
608  0.1578947368421053,
609  0.2631578947368421,
610  0.3684210526315789,
611  0.05263157894736842,
612  0.1578947368421053,
613  0.2631578947368421,
614  0.05263157894736842,
615  0.1578947368421053,
616  0.05263157894736842,
617  0.05882352941176471,
618  0.1764705882352941,
619  0.2941176470588235,
620  0.4117647058823529,
621  0.5294117647058824,
622  0.6470588235294118,
623  0.7647058823529411,
624  0.8823529411764706,
625  0.05882352941176471,
626  0.1764705882352941,
627  0.2941176470588235,
628  0.4117647058823529,
629  0.5294117647058824,
630  0.6470588235294118,
631  0.7647058823529411,
632  0.05882352941176471,
633  0.1764705882352941,
634  0.2941176470588235,
635  0.4117647058823529,
636  0.5294117647058824,
637  0.6470588235294118,
638  0.05882352941176471,
639  0.1764705882352941,
640  0.2941176470588235,
641  0.4117647058823529,
642  0.5294117647058824,
643  0.05882352941176471,
644  0.1764705882352941,
645  0.2941176470588235,
646  0.4117647058823529,
647  0.05882352941176471,
648  0.1764705882352941,
649  0.2941176470588235,
650  0.05882352941176471,
651  0.1764705882352941,
652  0.05882352941176471,
653  0.06666666666666667,
654  0.2,
655  0.3333333333333333,
656  0.4666666666666667,
657  0.6,
658  0.7333333333333333,
659  0.8666666666666667,
660  0.06666666666666667,
661  0.2,
662  0.3333333333333333,
663  0.4666666666666667,
664  0.6,
665  0.7333333333333333,
666  0.06666666666666667,
667  0.2,
668  0.3333333333333333,
669  0.4666666666666667,
670  0.6,
671  0.06666666666666667,
672  0.2,
673  0.3333333333333333,
674  0.4666666666666667,
675  0.06666666666666667,
676  0.2,
677  0.3333333333333333,
678  0.06666666666666667,
679  0.2,
680  0.06666666666666667,
681  0.07692307692307693,
682  0.2307692307692308,
683  0.3846153846153846,
684  0.5384615384615384,
685  0.6923076923076923,
686  0.8461538461538461,
687  0.07692307692307693,
688  0.2307692307692308,
689  0.3846153846153846,
690  0.5384615384615384,
691  0.6923076923076923,
692  0.07692307692307693,
693  0.2307692307692308,
694  0.3846153846153846,
695  0.5384615384615384,
696  0.07692307692307693,
697  0.2307692307692308,
698  0.3846153846153846,
699  0.07692307692307693,
700  0.2307692307692308,
701  0.07692307692307693,
702  0.09090909090909091,
703  0.2727272727272727,
704  0.4545454545454545,
705  0.6363636363636364,
706  0.8181818181818182,
707  0.09090909090909091,
708  0.2727272727272727,
709  0.4545454545454545,
710  0.6363636363636364,
711  0.09090909090909091,
712  0.2727272727272727,
713  0.4545454545454545,
714  0.09090909090909091,
715  0.2727272727272727,
716  0.09090909090909091,
717  0.1111111111111111,
718  0.3333333333333333,
719  0.5555555555555556,
720  0.7777777777777778,
721  0.1111111111111111,
722  0.3333333333333333,
723  0.5555555555555556,
724  0.1111111111111111,
725  0.3333333333333333,
726  0.1111111111111111,
727  0.1428571428571428,
728  0.4285714285714285,
729  0.7142857142857143,
730  0.1428571428571428,
731  0.4285714285714285,
732  0.1428571428571428,
733  0.2,
734  0.6,
735  0.2,
736  0.3333333333333333};
737 static const double G_TRI_Y286[] = {0.04347826086956522,
738  0.04347826086956522,
739  0.04347826086956522,
740  0.04347826086956522,
741  0.04347826086956522,
742  0.04347826086956522,
743  0.04347826086956522,
744  0.04347826086956522,
745  0.04347826086956522,
746  0.04347826086956522,
747  0.04347826086956522,
748  0.1304347826086956,
749  0.1304347826086956,
750  0.1304347826086956,
751  0.1304347826086956,
752  0.1304347826086956,
753  0.1304347826086956,
754  0.1304347826086956,
755  0.1304347826086956,
756  0.1304347826086956,
757  0.1304347826086956,
758  0.2173913043478261,
759  0.2173913043478261,
760  0.2173913043478261,
761  0.2173913043478261,
762  0.2173913043478261,
763  0.2173913043478261,
764  0.2173913043478261,
765  0.2173913043478261,
766  0.2173913043478261,
767  0.3043478260869565,
768  0.3043478260869565,
769  0.3043478260869565,
770  0.3043478260869565,
771  0.3043478260869565,
772  0.3043478260869565,
773  0.3043478260869565,
774  0.3043478260869565,
775  0.391304347826087,
776  0.391304347826087,
777  0.391304347826087,
778  0.391304347826087,
779  0.391304347826087,
780  0.391304347826087,
781  0.391304347826087,
782  0.4782608695652174,
783  0.4782608695652174,
784  0.4782608695652174,
785  0.4782608695652174,
786  0.4782608695652174,
787  0.4782608695652174,
788  0.5652173913043478,
789  0.5652173913043478,
790  0.5652173913043478,
791  0.5652173913043478,
792  0.5652173913043478,
793  0.6521739130434783,
794  0.6521739130434783,
795  0.6521739130434783,
796  0.6521739130434783,
797  0.7391304347826086,
798  0.7391304347826086,
799  0.7391304347826086,
800  0.8260869565217391,
801  0.8260869565217391,
802  0.9130434782608695,
803  0.04761904761904762,
804  0.04761904761904762,
805  0.04761904761904762,
806  0.04761904761904762,
807  0.04761904761904762,
808  0.04761904761904762,
809  0.04761904761904762,
810  0.04761904761904762,
811  0.04761904761904762,
812  0.04761904761904762,
813  0.1428571428571428,
814  0.1428571428571428,
815  0.1428571428571428,
816  0.1428571428571428,
817  0.1428571428571428,
818  0.1428571428571428,
819  0.1428571428571428,
820  0.1428571428571428,
821  0.1428571428571428,
822  0.2380952380952381,
823  0.2380952380952381,
824  0.2380952380952381,
825  0.2380952380952381,
826  0.2380952380952381,
827  0.2380952380952381,
828  0.2380952380952381,
829  0.2380952380952381,
830  0.3333333333333333,
831  0.3333333333333333,
832  0.3333333333333333,
833  0.3333333333333333,
834  0.3333333333333333,
835  0.3333333333333333,
836  0.3333333333333333,
837  0.4285714285714285,
838  0.4285714285714285,
839  0.4285714285714285,
840  0.4285714285714285,
841  0.4285714285714285,
842  0.4285714285714285,
843  0.5238095238095238,
844  0.5238095238095238,
845  0.5238095238095238,
846  0.5238095238095238,
847  0.5238095238095238,
848  0.6190476190476191,
849  0.6190476190476191,
850  0.6190476190476191,
851  0.6190476190476191,
852  0.7142857142857143,
853  0.7142857142857143,
854  0.7142857142857143,
855  0.8095238095238095,
856  0.8095238095238095,
857  0.9047619047619048,
858  0.05263157894736842,
859  0.05263157894736842,
860  0.05263157894736842,
861  0.05263157894736842,
862  0.05263157894736842,
863  0.05263157894736842,
864  0.05263157894736842,
865  0.05263157894736842,
866  0.05263157894736842,
867  0.1578947368421053,
868  0.1578947368421053,
869  0.1578947368421053,
870  0.1578947368421053,
871  0.1578947368421053,
872  0.1578947368421053,
873  0.1578947368421053,
874  0.1578947368421053,
875  0.2631578947368421,
876  0.2631578947368421,
877  0.2631578947368421,
878  0.2631578947368421,
879  0.2631578947368421,
880  0.2631578947368421,
881  0.2631578947368421,
882  0.3684210526315789,
883  0.3684210526315789,
884  0.3684210526315789,
885  0.3684210526315789,
886  0.3684210526315789,
887  0.3684210526315789,
888  0.4736842105263158,
889  0.4736842105263158,
890  0.4736842105263158,
891  0.4736842105263158,
892  0.4736842105263158,
893  0.5789473684210527,
894  0.5789473684210527,
895  0.5789473684210527,
896  0.5789473684210527,
897  0.6842105263157895,
898  0.6842105263157895,
899  0.6842105263157895,
900  0.7894736842105263,
901  0.7894736842105263,
902  0.8947368421052632,
903  0.05882352941176471,
904  0.05882352941176471,
905  0.05882352941176471,
906  0.05882352941176471,
907  0.05882352941176471,
908  0.05882352941176471,
909  0.05882352941176471,
910  0.05882352941176471,
911  0.1764705882352941,
912  0.1764705882352941,
913  0.1764705882352941,
914  0.1764705882352941,
915  0.1764705882352941,
916  0.1764705882352941,
917  0.1764705882352941,
918  0.2941176470588235,
919  0.2941176470588235,
920  0.2941176470588235,
921  0.2941176470588235,
922  0.2941176470588235,
923  0.2941176470588235,
924  0.4117647058823529,
925  0.4117647058823529,
926  0.4117647058823529,
927  0.4117647058823529,
928  0.4117647058823529,
929  0.5294117647058824,
930  0.5294117647058824,
931  0.5294117647058824,
932  0.5294117647058824,
933  0.6470588235294118,
934  0.6470588235294118,
935  0.6470588235294118,
936  0.7647058823529411,
937  0.7647058823529411,
938  0.8823529411764706,
939  0.06666666666666667,
940  0.06666666666666667,
941  0.06666666666666667,
942  0.06666666666666667,
943  0.06666666666666667,
944  0.06666666666666667,
945  0.06666666666666667,
946  0.2,
947  0.2,
948  0.2,
949  0.2,
950  0.2,
951  0.2,
952  0.3333333333333333,
953  0.3333333333333333,
954  0.3333333333333333,
955  0.3333333333333333,
956  0.3333333333333333,
957  0.4666666666666667,
958  0.4666666666666667,
959  0.4666666666666667,
960  0.4666666666666667,
961  0.6,
962  0.6,
963  0.6,
964  0.7333333333333333,
965  0.7333333333333333,
966  0.8666666666666667,
967  0.07692307692307693,
968  0.07692307692307693,
969  0.07692307692307693,
970  0.07692307692307693,
971  0.07692307692307693,
972  0.07692307692307693,
973  0.2307692307692308,
974  0.2307692307692308,
975  0.2307692307692308,
976  0.2307692307692308,
977  0.2307692307692308,
978  0.3846153846153846,
979  0.3846153846153846,
980  0.3846153846153846,
981  0.3846153846153846,
982  0.5384615384615384,
983  0.5384615384615384,
984  0.5384615384615384,
985  0.6923076923076923,
986  0.6923076923076923,
987  0.8461538461538461,
988  0.09090909090909091,
989  0.09090909090909091,
990  0.09090909090909091,
991  0.09090909090909091,
992  0.09090909090909091,
993  0.2727272727272727,
994  0.2727272727272727,
995  0.2727272727272727,
996  0.2727272727272727,
997  0.4545454545454545,
998  0.4545454545454545,
999  0.4545454545454545,
1000  0.6363636363636364,
1001  0.6363636363636364,
1002  0.8181818181818182,
1003  0.1111111111111111,
1004  0.1111111111111111,
1005  0.1111111111111111,
1006  0.1111111111111111,
1007  0.3333333333333333,
1008  0.3333333333333333,
1009  0.3333333333333333,
1010  0.5555555555555556,
1011  0.5555555555555556,
1012  0.7777777777777778,
1013  0.1428571428571428,
1014  0.1428571428571428,
1015  0.1428571428571428,
1016  0.4285714285714285,
1017  0.4285714285714285,
1018  0.7142857142857143,
1019  0.2,
1020  0.2,
1021  0.6,
1022  0.3333333333333333};
1023 static const double G_TRI_W286[] = {
1024  2.912193380035668, 2.912193380035668, 2.912193380035668,
1025  2.912193380035668, 2.912193380035668, 2.912193380035668,
1026  2.912193380035668, 2.912193380035668, 2.912193380035668,
1027  2.912193380035668, 2.912193380035668, 2.912193380035668,
1028  2.912193380035668, 2.912193380035668, 2.912193380035668,
1029  2.912193380035668, 2.912193380035668, 2.912193380035668,
1030  2.912193380035668, 2.912193380035668, 2.912193380035668,
1031  2.912193380035668, 2.912193380035668, 2.912193380035668,
1032  2.912193380035668, 2.912193380035668, 2.912193380035668,
1033  2.912193380035668, 2.912193380035668, 2.912193380035668,
1034  2.912193380035668, 2.912193380035668, 2.912193380035668,
1035  2.912193380035668, 2.912193380035668, 2.912193380035668,
1036  2.912193380035668, 2.912193380035668, 2.912193380035668,
1037  2.912193380035668, 2.912193380035668, 2.912193380035668,
1038  2.912193380035668, 2.912193380035668, 2.912193380035668,
1039  2.912193380035668, 2.912193380035668, 2.912193380035668,
1040  2.912193380035668, 2.912193380035668, 2.912193380035668,
1041  2.912193380035668, 2.912193380035668, 2.912193380035668,
1042  2.912193380035668, 2.912193380035668, 2.912193380035668,
1043  2.912193380035668, 2.912193380035668, 2.912193380035668,
1044  2.912193380035668, 2.912193380035668, 2.912193380035668,
1045  2.912193380035668, 2.912193380035668, 2.912193380035668,
1046  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1047  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1048  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1049  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1050  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1051  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1052  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1053  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1054  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1055  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1056  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1057  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1058  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1059  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1060  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1061  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1062  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1063  -9.914451197589852, -9.914451197589852, -9.914451197589852,
1064  -9.914451197589852, 13.33158527957992, 13.33158527957992,
1065  13.33158527957992, 13.33158527957992, 13.33158527957992,
1066  13.33158527957992, 13.33158527957992, 13.33158527957992,
1067  13.33158527957992, 13.33158527957992, 13.33158527957992,
1068  13.33158527957992, 13.33158527957992, 13.33158527957992,
1069  13.33158527957992, 13.33158527957992, 13.33158527957992,
1070  13.33158527957992, 13.33158527957992, 13.33158527957992,
1071  13.33158527957992, 13.33158527957992, 13.33158527957992,
1072  13.33158527957992, 13.33158527957992, 13.33158527957992,
1073  13.33158527957992, 13.33158527957992, 13.33158527957992,
1074  13.33158527957992, 13.33158527957992, 13.33158527957992,
1075  13.33158527957992, 13.33158527957992, 13.33158527957992,
1076  13.33158527957992, 13.33158527957992, 13.33158527957992,
1077  13.33158527957992, 13.33158527957992, 13.33158527957992,
1078  13.33158527957992, 13.33158527957992, 13.33158527957992,
1079  13.33158527957992, -9.027792408986382, -9.027792408986382,
1080  -9.027792408986382, -9.027792408986382, -9.027792408986382,
1081  -9.027792408986382, -9.027792408986382, -9.027792408986382,
1082  -9.027792408986382, -9.027792408986382, -9.027792408986382,
1083  -9.027792408986382, -9.027792408986382, -9.027792408986382,
1084  -9.027792408986382, -9.027792408986382, -9.027792408986382,
1085  -9.027792408986382, -9.027792408986382, -9.027792408986382,
1086  -9.027792408986382, -9.027792408986382, -9.027792408986382,
1087  -9.027792408986382, -9.027792408986382, -9.027792408986382,
1088  -9.027792408986382, -9.027792408986382, -9.027792408986382,
1089  -9.027792408986382, -9.027792408986382, -9.027792408986382,
1090  -9.027792408986382, -9.027792408986382, -9.027792408986382,
1091  -9.027792408986382, 3.258672079964582, 3.258672079964582,
1092  3.258672079964582, 3.258672079964582, 3.258672079964582,
1093  3.258672079964582, 3.258672079964582, 3.258672079964582,
1094  3.258672079964582, 3.258672079964582, 3.258672079964582,
1095  3.258672079964582, 3.258672079964582, 3.258672079964582,
1096  3.258672079964582, 3.258672079964582, 3.258672079964582,
1097  3.258672079964582, 3.258672079964582, 3.258672079964582,
1098  3.258672079964582, 3.258672079964582, 3.258672079964582,
1099  3.258672079964582, 3.258672079964582, 3.258672079964582,
1100  3.258672079964582, 3.258672079964582, -0.6133639040302452,
1101  -0.6133639040302452, -0.6133639040302452, -0.6133639040302452,
1102  -0.6133639040302452, -0.6133639040302452, -0.6133639040302452,
1103  -0.6133639040302452, -0.6133639040302452, -0.6133639040302452,
1104  -0.6133639040302452, -0.6133639040302452, -0.6133639040302452,
1105  -0.6133639040302452, -0.6133639040302452, -0.6133639040302452,
1106  -0.6133639040302452, -0.6133639040302452, -0.6133639040302452,
1107  -0.6133639040302452, -0.6133639040302452, 0.05511571669513555,
1108  0.05511571669513555, 0.05511571669513555, 0.05511571669513555,
1109  0.05511571669513555, 0.05511571669513555, 0.05511571669513555,
1110  0.05511571669513555, 0.05511571669513555, 0.05511571669513555,
1111  0.05511571669513555, 0.05511571669513555, 0.05511571669513555,
1112  0.05511571669513555, 0.05511571669513555, -0.001979122382447095,
1113  -0.001979122382447095, -0.001979122382447095, -0.001979122382447095,
1114  -0.001979122382447095, -0.001979122382447095, -0.001979122382447095,
1115  -0.001979122382447095, -0.001979122382447095, -0.001979122382447095,
1116  2.02054621415273e-05, 2.02054621415273e-05, 2.02054621415273e-05,
1117  2.02054621415273e-05, 2.02054621415273e-05, 2.02054621415273e-05,
1118  -2.874940020535803e-08, -2.874940020535803e-08, -2.874940020535803e-08,
1119  8.829438425435718e-13};
1120 
1121 // TET
1122 static const double G_TET_X1[] = {0.25};
1123 static const double G_TET_Y1[] = {0.25};
1124 static const double G_TET_Z1[] = {0.25};
1125 static const double G_TET_W1[] = {1.};
1126 static const double G_TET_X4[] = {0.1757281246520584, 0.2445310270213291,
1127  0.5556470949048655, 0.0240937534217468};
1128 static const double G_TET_Y4[] = {0.5656137776620919, 0.0501800797762026,
1129  0.1487681308666864, 0.2354380116950194};
1130 static const double G_TET_Z4[] = {0.2180665126782654, 0.5635595064952189,
1131  0.0350112499848832, 0.1833627308416330};
1132 static const double G_TET_W4[] = {0.25, 0.25, 0.25, 0.25};
1133 static const double G_TET_X5[] = {0.25000000000000000, 0.50000000000000000,
1134  0.16666666666666667, 0.16666666666666667,
1135  0.16666666666666667};
1136 static const double G_TET_Y5[] = {0.25000000000000000, 0.16666666666666667,
1137  0.50000000000000000, 0.16666666666666667,
1138  0.16666666666666667};
1139 static const double G_TET_Z5[] = {0.25000000000000000, 0.16666666666666667,
1140  0.16666666666666667, 0.50000000000000000,
1141  0.16666666666666667};
1142 static const double G_TET_W5[] = {-0.80000000000000000, 0.45000000000000000,
1143  0.45000000000000000, 0.45000000000000000,
1144  0.45000000000000000};
1145 static const double G_TET_X10[] = {0.5684305841968444, 0.1438564719343852,
1146  0.1438564719343852, 0.1438564719343852,
1147  0.0000000000000000, 0.5000000000000000,
1148  0.5000000000000000, 0.5000000000000000,
1149  0.0000000000000000, 0.0000000000000000};
1150 static const double G_TET_Y10[] = {0.1438564719343852, 0.1438564719343852,
1151  0.1438564719343852, 0.5684305841968444,
1152  0.5000000000000000, 0.0000000000000000,
1153  0.5000000000000000, 0.0000000000000000,
1154  0.5000000000000000, 0.0000000000000000};
1155 static const double G_TET_Z10[] = {0.1438564719343852, 0.1438564719343852,
1156  0.5684305841968444, 0.1438564719343852,
1157  0.5000000000000000, 0.5000000000000000,
1158  0.0000000000000000, 0.0000000000000000,
1159  0.0000000000000000, 0.5000000000000000};
1160 static const double G_TET_W10[] = {0.2177650698804054, 0.2177650698804054,
1161  0.2177650698804054, 0.2177650698804054,
1162  0.0214899534130631, 0.0214899534130631,
1163  0.0214899534130631, 0.0214899534130631,
1164  0.0214899534130631, 0.0214899534130631};
1165 
1166 static const double G_TET_X45[] = {
1167  0.2500000000000000, 0.6175871903000830, 0.1274709365666390,
1168  0.1274709365666390, 0.1274709365666390, 0.9037635088221031,
1169  0.0320788303926323, 0.0320788303926323, 0.0320788303926323,
1170  0.4502229043567190, 0.0497770956432810, 0.0497770956432810,
1171  0.0497770956432810, 0.4502229043567190, 0.4502229043567190,
1172  0.3162695526014501, 0.1837304473985499, 0.1837304473985499,
1173  0.1837304473985499, 0.3162695526014501, 0.3162695526014501,
1174  0.0229177878448171, 0.2319010893971509, 0.2319010893971509,
1175  0.5132800333608811, 0.2319010893971509, 0.2319010893971509,
1176  0.2319010893971509, 0.0229177878448171, 0.5132800333608811,
1177  0.2319010893971509, 0.0229177878448171, 0.5132800333608811,
1178  0.7303134278075384, 0.0379700484718286, 0.0379700484718286,
1179  0.1937464752488044, 0.0379700484718286, 0.0379700484718286,
1180  0.0379700484718286, 0.7303134278075384, 0.1937464752488044,
1181  0.0379700484718286, 0.7303134278075384, 0.1937464752488044};
1182 static const double G_TET_Y45[] = {
1183  0.2500000000000000, 0.1274709365666390, 0.1274709365666390,
1184  0.1274709365666390, 0.6175871903000830, 0.0320788303926323,
1185  0.0320788303926323, 0.0320788303926323, 0.9037635088221031,
1186  0.0497770956432810, 0.4502229043567190, 0.0497770956432810,
1187  0.4502229043567190, 0.0497770956432810, 0.4502229043567190,
1188  0.1837304473985499, 0.3162695526014501, 0.1837304473985499,
1189  0.3162695526014501, 0.1837304473985499, 0.3162695526014501,
1190  0.2319010893971509, 0.0229177878448171, 0.2319010893971509,
1191  0.2319010893971509, 0.5132800333608811, 0.2319010893971509,
1192  0.0229177878448171, 0.5132800333608811, 0.2319010893971509,
1193  0.5132800333608811, 0.2319010893971509, 0.0229177878448171,
1194  0.0379700484718286, 0.7303134278075384, 0.0379700484718286,
1195  0.0379700484718286, 0.1937464752488044, 0.0379700484718286,
1196  0.7303134278075384, 0.1937464752488044, 0.0379700484718286,
1197  0.1937464752488044, 0.0379700484718286, 0.7303134278075384};
1198 static const double G_TET_Z45[] = {
1199  0.2500000000000000, 0.1274709365666390, 0.1274709365666390,
1200  0.6175871903000830, 0.1274709365666390, 0.0320788303926323,
1201  0.0320788303926323, 0.9037635088221031, 0.0320788303926323,
1202  0.0497770956432810, 0.0497770956432810, 0.4502229043567190,
1203  0.4502229043567190, 0.4502229043567190, 0.0497770956432810,
1204  0.1837304473985499, 0.1837304473985499, 0.3162695526014501,
1205  0.3162695526014501, 0.3162695526014501, 0.1837304473985499,
1206  0.2319010893971509, 0.2319010893971509, 0.0229177878448171,
1207  0.2319010893971509, 0.2319010893971509, 0.5132800333608811,
1208  0.5132800333608811, 0.2319010893971509, 0.0229177878448171,
1209  0.0229177878448171, 0.5132800333608811, 0.2319010893971509,
1210  0.0379700484718286, 0.0379700484718286, 0.7303134278075384,
1211  0.0379700484718286, 0.0379700484718286, 0.1937464752488044,
1212  0.1937464752488044, 0.0379700484718286, 0.7303134278075384,
1213  0.7303134278075384, 0.1937464752488044, 0.0379700484718286};
1214 static const double G_TET_W45[] = {
1215  -0.2359620398477557, 0.0244878963560562, 0.0244878963560562,
1216  0.0244878963560562, 0.0244878963560562, 0.0039485206398261,
1217  0.0039485206398261, 0.0039485206398261, 0.0039485206398261,
1218  0.0263055529507371, 0.0263055529507371, 0.0263055529507371,
1219  0.0263055529507371, 0.0263055529507371, 0.0263055529507371,
1220  0.0829803830550589, 0.0829803830550589, 0.0829803830550589,
1221  0.0829803830550589, 0.0829803830550589, 0.0829803830550589,
1222  0.0254426245481023, 0.0254426245481023, 0.0254426245481023,
1223  0.0254426245481023, 0.0254426245481023, 0.0254426245481023,
1224  0.0254426245481023, 0.0254426245481023, 0.0254426245481023,
1225  0.0254426245481023, 0.0254426245481023, 0.0254426245481023,
1226  0.0134324384376852, 0.0134324384376852, 0.0134324384376852,
1227  0.0134324384376852, 0.0134324384376852, 0.0134324384376852,
1228  0.0134324384376852, 0.0134324384376852, 0.0134324384376852,
1229  0.0134324384376852, 0.0134324384376852, 0.0134324384376852};
1230 static const double NC_TET_X84[] = {
1231  0.1000000000000000, 0.1000000000000000, 0.1000000000000000,
1232  0.7000000000000000, 0.1000000000000000, 0.1000000000000000,
1233  0.1000000000000000, 0.1000000000000000, 0.1000000000000000,
1234  0.1000000000000000, 0.6000000000000000, 0.6000000000000000,
1235  0.6000000000000000, 0.2000000000000000, 0.2000000000000000,
1236  0.2000000000000000, 0.1000000000000000, 0.1000000000000000,
1237  0.1000000000000000, 0.1000000000000000, 0.1000000000000000,
1238  0.1000000000000000, 0.5000000000000000, 0.5000000000000000,
1239  0.5000000000000000, 0.3000000000000000, 0.3000000000000000,
1240  0.3000000000000000, 0.2000000000000000, 0.2000000000000000,
1241  0.2000000000000000, 0.2000000000000000, 0.2000000000000000,
1242  0.2000000000000000, 0.1000000000000000, 0.1000000000000000,
1243  0.1000000000000000, 0.5000000000000000, 0.5000000000000000,
1244  0.5000000000000000, 0.4000000000000000, 0.4000000000000000,
1245  0.4000000000000000, 0.1000000000000000, 0.1000000000000000,
1246  0.1000000000000000, 0.4000000000000000, 0.4000000000000000,
1247  0.4000000000000000, 0.4000000000000000, 0.4000000000000000,
1248  0.4000000000000000, 0.3000000000000000, 0.3000000000000000,
1249  0.3000000000000000, 0.3000000000000000, 0.3000000000000000,
1250  0.3000000000000000, 0.2000000000000000, 0.2000000000000000,
1251  0.2000000000000000, 0.2000000000000000, 0.2000000000000000,
1252  0.2000000000000000, 0.1000000000000000, 0.1000000000000000,
1253  0.1000000000000000, 0.1000000000000000, 0.1000000000000000,
1254  0.1000000000000000, 0.2000000000000000, 0.2000000000000000,
1255  0.2000000000000000, 0.4000000000000000, 0.3000000000000000,
1256  0.3000000000000000, 0.3000000000000000, 0.1000000000000000,
1257  0.3000000000000000, 0.3000000000000000, 0.3000000000000000,
1258  0.2000000000000000, 0.2000000000000000, 0.2000000000000000};
1259 static const double NC_TET_Y84[] = {
1260  0.1000000000000000, 0.1000000000000000, 0.7000000000000000,
1261  0.1000000000000000, 0.1000000000000000, 0.1000000000000000,
1262  0.6000000000000000, 0.6000000000000000, 0.2000000000000000,
1263  0.2000000000000000, 0.1000000000000000, 0.1000000000000000,
1264  0.2000000000000000, 0.1000000000000000, 0.1000000000000000,
1265  0.6000000000000000, 0.1000000000000000, 0.1000000000000000,
1266  0.5000000000000000, 0.5000000000000000, 0.3000000000000000,
1267  0.3000000000000000, 0.1000000000000000, 0.1000000000000000,
1268  0.3000000000000000, 0.1000000000000000, 0.1000000000000000,
1269  0.5000000000000000, 0.2000000000000000, 0.2000000000000000,
1270  0.1000000000000000, 0.1000000000000000, 0.5000000000000000,
1271  0.5000000000000000, 0.2000000000000000, 0.2000000000000000,
1272  0.5000000000000000, 0.2000000000000000, 0.2000000000000000,
1273  0.1000000000000000, 0.4000000000000000, 0.1000000000000000,
1274  0.1000000000000000, 0.4000000000000000, 0.4000000000000000,
1275  0.1000000000000000, 0.3000000000000000, 0.3000000000000000,
1276  0.2000000000000000, 0.2000000000000000, 0.1000000000000000,
1277  0.1000000000000000, 0.4000000000000000, 0.4000000000000000,
1278  0.2000000000000000, 0.2000000000000000, 0.1000000000000000,
1279  0.1000000000000000, 0.4000000000000000, 0.4000000000000000,
1280  0.3000000000000000, 0.3000000000000000, 0.1000000000000000,
1281  0.1000000000000000, 0.4000000000000000, 0.4000000000000000,
1282  0.3000000000000000, 0.3000000000000000, 0.2000000000000000,
1283  0.2000000000000000, 0.2000000000000000, 0.2000000000000000,
1284  0.4000000000000000, 0.2000000000000000, 0.3000000000000000,
1285  0.3000000000000000, 0.1000000000000000, 0.3000000000000000,
1286  0.3000000000000000, 0.2000000000000000, 0.2000000000000000,
1287  0.3000000000000000, 0.3000000000000000, 0.2000000000000000};
1288 static const double NC_TET_Z84[] = {
1289  0.1000000000000000, 0.7000000000000000, 0.1000000000000000,
1290  0.1000000000000000, 0.6000000000000000, 0.2000000000000000,
1291  0.1000000000000000, 0.2000000000000000, 0.1000000000000000,
1292  0.6000000000000000, 0.1000000000000000, 0.2000000000000000,
1293  0.1000000000000000, 0.1000000000000000, 0.6000000000000000,
1294  0.1000000000000000, 0.5000000000000000, 0.3000000000000000,
1295  0.1000000000000000, 0.3000000000000000, 0.1000000000000000,
1296  0.5000000000000000, 0.1000000000000000, 0.3000000000000000,
1297  0.1000000000000000, 0.1000000000000000, 0.5000000000000000,
1298  0.1000000000000000, 0.1000000000000000, 0.5000000000000000,
1299  0.2000000000000000, 0.5000000000000000, 0.2000000000000000,
1300  0.1000000000000000, 0.2000000000000000, 0.5000000000000000,
1301  0.2000000000000000, 0.2000000000000000, 0.1000000000000000,
1302  0.2000000000000000, 0.1000000000000000, 0.4000000000000000,
1303  0.1000000000000000, 0.4000000000000000, 0.1000000000000000,
1304  0.4000000000000000, 0.2000000000000000, 0.1000000000000000,
1305  0.3000000000000000, 0.1000000000000000, 0.3000000000000000,
1306  0.2000000000000000, 0.2000000000000000, 0.1000000000000000,
1307  0.4000000000000000, 0.1000000000000000, 0.4000000000000000,
1308  0.2000000000000000, 0.3000000000000000, 0.1000000000000000,
1309  0.4000000000000000, 0.1000000000000000, 0.4000000000000000,
1310  0.3000000000000000, 0.3000000000000000, 0.2000000000000000,
1311  0.4000000000000000, 0.2000000000000000, 0.4000000000000000,
1312  0.3000000000000000, 0.2000000000000000, 0.4000000000000000,
1313  0.2000000000000000, 0.2000000000000000, 0.3000000000000000,
1314  0.1000000000000000, 0.3000000000000000, 0.3000000000000000,
1315  0.2000000000000000, 0.3000000000000000, 0.2000000000000000,
1316  0.3000000000000000, 0.2000000000000000, 0.3000000000000000};
1317 static const double NC_TET_W84[] = {
1318  0.2843915343915344, 0.2843915343915344, 0.2843915343915344,
1319  0.2843915343915344, -0.3882275132275133, -0.3882275132275133,
1320  -0.3882275132275133, -0.3882275132275133, -0.3882275132275133,
1321  -0.3882275132275133, -0.3882275132275133, -0.3882275132275133,
1322  -0.3882275132275133, -0.3882275132275133, -0.3882275132275133,
1323  -0.3882275132275133, 0.8776455026455027, 0.8776455026455027,
1324  0.8776455026455027, 0.8776455026455027, 0.8776455026455027,
1325  0.8776455026455027, 0.8776455026455027, 0.8776455026455027,
1326  0.8776455026455027, 0.8776455026455027, 0.8776455026455027,
1327  0.8776455026455027, 0.1236772486772487, 0.1236772486772487,
1328  0.1236772486772487, 0.1236772486772487, 0.1236772486772487,
1329  0.1236772486772487, 0.1236772486772487, 0.1236772486772487,
1330  0.1236772486772487, 0.1236772486772487, 0.1236772486772487,
1331  0.1236772486772487, -0.8584656084656085, -0.8584656084656085,
1332  -0.8584656084656085, -0.8584656084656085, -0.8584656084656085,
1333  -0.8584656084656085, -0.2632275132275133, -0.2632275132275133,
1334  -0.2632275132275133, -0.2632275132275133, -0.2632275132275133,
1335  -0.2632275132275133, -0.2632275132275133, -0.2632275132275133,
1336  -0.2632275132275133, -0.2632275132275133, -0.2632275132275133,
1337  -0.2632275132275133, -0.2632275132275133, -0.2632275132275133,
1338  -0.2632275132275133, -0.2632275132275133, -0.2632275132275133,
1339  -0.2632275132275133, -0.2632275132275133, -0.2632275132275133,
1340  -0.2632275132275133, -0.2632275132275133, -0.2632275132275133,
1341  -0.2632275132275133, 0.0145502645502645, 0.0145502645502645,
1342  0.0145502645502645, 0.0145502645502645, 1.0165343915343916,
1343  1.0165343915343916, 1.0165343915343916, 1.0165343915343916,
1344  -0.0251322751322751, -0.0251322751322751, -0.0251322751322751,
1345  -0.0251322751322751, -0.0251322751322751, -0.0251322751322751};
1346 
1347 #ifdef __cplusplus
1348 }
1349 #endif
1350 
1351 #endif
static Index< 'M', 3 > M
static Index< 'n', 3 > n
MatrixDouble invJac
static const double eps
double ShapeVolumeMBTETQ(const double *diffN, const double *coords, int G_DIM, double *G_TET_W)
Definition: fem_tools.c:1005
PetscErrorCode MakeComplexTensor(double *reA, double *imA, __CLPK_doublecomplex *xA)
Definition: fem_tools.c:499
static const double G_TRI_Y28[]
Definition: fem_tools.h:398
void ShapeDiffMBTRIinvJ(double *diffN, double *invJac, double *diffNinvJac)
calculate derivatives of shape functions in space
static const double G_TRI_W3[]
Definition: fem_tools.h:327
PetscErrorCode ShapeDiffMBTET(double *diffN)
calculate derivatives of shape functions
Definition: fem_tools.c:331
static const double G_TRI_W13[]
Definition: fem_tools.h:357
static const double G_TET_X10[]
Definition: fem_tools.h:1145
PetscErrorCode ShapeJacMBTET(double *diffN, const double *coords, double *jac)
calculate jacobian
Definition: fem_tools.c:300
static const double G_TRI_W19[]
Definition: fem_tools.h:375
static const double G_TRI_Y3[]
Definition: fem_tools.h:326
PetscErrorCode GradientOfDeformation(double *diffN, double *dofs, double *F)
calculate gradient of deformation
Definition: fem_tools.c:437
static const double G_TRI_W37[]
Definition: fem_tools.h:440
void print_mat_sym_upper(double *M, int m, int n)
print upper part of the symmetric matrix
PetscErrorCode Grundmann_Moeller_integration_points_2D_TRI(int rule, double *G_TRI_X, double *G_TRI_Y, double *G_TRI_W)
Definition: fem_tools.c:110
static const double G_TRI_Y1[]
Definition: fem_tools.h:323
static const double G_TET_X1[]
Definition: fem_tools.h:1122
static const double G_TRI_X28[]
Definition: fem_tools.h:387
static const double NC_TET_Y84[]
Definition: fem_tools.h:1259
static const double G_TRI_W4[]
Definition: fem_tools.h:335
PetscErrorCode ShapeMBTET_inverse(double *N, double *diffN, const double *elem_coords, const double *glob_coords, double *loc_coords)
calculate local coordinates for given global coordinates
Definition: fem_tools.c:347
static const double G_TRI_Y19[]
Definition: fem_tools.h:369
PetscErrorCode ShapeDiffMBTRIQ(double *diffN, const double *X, const double *Y, const int G_DIM)
Definition: fem_tools.c:807
static const double G_TET_W4[]
Definition: fem_tools.h:1132
PetscErrorCode ShapeFaceDiffNormalMBTRI(double *diffN, const double *coords, double *diff_normal)
calculate derivative of normal in respect to nodal positions
Definition: fem_tools.c:249
PetscErrorCode ShapeFaceBaseMBTRI(double *diffN, const double *coords, double *normal, double *s1, double *s2)
Definition: fem_tools.c:216
static const double G_TET_W45[]
Definition: fem_tools.h:1214
PetscErrorCode ShapeFaceNormalMBTRI_complex(double *diffN, __CLPK_doublecomplex *xcoords, __CLPK_doublecomplex *xnormal)
Definition: fem_tools.c:464
void ShapeJacMBTRI(double *diffN, const double *coords, double *Jac)
calculate jacobioan
PetscErrorCode ShapeFaceNormalMBTRI(double *diffN, const double *coords, double *normal)
Definition: fem_tools.c:241
static const double NC_TET_X84[]
Definition: fem_tools.h:1230
static const double G_TRI_Y37[]
Definition: fem_tools.h:429
PetscErrorCode ShapeDiffMBTETQ(double *diffN, const double x, const double y, const double z)
Definition: fem_tools.c:885
static const double G_TET_W1[]
Definition: fem_tools.h:1125
PetscErrorCode Grundmann_Moeller_integration_points_3D_TET(int rule, double *G_TET_X, double *G_TET_Y, double *G_TET_Z, double *G_TET_W)
Definition: fem_tools.c:155
PetscErrorCode ShapeInvJacVolume(double *jac)
Definition: fem_tools.c:51
static const double G_TRI_Y286[]
Definition: fem_tools.h:737
static const double G_TET_Y1[]
Definition: fem_tools.h:1123
PetscErrorCode ShapeDiffMBTETinvJ(double *diffN, double *invJac, double *diffNinvJac)
calculate shape functions derivatives in space
Definition: fem_tools.c:427
static const double G_TRI_W28[]
Definition: fem_tools.h:406
static const double G_TET_W10[]
Definition: fem_tools.h:1160
static const double G_TET_Z45[]
Definition: fem_tools.h:1198
PetscErrorCode ShapeDiffMBTRI(double *diffN)
calculate derivatives of shape functions
Definition: fem_tools.c:206
static const double G_TRI_X19[]
Definition: fem_tools.h:363
static const double G_TET_Z4[]
Definition: fem_tools.h:1130
static const double G_TRI_X4[]
Definition: fem_tools.h:329
PetscErrorCode ShapeMBTET(double *N, const double *G_X, const double *G_Y, const double *G_Z, int DIM)
calculate shape functions
Definition: fem_tools.c:318
void print_mat(double *M, int m, int n)
print matric M
PetscErrorCode ShapeDiffMBEDGE(double *diffN)
Definition: fem_tools.c:783
static const double G_TRI_W1[]
Definition: fem_tools.h:324
static const double G_TRI_Y4[]
Definition: fem_tools.h:332
PetscErrorCode ShapeMBTETQ_GAUSS(double *N, const double *X, const double *Y, const double *Z, const int G_DIM)
Definition: fem_tools.c:920
static const double G_TET_X45[]
Definition: fem_tools.h:1166
PetscErrorCode DeterminantComplexGradient(__CLPK_doublecomplex *xF, __CLPK_doublecomplex *det_xF)
Definition: fem_tools.c:538
PetscErrorCode Spin(double *spinOmega, double *vecOmega)
calculate spin matrix from vector
Definition: fem_tools.c:558
PetscErrorCode Grundmann_Moeller_integration_points_1D_EDGE(int rule, double *G_TRI_X, double *G_TRI_W)
Compute weights and integration points for edge using Grundmann_Moeller rule.
Definition: fem_tools.c:67
static const double G_TET_Z1[]
Definition: fem_tools.h:1124
static const double G_TET_Y4[]
Definition: fem_tools.h:1128
static const double G_TET_X4[]
Definition: fem_tools.h:1126
PetscErrorCode InvertComplexGradient(__CLPK_doublecomplex *xF)
Definition: fem_tools.c:511
static const double NC_TET_Z84[]
Definition: fem_tools.h:1288
PetscErrorCode ShapeMBTRIQ(double *N, const double *X, const double *Y, const int G_DIM)
Definition: fem_tools.c:792
PetscErrorCode ShapeMBTETQ(double *N, const double x, const double y, const double z)
Definition: fem_tools.c:870
static const double G_TET_Y45[]
Definition: fem_tools.h:1182
static const double NC_TET_W84[]
Definition: fem_tools.h:1317
static const double G_TET_Z10[]
Definition: fem_tools.h:1155
PetscErrorCode ShapeMBEDGE(double *N, const double *G_X, int DIM)
Definition: fem_tools.c:773
double ShapeVolumeMBTET(double *diffN, const double *coords)
calculate TET volume
Definition: fem_tools.c:310
PetscErrorCode ShapeMBTETQ_inverse(double *N, double *diffN, const double *elem_coords, const double *glob_coords, double *loc_coords, const double eps)
Definition: fem_tools.c:1019
PetscErrorCode make_complex_matrix(double *reA, double *imA, __CLPK_doublecomplex *xA)
Compose complex matrix (3x3) from two real matrices.
Definition: fem_tools.c:569
static const double G_TET_Z5[]
Definition: fem_tools.h:1139
PetscErrorCode ShapeDiffMBTETQ_GAUSS(double *diffN, const double *X, const double *Y, const double *Z, const int G_DIM)
Definition: fem_tools.c:939
PetscErrorCode ShapeMBTETQ_detJac_at_Gauss_Points(double *detJac_at_Gauss_Points, const double *diffN, const double *coords, int G_DIM)
Definition: fem_tools.c:991
static const double G_TRI_X13[]
Definition: fem_tools.h:347
double ShapeDetJacVolume(double *jac)
determined of jacobian
Definition: fem_tools.c:34
static const double G_TET_X5[]
Definition: fem_tools.h:1133
static const double G_TET_W5[]
Definition: fem_tools.h:1142
static const double G_TRI_X37[]
Definition: fem_tools.h:418
void ShapeDiffMBTETinvJ_complex(double *diffN, __CLPK_doublecomplex *invJac, __CLPK_doublecomplex *diffNinvJac, CBLAS_TRANSPOSE Trans)
Definition: fem_tools.c:449
PetscErrorCode Base_scale(__CLPK_doublecomplex *xnormal, __CLPK_doublecomplex *xs1, __CLPK_doublecomplex *xs2)
Definition: fem_tools.c:753
static const double G_TRI_X286[]
Definition: fem_tools.h:451
PetscErrorCode ShapeJacMBTETQ(const double *diffN, const double *coords, double *Jac)
Definition: fem_tools.c:979
PetscErrorCode ShapeMBTRI_inverse(double *N, double *diffN, const double *elem_coords, const double *glob_coords, double *loc_coords)
calculate local coordinates of triangle element for given global coordinates in 2D (Assume e....
Definition: fem_tools.c:392
static const double G_TET_Y5[]
Definition: fem_tools.h:1136
static const double G_TRI_W286[]
Definition: fem_tools.h:1023
static const double G_TRI_X7[]
Definition: fem_tools.h:338
void print_mat_complex(__CLPK_doublecomplex *M, int m, int n)
priint complex matrix
PetscErrorCode ShapeMBTRI(double *N, const double *X, const double *Y, const int G_DIM)
calculate shape functions on triangle
Definition: fem_tools.c:194
PetscErrorCode Normal_hierarchical(int order_approx, int *order_edge_approx, int order, int *order_edge, double *diffN, double *diffN_face, double *diffN_edge[], double *dofs, double *dofs_edge[], double *dofs_face, double *idofs, double *idofs_edge[], double *idofs_face, __CLPK_doublecomplex *xnormal, __CLPK_doublecomplex *s1, __CLPK_doublecomplex *s2, int gg)
Complex normal.
Definition: fem_tools.c:581
static const double G_TRI_Y7[]
Definition: fem_tools.h:341
static const double G_TRI_X3[]
Definition: fem_tools.h:325
static const double G_TRI_Y13[]
Definition: fem_tools.h:352
static const double G_TRI_X1[]
Definition: fem_tools.h:322
static const double G_TRI_W7[]
Definition: fem_tools.h:344
static const double G_TET_Y10[]
Definition: fem_tools.h:1150
FTensor::Tensor1< double, SPACE_DIM > normal(FTensor::Tensor1< T1, 3 > &t_coords, FTensor::Tensor1< T2, SPACE_DIM > &t_disp)
Definition: ContactOps.hpp:77
FTensor::Index< 'm', 3 > m
const int N
Definition: speed_test.cpp:3