v0.14.0
tetrahedron_ncc_rule.c
Go to the documentation of this file.
1 # include "tetrahedron_ncc_rule.h"
2 
3 /******************************************************************************/
4 
5 double r8mat_det_4d ( double a[] )
6 
7 /******************************************************************************/
8 /*
9  Purpose:
10 
11  R8MAT_DET_4D computes the determinant of a 4 by 4 R8MAT.
12 
13  Discussion:
14 
15  An R8MAT is a doubly dimensioned array of R8 values, stored as a vector
16  in column-major order.
17 
18  Licensing:
19 
20  This code is distributed under the GNU LGPL license.
21 
22  Modified:
23 
24  15 May 2010
25 
26  Author:
27 
28  John Burkardt
29 
30  Parameters:
31 
32  Input, double A[4*4], the matrix whose determinant is desired.
33 
34  Output, double R8MAT_DET_4D, the determinant of the matrix.
35 */
36 {
37  double det;
38 
39  det =
40  a[0+0*4] * (
41  a[1+1*4] * ( a[2+2*4] * a[3+3*4] - a[2+3*4] * a[3+2*4] )
42  - a[1+2*4] * ( a[2+1*4] * a[3+3*4] - a[2+3*4] * a[3+1*4] )
43  + a[1+3*4] * ( a[2+1*4] * a[3+2*4] - a[2+2*4] * a[3+1*4] ) )
44  - a[0+1*4] * (
45  a[1+0*4] * ( a[2+2*4] * a[3+3*4] - a[2+3*4] * a[3+2*4] )
46  - a[1+2*4] * ( a[2+0*4] * a[3+3*4] - a[2+3*4] * a[3+0*4] )
47  + a[1+3*4] * ( a[2+0*4] * a[3+2*4] - a[2+2*4] * a[3+0*4] ) )
48  + a[0+2*4] * (
49  a[1+0*4] * ( a[2+1*4] * a[3+3*4] - a[2+3*4] * a[3+1*4] )
50  - a[1+1*4] * ( a[2+0*4] * a[3+3*4] - a[2+3*4] * a[3+0*4] )
51  + a[1+3*4] * ( a[2+0*4] * a[3+1*4] - a[2+1*4] * a[3+0*4] ) )
52  - a[0+3*4] * (
53  a[1+0*4] * ( a[2+1*4] * a[3+2*4] - a[2+2*4] * a[3+1*4] )
54  - a[1+1*4] * ( a[2+0*4] * a[3+2*4] - a[2+2*4] * a[3+0*4] )
55  + a[1+2*4] * ( a[2+0*4] * a[3+1*4] - a[2+1*4] * a[3+0*4] ) );
56 
57  return det;
58 }
59 /******************************************************************************/
60 
61 void reference_to_physical_t4 ( double t[], int n, double ref[], double phy[] )
62 
63 /******************************************************************************/
64 /*
65  Purpose:
66 
67  REFERENCE_TO_PHYSICAL_T4 maps T4 reference points to physical points.
68 
69  Discussion:
70 
71  Given the vertices of an order 4 physical tetrahedron and a point
72  (R,S,T) in the reference tetrahedron, the routine computes the value
73  of the corresponding image point (X,Y,Z) in physical space.
74 
75  This routine will also be correct for an order 10 tetrahedron,
76  if the mapping between reference and physical space
77  is linear. This implies, in particular, that the sides of the
78  image tetrahedron are straight, the faces are flat, and
79  the "midside" nodes in the physical tetrahedron are
80  halfway along the edges of the physical tetrahedron.
81 
82  Licensing:
83 
84  This code is distributed under the GNU LGPL license.
85 
86  Modified:
87 
88  31 January 2007
89 
90  Author:
91 
92  John Burkardt
93 
94  Parameters:
95 
96  Input, double T[3*4], the coordinates of the vertices.
97  The vertices are assumed to be the images of (0,0,0), (1,0,0),
98  (0,1,0) and (0,0,1) respectively.
99 
100  Input, int N, the number of objects to transform.
101 
102  Input, double REF[3*N], points in the reference tetrahedron.
103 
104  Output, double PHY[3*N], corresponding points in the
105  physical tetrahedron.
106 */
107 {
108  int i;
109  int j;
110 
111  for ( i = 0; i < 3; i++ )
112  {
113  for ( j = 0; j < n; j++ )
114  {
115  phy[i+j*3] = t[i+0*3] * ( 1.0 - ref[0+j*3] - ref[1+j*3] - ref[2+j*3] )
116  + t[i+1*3] * + ref[0+j*3]
117  + t[i+2*3] * + ref[1+j*3]
118  + t[i+3*3] * + ref[2+j*3];
119  }
120  }
121 
122  return;
123 }
124 /******************************************************************************/
125 
126 int tetrahedron_ncc_degree ( int rule )
127 
128 /******************************************************************************/
129 /*
130  Purpose:
131 
132  TETRAHEDRON_NCC_DEGREE: degree of an NCC rule for the tetrahedron.
133 
134  Licensing:
135 
136  This code is distributed under the GNU LGPL license.
137 
138  Modified:
139 
140  31 January 2007
141 
142  Author:
143 
144  John Burkardt
145 
146  Reference:
147 
148  Peter Silvester,
149  Symmetric Quadrature Formulae for Simplexes,
150  Mathematics of Computation,
151  Volume 24, Number 109, January 1970, pages 95-100.
152 
153  Parameters:
154 
155  Input, int RULE, the index of the rule.
156 
157  Output, int TETRAHEDRON_NCC_DEGREE, the polynomial degree of exactness of
158  the rule.
159 */
160 {
161  int degree;
162 
163  if ( 1 <= rule && rule <= 7 )
164  {
165  degree = rule - 1;
166  }
167  else
168  {
169  degree = -1;
170  fprintf ( stderr, "\n" );
171  fprintf ( stderr, "TETRAHEDRON_NCC_DEGREE - Fatal error!\n" );
172  fprintf ( stderr, " Illegal RULE = %d\n", rule );
173  exit ( 1 );
174  }
175 
176  return degree;
177 }
178 /******************************************************************************/
179 
180 int tetrahedron_ncc_order_num ( int rule )
181 
182 /******************************************************************************/
183 /*
184  Purpose:
185 
186  TETRAHEDRON_NCC_ORDER_NUM: order of an NCC rule for the tetrahedron.
187 
188  Licensing:
189 
190  This code is distributed under the GNU LGPL license.
191 
192  Modified:
193 
194  31 January 2007
195 
196  Author:
197 
198  John Burkardt
199 
200  Reference:
201 
202  Peter Silvester,
203  Symmetric Quadrature Formulae for Simplexes,
204  Mathematics of Computation,
205  Volume 24, Number 109, January 1970, pages 95-100.
206 
207  Parameters:
208 
209  Input, int RULE, the index of the rule.
210 
211  Output, int TETRAHEDRON_NCC_ORDER_NUM, the order (number of points)
212  of the rule.
213 */
214 {
215  int order;
216  int order_num;
217  int *suborder;
218  int suborder_num;
219 
220  suborder_num = tetrahedron_ncc_suborder_num ( rule );
221 
222  suborder = tetrahedron_ncc_suborder ( rule, suborder_num );
223 
224  order_num = 0;
225  for ( order = 0; order < suborder_num; order++ )
226  {
227  order_num = order_num + suborder[order];
228  }
229 
230  free ( suborder );
231 
232  return order_num;
233 }
234 /******************************************************************************/
235 
236 void tetrahedron_ncc_rule ( int rule, int order_num, double xyz[], double w[] )
237 
238 /******************************************************************************/
239 /*
240  Purpose:
241 
242  TETRAHEDRON_NCC_RULE returns the points and weights of an NCC rule.
243 
244  Licensing:
245 
246  This code is distributed under the GNU LGPL license.
247 
248  Modified:
249 
250  31 January 2007
251 
252  Author:
253 
254  John Burkardt
255 
256  Reference:
257 
258  Peter Silvester,
259  Symmetric Quadrature Formulae for Simplexes,
260  Mathematics of Computation,
261  Volume 24, Number 109, January 1970, pages 95-100.
262 
263  Parameters:
264 
265  Input, int RULE, the index of the rule.
266 
267  Input, int ORDER_NUM, the order (number of points) of the rule.
268 
269  Output, double XYZ[3*ORDER_NUM], the points of the rule.
270 
271  Output, double W[ORDER_NUM], the weights of the rule.
272 */
273 {
274  int o;
275  int s;
276  int *suborder;
277  int suborder_num;
278  double *suborder_w;
279  double *suborder_xyz;
280 /*
281  Get the suborder information.
282 */
283  suborder_num = tetrahedron_ncc_suborder_num ( rule );
284 
285  suborder_xyz = ( double * ) malloc ( 4 * suborder_num * sizeof ( double ) );
286  suborder_w = ( double * ) malloc ( suborder_num * sizeof ( double ) );
287 
288  suborder = tetrahedron_ncc_suborder ( rule, suborder_num );
289 
290  tetrahedron_ncc_subrule ( rule, suborder_num, suborder_xyz, suborder_w );
291 /*
292  Expand the suborder information to a full order rule.
293 */
294  o = 0;
295 
296  for ( s = 0; s < suborder_num; s++ )
297  {
298  if ( suborder[s] == 1 )
299  {
300  xyz[0+o*3] = suborder_xyz[0+s*4];
301  xyz[1+o*3] = suborder_xyz[1+s*4];
302  xyz[2+o*3] = suborder_xyz[2+s*4];
303  w[o] = suborder_w[s];
304  o = o + 1;
305  }
306 /*
307  Fourfold symmetry on (A,A,A,B)
308 
309  123 AAA
310  124 AAB
311  142 ABA
312  412 BAA
313 */
314  else if ( suborder[s] == 4 )
315  {
316  xyz[0+o*3] = suborder_xyz[0+s*4];
317  xyz[1+o*3] = suborder_xyz[1+s*4];
318  xyz[2+o*3] = suborder_xyz[2+s*4];
319  w[o] = suborder_w[s];
320  o = o + 1;
321 
322  xyz[0+o*3] = suborder_xyz[0+s*4];
323  xyz[1+o*3] = suborder_xyz[1+s*4];
324  xyz[2+o*3] = suborder_xyz[3+s*4];
325  w[o] = suborder_w[s];
326  o = o + 1;
327 
328  xyz[0+o*3] = suborder_xyz[0+s*4];
329  xyz[1+o*3] = suborder_xyz[3+s*4];
330  xyz[2+o*3] = suborder_xyz[1+s*4];
331  w[o] = suborder_w[s];
332  o = o + 1;
333 
334  xyz[0+o*3] = suborder_xyz[3+s*4];
335  xyz[1+o*3] = suborder_xyz[0+s*4];
336  xyz[2+o*3] = suborder_xyz[1+s*4];
337  w[o] = suborder_w[s];
338  o = o + 1;
339  }
340 /*
341  Sixfold symmetry on (A,A,B,B):
342 
343  123 (A,A,B)
344  132 (A,B,A),
345  134 (A,B,B)
346  312 (B,A,A)
347  314 (B,A,B)
348  341 (B,B,A)
349 */
350  else if ( suborder[s] == 6 )
351  {
352  xyz[0+o*3] = suborder_xyz[0+s*4];
353  xyz[1+o*3] = suborder_xyz[1+s*4];
354  xyz[2+o*3] = suborder_xyz[2+s*4];
355  w[o] = suborder_w[s];
356  o = o + 1;
357 
358  xyz[0+o*3] = suborder_xyz[0+s*4];
359  xyz[1+o*3] = suborder_xyz[2+s*4];
360  xyz[2+o*3] = suborder_xyz[1+s*4];
361  w[o] = suborder_w[s];
362  o = o + 1;
363 
364  xyz[0+o*3] = suborder_xyz[0+s*4];
365  xyz[1+o*3] = suborder_xyz[2+s*4];
366  xyz[2+o*3] = suborder_xyz[3+s*4];
367  w[o] = suborder_w[s];
368  o = o + 1;
369 
370  xyz[0+o*3] = suborder_xyz[2+s*4];
371  xyz[1+o*3] = suborder_xyz[0+s*4];
372  xyz[2+o*3] = suborder_xyz[1+s*4];
373  w[o] = suborder_w[s];
374  o = o + 1;
375 
376  xyz[0+o*3] = suborder_xyz[2+s*4];
377  xyz[1+o*3] = suborder_xyz[0+s*4];
378  xyz[2+o*3] = suborder_xyz[3+s*4];
379  w[o] = suborder_w[s];
380  o = o + 1;
381 
382  xyz[0+o*3] = suborder_xyz[2+s*4];
383  xyz[1+o*3] = suborder_xyz[3+s*4];
384  xyz[2+o*3] = suborder_xyz[0+s*4];
385  w[o] = suborder_w[s];
386  o = o + 1;
387  }
388 /*
389  Twelvefold symmetry on (A,A,B,C):
390 
391  123 (A,A,B)
392  124 (A,A,C)
393  132 (A,B,A)
394  134 (A,B,C)
395  142 (A,C,A)
396  143 (A,C,B)
397  312 (B,A,A)
398  314 (B,A,C)
399  341 (B,C,A)
400  412 (C,A,A)
401  413 (C,A,B)
402  431 (C,B,A)
403 */
404  else if ( suborder[s] == 12 )
405  {
406  xyz[0+o*3] = suborder_xyz[0+s*4];
407  xyz[1+o*3] = suborder_xyz[1+s*4];
408  xyz[2+o*3] = suborder_xyz[2+s*4];
409  w[o] = suborder_w[s];
410  o = o + 1;
411 
412  xyz[0+o*3] = suborder_xyz[0+s*4];
413  xyz[1+o*3] = suborder_xyz[1+s*4];
414  xyz[2+o*3] = suborder_xyz[3+s*4];
415  w[o] = suborder_w[s];
416  o = o + 1;
417 
418  xyz[0+o*3] = suborder_xyz[0+s*4];
419  xyz[1+o*3] = suborder_xyz[2+s*4];
420  xyz[2+o*3] = suborder_xyz[1+s*4];
421  w[o] = suborder_w[s];
422  o = o + 1;
423 
424  xyz[0+o*3] = suborder_xyz[0+s*4];
425  xyz[1+o*3] = suborder_xyz[2+s*4];
426  xyz[2+o*3] = suborder_xyz[3+s*4];
427  w[o] = suborder_w[s];
428  o = o + 1;
429 
430  xyz[0+o*3] = suborder_xyz[0+s*4];
431  xyz[1+o*3] = suborder_xyz[3+s*4];
432  xyz[2+o*3] = suborder_xyz[1+s*4];
433  w[o] = suborder_w[s];
434  o = o + 1;
435 
436  xyz[0+o*3] = suborder_xyz[0+s*4];
437  xyz[1+o*3] = suborder_xyz[3+s*4];
438  xyz[2+o*3] = suborder_xyz[2+s*4];
439  w[o] = suborder_w[s];
440  o = o + 1;
441 
442  xyz[0+o*3] = suborder_xyz[2+s*4];
443  xyz[1+o*3] = suborder_xyz[0+s*4];
444  xyz[2+o*3] = suborder_xyz[1+s*4];
445  w[o] = suborder_w[s];
446  o = o + 1;
447 
448  xyz[0+o*3] = suborder_xyz[2+s*4];
449  xyz[1+o*3] = suborder_xyz[0+s*4];
450  xyz[2+o*3] = suborder_xyz[3+s*4];
451  w[o] = suborder_w[s];
452  o = o + 1;
453 
454  xyz[0+o*3] = suborder_xyz[2+s*4];
455  xyz[1+o*3] = suborder_xyz[3+s*4];
456  xyz[2+o*3] = suborder_xyz[1+s*4];
457  w[o] = suborder_w[s];
458  o = o + 1;
459 
460  xyz[0+o*3] = suborder_xyz[3+s*4];
461  xyz[1+o*3] = suborder_xyz[0+s*4];
462  xyz[2+o*3] = suborder_xyz[1+s*4];
463  w[o] = suborder_w[s];
464  o = o + 1;
465 
466  xyz[0+o*3] = suborder_xyz[3+s*4];
467  xyz[1+o*3] = suborder_xyz[0+s*4];
468  xyz[2+o*3] = suborder_xyz[2+s*4];
469  w[o] = suborder_w[s];
470  o = o + 1;
471 
472  xyz[0+o*3] = suborder_xyz[3+s*4];
473  xyz[1+o*3] = suborder_xyz[2+s*4];
474  xyz[2+o*3] = suborder_xyz[0+s*4];
475  w[o] = suborder_w[s];
476  o = o + 1;
477  }
478 /*
479  24 fold symmetry on (A,B,C,D):
480 
481  123 (A,B,C)
482  124 (A,B,D)
483  132 (A,C,B)
484  134 (A,C,D)
485  142 (A,D,B)
486  143 (A,D,C)
487  213 (B,A,C)
488  214 (B,A,D)
489  231 (B,C,A)
490  234 (B,C,D)
491  241 (B,D,A)
492  243 (B,D,C)
493  312 (C,A,B)
494  314 (C,A,D)
495  321 (C,B,A)
496  324 (C,B,D)
497  341 (C,D,A)
498  342 (C,D,B)
499  412 (D,A,B)
500  413 (D,A,C)
501  421 (D,B,A)
502  423 (D,B,C)
503  431 (D,C,A)
504  432 (D,C,B)
505 */
506  else if ( suborder[s] == 24 )
507  {
508  xyz[0+o*3] = suborder_xyz[0+s*4];
509  xyz[1+o*3] = suborder_xyz[1+s*4];
510  xyz[2+o*3] = suborder_xyz[2+s*4];
511  w[o] = suborder_w[s];
512  o = o + 1;
513 
514  xyz[0+o*3] = suborder_xyz[0+s*4];
515  xyz[1+o*3] = suborder_xyz[1+s*4];
516  xyz[2+o*3] = suborder_xyz[3+s*4];
517  w[o] = suborder_w[s];
518  o = o + 1;
519 
520  xyz[0+o*3] = suborder_xyz[0+s*4];
521  xyz[1+o*3] = suborder_xyz[2+s*4];
522  xyz[2+o*3] = suborder_xyz[1+s*4];
523  w[o] = suborder_w[s];
524  o = o + 1;
525 
526  xyz[0+o*3] = suborder_xyz[0+s*4];
527  xyz[1+o*3] = suborder_xyz[2+s*4];
528  xyz[2+o*3] = suborder_xyz[3+s*4];
529  w[o] = suborder_w[s];
530  o = o + 1;
531 
532  xyz[0+o*3] = suborder_xyz[0+s*4];
533  xyz[1+o*3] = suborder_xyz[3+s*4];
534  xyz[2+o*3] = suborder_xyz[1+s*4];
535  w[o] = suborder_w[s];
536  o = o + 1;
537 
538  xyz[0+o*3] = suborder_xyz[0+s*4];
539  xyz[1+o*3] = suborder_xyz[3+s*4];
540  xyz[2+o*3] = suborder_xyz[2+s*4];
541  w[o] = suborder_w[s];
542  o = o + 1;
543 
544  xyz[0+o*3] = suborder_xyz[1+s*4];
545  xyz[1+o*3] = suborder_xyz[0+s*4];
546  xyz[2+o*3] = suborder_xyz[3+s*4];
547  w[o] = suborder_w[s];
548  o = o + 1;
549 
550  xyz[0+o*3] = suborder_xyz[1+s*4];
551  xyz[1+o*3] = suborder_xyz[0+s*4];
552  xyz[2+o*3] = suborder_xyz[4+s*4];
553  w[o] = suborder_w[s];
554  o = o + 1;
555 
556  xyz[0+o*3] = suborder_xyz[1+s*4];
557  xyz[1+o*3] = suborder_xyz[2+s*4];
558  xyz[2+o*3] = suborder_xyz[0+s*4];
559  w[o] = suborder_w[s];
560  o = o + 1;
561 
562  xyz[0+o*3] = suborder_xyz[1+s*4];
563  xyz[1+o*3] = suborder_xyz[2+s*4];
564  xyz[2+o*3] = suborder_xyz[3+s*4];
565  w[o] = suborder_w[s];
566  o = o + 1;
567 
568  xyz[0+o*3] = suborder_xyz[1+s*4];
569  xyz[1+o*3] = suborder_xyz[3+s*4];
570  xyz[2+o*3] = suborder_xyz[0+s*4];
571  w[o] = suborder_w[s];
572  o = o + 1;
573 
574  xyz[0+o*3] = suborder_xyz[1+s*4];
575  xyz[1+o*3] = suborder_xyz[3+s*4];
576  xyz[2+o*3] = suborder_xyz[2+s*4];
577  w[o] = suborder_w[s];
578  o = o + 1;
579 
580  xyz[0+o*3] = suborder_xyz[2+s*4];
581  xyz[1+o*3] = suborder_xyz[0+s*4];
582  xyz[2+o*3] = suborder_xyz[1+s*4];
583  w[o] = suborder_w[s];
584  o = o + 1;
585 
586  xyz[0+o*3] = suborder_xyz[2+s*4];
587  xyz[1+o*3] = suborder_xyz[0+s*4];
588  xyz[2+o*3] = suborder_xyz[3+s*4];
589  w[o] = suborder_w[s];
590  o = o + 1;
591 
592  xyz[0+o*3] = suborder_xyz[2+s*4];
593  xyz[1+o*3] = suborder_xyz[1+s*4];
594  xyz[2+o*3] = suborder_xyz[0+s*4];
595  w[o] = suborder_w[s];
596  o = o + 1;
597 
598  xyz[0+o*3] = suborder_xyz[2+s*4];
599  xyz[1+o*3] = suborder_xyz[1+s*4];
600  xyz[2+o*3] = suborder_xyz[3+s*4];
601  w[o] = suborder_w[s];
602  o = o + 1;
603 
604  xyz[0+o*3] = suborder_xyz[2+s*4];
605  xyz[1+o*3] = suborder_xyz[3+s*4];
606  xyz[2+o*3] = suborder_xyz[0+s*4];
607  w[o] = suborder_w[s];
608  o = o + 1;
609 
610  xyz[0+o*3] = suborder_xyz[2+s*4];
611  xyz[1+o*3] = suborder_xyz[3+s*4];
612  xyz[2+o*3] = suborder_xyz[1+s*4];
613  w[o] = suborder_w[s];
614  o = o + 1;
615 
616  xyz[0+o*3] = suborder_xyz[3+s*4];
617  xyz[1+o*3] = suborder_xyz[0+s*4];
618  xyz[2+o*3] = suborder_xyz[1+s*4];
619  w[o] = suborder_w[s];
620  o = o + 1;
621 
622  xyz[0+o*3] = suborder_xyz[3+s*4];
623  xyz[1+o*3] = suborder_xyz[0+s*4];
624  xyz[2+o*3] = suborder_xyz[2+s*4];
625  w[o] = suborder_w[s];
626  o = o + 1;
627 
628  xyz[0+o*3] = suborder_xyz[3+s*4];
629  xyz[1+o*3] = suborder_xyz[1+s*4];
630  xyz[2+o*3] = suborder_xyz[0+s*4];
631  w[o] = suborder_w[s];
632  o = o + 1;
633 
634  xyz[0+o*3] = suborder_xyz[3+s*4];
635  xyz[1+o*3] = suborder_xyz[1+s*4];
636  xyz[2+o*3] = suborder_xyz[2+s*4];
637  w[o] = suborder_w[s];
638  o = o + 1;
639 
640  xyz[0+o*3] = suborder_xyz[3+s*4];
641  xyz[1+o*3] = suborder_xyz[2+s*4];
642  xyz[2+o*3] = suborder_xyz[0+s*4];
643  w[o] = suborder_w[s];
644  o = o + 1;
645 
646  xyz[0+o*3] = suborder_xyz[3+s*4];
647  xyz[1+o*3] = suborder_xyz[2+s*4];
648  xyz[2+o*3] = suborder_xyz[1+s*4];
649  w[o] = suborder_w[s];
650  o = o + 1;
651  }
652  else
653  {
654  fprintf ( stderr, "\n" );
655  fprintf ( stderr, "TETRAHEDRON_NCC_RULE - Fatal error!\n" );
656  fprintf ( stderr, " Illegal SUBORDER(%d) = %d\n", s, suborder[s] );
657  exit ( 1 );
658  }
659  }
660 
661  free ( suborder );
662  free ( suborder_xyz );
663  free ( suborder_w );
664 
665  return;
666 }
667 /******************************************************************************/
668 
669 int tetrahedron_ncc_rule_num ( )
670 
671 /******************************************************************************/
672 /*
673  Purpose:
674 
675  TETRAHEDRON_NCC_RULE_NUM returns the number of NCC rules available.
676 
677  Licensing:
678 
679  This code is distributed under the GNU LGPL license.
680 
681  Modified:
682 
683  31 January 2007
684 
685  Author:
686 
687  John Burkardt
688 
689  Reference:
690 
691  Peter Silvester,
692  Symmetric Quadrature Formulae for Simplexes,
693  Mathematics of Computation,
694  Volume 24, Number 109, January 1970, pages 95-100.
695 
696  Parameters:
697 
698  Output, int TETRAHEDRON_NCC_RULE_NUM, the number of rules available.
699 */
700 {
701  int rule_num;
702 
703  rule_num = 7;
704 
705  return rule_num;
706 }
707 /******************************************************************************/
708 
709 int *tetrahedron_ncc_suborder ( int rule, int suborder_num )
710 
711 /******************************************************************************/
712 /*
713  Purpose:
714 
715  TETRAHEDRON_NCC_SUBORDER returns the suborders for an NCC rule.
716 
717  Licensing:
718 
719  This code is distributed under the GNU LGPL license.
720 
721  Modified:
722 
723  31 January 2007
724 
725  Author:
726 
727  John Burkardt
728 
729  Reference:
730 
731  Peter Silvester,
732  Symmetric Quadrature Formulae for Simplexes,
733  Mathematics of Computation,
734  Volume 24, Number 109, January 1970, pages 95-100.
735 
736  Parameters:
737 
738  Input, int RULE, the index of the rule.
739 
740  Input, int SUBORDER_NUM, the number of suborders of the rule.
741 
742  Output, int TETRAHEDRON_NCC_SUBORDER[SUBORDER_NUM],
743  the suborders of the rule.
744 */
745 {
746  int *suborder;
747 
748  suborder = ( int * ) malloc ( suborder_num * sizeof ( int ) );
749 
750  if ( rule == 1 )
751  {
752  suborder[0] = 1;
753  }
754  else if ( rule == 2 )
755  {
756  suborder[0] = 4;
757  }
758  else if ( rule == 3 )
759  {
760  suborder[0] = 4;
761  suborder[1] = 6;
762  }
763  else if ( rule == 4 )
764  {
765  suborder[0] = 4;
766  suborder[1] = 12;
767  suborder[2] = 4;
768  }
769  else if ( rule == 5 )
770  {
771  suborder[0] = 4;
772  suborder[1] = 12;
773  suborder[2] = 6;
774  suborder[3] = 12;
775  suborder[4] = 1;
776  }
777  else if ( rule == 6 )
778  {
779  suborder[0] = 4;
780  suborder[1] = 12;
781  suborder[2] = 12;
782  suborder[3] = 12;
783  suborder[4] = 12;
784  suborder[5] = 4;
785  }
786  else if ( rule == 7 )
787  {
788  suborder[0] = 4;
789  suborder[1] = 12;
790  suborder[2] = 12;
791  suborder[3] = 12;
792  suborder[4] = 6;
793  suborder[5] = 24;
794  suborder[6] = 4;
795  suborder[7] = 4;
796  suborder[8] = 6;
797  }
798  else
799  {
800  fprintf ( stderr, "\n" );
801  fprintf ( stderr, "TETRAHEDRON_NCC_SUBORDER - Fatal error!\n" );
802  fprintf ( stderr, " Illegal RULE = %d\n", rule );
803  exit ( 1 );
804  }
805 
806  return suborder;
807 }
808 /******************************************************************************/
809 
810 int tetrahedron_ncc_suborder_num ( int rule )
811 
812 /******************************************************************************/
813 /*
814  Purpose:
815 
816  TETRAHEDRON_NCC_SUBORDER_NUM returns the number of suborders for an NCC rule.
817 
818  Licensing:
819 
820  This code is distributed under the GNU LGPL license.
821 
822  Modified:
823 
824  31 January 2007
825 
826  Author:
827 
828  John Burkardt
829 
830  Reference:
831 
832  Peter Silvester,
833  Symmetric Quadrature Formulae for Simplexes,
834  Mathematics of Computation,
835  Volume 24, Number 109, January 1970, pages 95-100.
836 
837  Parameters:
838 
839  Input, int RULE, the index of the rule.
840 
841  Output, int TETRAHEDRON_NCC_SUBORDER_NUM, the number of suborders
842  of the rule.
843 */
844 {
845  int suborder_num;
846 
847  if ( rule == 1 )
848  {
849  suborder_num = 1;
850  }
851  else if ( rule == 2 )
852  {
853  suborder_num = 1;
854  }
855  else if ( rule == 3 )
856  {
857  suborder_num = 2;
858  }
859  else if ( rule == 4 )
860  {
861  suborder_num = 3;
862  }
863  else if ( rule == 5 )
864  {
865  suborder_num = 5;
866  }
867  else if ( rule == 6 )
868  {
869  suborder_num = 6;
870  }
871  else if ( rule == 7 )
872  {
873  suborder_num = 9;
874  }
875  else
876  {
877  suborder_num = -1;
878  fprintf ( stderr, "\n" );
879  fprintf ( stderr, "TETRAHEDRON_NCC_SUBORDER_NUM - Fatal error!\n" );
880  fprintf ( stderr, " Illegal RULE = %d\n", rule );
881  exit ( 1 );
882  }
883 
884  return suborder_num;
885 }
886 /******************************************************************************/
887 
888 void tetrahedron_ncc_subrule ( int rule, int suborder_num,
889  double suborder_xyz[], double suborder_w[] )
890 
891 /******************************************************************************/
892 /*
893  Purpose:
894 
895  TETRAHEDRON_NCC_SUBRULE returns a compressed NCC rule.
896 
897  Licensing:
898 
899  This code is distributed under the GNU LGPL license.
900 
901  Modified:
902 
903  31 January 2007
904 
905  Author:
906 
907  John Burkardt
908 
909  Reference:
910 
911  Peter Silvester,
912  Symmetric Quadrature Formulae for Simplexes,
913  Mathematics of Computation,
914  Volume 24, Number 109, January 1970, pages 95-100.
915 
916  Parameters:
917 
918  Input, int RULE, the index of the rule.
919 
920  Input, int SUBORDER_NUM, the number of suborders of the rule.
921 
922  Output, double SUBORDER_XYZ[4*SUBORDER_NUM],
923  the barycentric coordinates of the abscissas.
924 
925  Output, double SUBORDER_W[SUBORDER_NUM], the suborder weights.
926 */
927 {
928  int i;
929  int s;
930  int suborder_w_d;
931  int *suborder_w_n;
932  int suborder_xyz_d;
933  int *suborder_xyz_n;
934 
935  suborder_xyz_n = ( int * ) malloc ( 4 * suborder_num * sizeof ( int ) );
936  suborder_w_n = ( int * ) malloc ( suborder_num * sizeof ( int ) );
937 
938  if ( rule == 1 )
939  {
940  tetrahedron_ncc_subrule_01 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
941  suborder_w_n, &suborder_w_d );
942  }
943  else if ( rule == 2 )
944  {
945  tetrahedron_ncc_subrule_02 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
946  suborder_w_n, &suborder_w_d );
947  }
948  else if ( rule == 3 )
949  {
950  tetrahedron_ncc_subrule_03 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
951  suborder_w_n, &suborder_w_d );
952  }
953  else if ( rule == 4 )
954  {
955  tetrahedron_ncc_subrule_04 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
956  suborder_w_n, &suborder_w_d );
957  }
958  else if ( rule == 5 )
959  {
960  tetrahedron_ncc_subrule_05 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
961  suborder_w_n, &suborder_w_d );
962  }
963  else if ( rule == 6 )
964  {
965  tetrahedron_ncc_subrule_06 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
966  suborder_w_n, &suborder_w_d );
967  }
968  else if ( rule == 7 )
969  {
970  tetrahedron_ncc_subrule_07 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
971  suborder_w_n, &suborder_w_d );
972  }
973  else
974  {
975  fprintf ( stderr, "\n" );
976  fprintf ( stderr, "TETRAHEDRON_NCC_SUBRULE - Fatal error!\n" );
977  fprintf ( stderr, " Illegal RULE = %d\n", rule );
978  exit ( 1 );
979  }
980 
981  for ( s = 0; s < suborder_num; s++ )
982  {
983  for ( i = 0; i < 4; i++ )
984  {
985  suborder_xyz[i+s*4] =
986  ( double ) ( suborder_xyz_n[i+s*4] )
987  / ( double ) ( suborder_xyz_d );
988  }
989  }
990  for ( s = 0; s < suborder_num; s++ )
991  {
992  suborder_w[s] = ( double ) suborder_w_n[s] / ( double ) suborder_w_d;
993  }
994 
995  free ( suborder_w_n );
996  free ( suborder_xyz_n );
997 
998  return;
999 }
1000 /******************************************************************************/
1001 
1002 void tetrahedron_ncc_subrule_01 ( int suborder_num, int suborder_xyz_n[],
1003  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1004 
1005 /******************************************************************************/
1006 /*
1007  Purpose:
1008 
1009  TETRAHEDRON_NCC_SUBRULE_01 returns a compressed NCC rule 1.
1010 
1011  Licensing:
1012 
1013  This code is distributed under the GNU LGPL license.
1014 
1015  Modified:
1016 
1017  31 January 2007
1018 
1019  Author:
1020 
1021  John Burkardt
1022 
1023  Reference:
1024 
1025  Peter Silvester,
1026  Symmetric Quadrature Formulae for Simplexes,
1027  Mathematics of Computation,
1028  Volume 24, Number 109, January 1970, pages 95-100.
1029 
1030  Parameters:
1031 
1032  Input, int SUBORDER_NUM, the number of suborders of the rule.
1033 
1034  Output, int SUBORDER_XYZ_N[4*SUBORDER_NUM],
1035  the numerators of the barycentric coordinates of the abscissas.
1036 
1037  Output, int *SUBORDER_XYZ_D,
1038  the denominator of the barycentric coordinates of the abscissas.
1039 
1040  Output, int SUBORDER_W_N[SUBORDER_NUM],
1041  the numerator of the suborder weights.
1042 
1043  Output, int SUBORDER_W_D,
1044  the denominator of the suborder weights.
1045 */
1046 {
1047  int i;
1048  int s;
1049  int suborder_xyz_n_01[4*1] = {
1050  1, 1, 1, 1
1051  };
1052  int suborder_xyz_d_01 = 4;
1053  int suborder_w_n_01[1] = { 1 };
1054  int suborder_w_d_01 = 1;
1055 
1056  for ( s = 0; s < suborder_num; s++ )
1057  {
1058  for ( i = 0; i < 4; i++ )
1059  {
1060  suborder_xyz_n[i+s*4] = suborder_xyz_n_01[i+s*4];
1061  }
1062  }
1063  *suborder_xyz_d = suborder_xyz_d_01;
1064 
1065  for ( s = 0; s < suborder_num; s++ )
1066  {
1067  suborder_w_n[s] = suborder_w_n_01[s];
1068  }
1069  *suborder_w_d = suborder_w_d_01;
1070 
1071  return;
1072 }
1073 /******************************************************************************/
1074 
1075 void tetrahedron_ncc_subrule_02 ( int suborder_num, int suborder_xyz_n[],
1076  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1077 
1078 /******************************************************************************/
1079 /*
1080  Purpose:
1081 
1082  TETRAHEDRON_NCC_SUBRULE_02 returns a compressed NCC rule 2.
1083 
1084  Licensing:
1085 
1086  This code is distributed under the GNU LGPL license.
1087 
1088  Modified:
1089 
1090  31 January 2007
1091 
1092  Author:
1093 
1094  John Burkardt
1095 
1096  Reference:
1097 
1098  Peter Silvester,
1099  Symmetric Quadrature Formulae for Simplexes,
1100  Mathematics of Computation,
1101  Volume 24, Number 109, January 1970, pages 95-100.
1102 
1103  Parameters:
1104 
1105  Input, int SUBORDER_NUM, the number of suborders of the rule.
1106 
1107  Output, int SUBORDER_XYZ_N[4*SUBORDER_NUM],
1108  the numerators of the barycentric coordinates of the abscissas.
1109 
1110  Output, int *SUBORDER_XYZ_D,
1111  the denominator of the barycentric coordinates of the abscissas.
1112 
1113  Output, int SUBORDER_W_N[SUBORDER_NUM],
1114  the numerator of the suborder weights.
1115 
1116  Output, int SUBORDER_W_D,
1117  the denominator of the suborder weights.
1118 */
1119 {
1120  int i;
1121  int s;
1122  int suborder_xyz_n_02[4*1] = {
1123  0, 0, 0, 1
1124  };
1125  int suborder_xyz_d_02 = 1;
1126  int suborder_w_n_02[1] = { 1 };
1127  int suborder_w_d_02 = 4;
1128 
1129  for ( s = 0; s < suborder_num; s++ )
1130  {
1131  for ( i = 0; i < 4; i++ )
1132  {
1133  suborder_xyz_n[i+s*4] = suborder_xyz_n_02[i+s*4];
1134  }
1135  }
1136  *suborder_xyz_d = suborder_xyz_d_02;
1137 
1138  for ( s = 0; s < suborder_num; s++ )
1139  {
1140  suborder_w_n[s] = suborder_w_n_02[s];
1141  }
1142  *suborder_w_d = suborder_w_d_02;
1143 
1144  return;
1145 }
1146 /******************************************************************************/
1147 
1148 void tetrahedron_ncc_subrule_03 ( int suborder_num, int suborder_xyz_n[],
1149  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1150 
1151 /******************************************************************************/
1152 /*
1153  Purpose:
1154 
1155  TETRAHEDRON_NCC_SUBRULE_03 returns a compressed NCC rule 3.
1156 
1157  Licensing:
1158 
1159  This code is distributed under the GNU LGPL license.
1160 
1161  Modified:
1162 
1163  31 January 2007
1164 
1165  Author:
1166 
1167  John Burkardt
1168 
1169  Reference:
1170 
1171  Peter Silvester,
1172  Symmetric Quadrature Formulae for Simplexes,
1173  Mathematics of Computation,
1174  Volume 24, Number 109, January 1970, pages 95-100.
1175 
1176  Parameters:
1177 
1178  Input, int SUBORDER_NUM, the number of suborders of the rule.
1179 
1180  Output, int SUBORDER_XYZ_N[4*SUBORDER_NUM],
1181  the numerators of the barycentric coordinates of the abscissas.
1182 
1183  Output, int *SUBORDER_XYZ_D,
1184  the denominator of the barycentric coordinates of the abscissas.
1185 
1186  Output, int SUBORDER_W_N[SUBORDER_NUM],
1187  the numerator of the suborder weights.
1188 
1189  Output, int SUBORDER_W_D,
1190  the denominator of the suborder weights.
1191 */
1192 {
1193  int i;
1194  int s;
1195  int suborder_xyz_n_03[4*2] = {
1196  0, 0, 0, 2,
1197  1, 1, 0, 0
1198  };
1199  int suborder_xyz_d_03 = 2;
1200  int suborder_w_n_03[2] = { -1, 4 };
1201  int suborder_w_d_03 = 20;
1202 
1203  for ( s = 0; s < suborder_num; s++ )
1204  {
1205  for ( i = 0; i < 4; i++ )
1206  {
1207  suborder_xyz_n[i+s*4] = suborder_xyz_n_03[i+s*4];
1208  }
1209  }
1210  *suborder_xyz_d = suborder_xyz_d_03;
1211 
1212  for ( s = 0; s < suborder_num; s++ )
1213  {
1214  suborder_w_n[s] = suborder_w_n_03[s];
1215  }
1216  *suborder_w_d = suborder_w_d_03;
1217 
1218  return;
1219 }
1220 /******************************************************************************/
1221 
1222 void tetrahedron_ncc_subrule_04 ( int suborder_num, int suborder_xyz_n[],
1223  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1224 
1225 /******************************************************************************/
1226 /*
1227  Purpose:
1228 
1229  TETRAHEDRON_NCC_SUBRULE_04 returns a compressed NCC rule 4.
1230 
1231  Licensing:
1232 
1233  This code is distributed under the GNU LGPL license.
1234 
1235  Modified:
1236 
1237  31 January 2007
1238 
1239  Author:
1240 
1241  John Burkardt
1242 
1243  Reference:
1244 
1245  Peter Silvester,
1246  Symmetric Quadrature Formulae for Simplexes,
1247  Mathematics of Computation,
1248  Volume 24, Number 109, January 1970, pages 95-100.
1249 
1250  Parameters:
1251 
1252  Input, int SUBORDER_NUM, the number of suborders of the rule.
1253 
1254  Output, int SUBORDER_XYZ_N[4*SUBORDER_NUM],
1255  the numerators of the barycentric coordinates of the abscissas.
1256 
1257  Output, int *SUBORDER_XYZ_D,
1258  the denominator of the barycentric coordinates of the abscissas.
1259 
1260  Output, int SUBORDER_W_N[SUBORDER_NUM],
1261  the numerator of the suborder weights.
1262 
1263  Output, int SUBORDER_W_D,
1264  the denominator of the suborder weights.
1265 */
1266 {
1267  int i;
1268  int s;
1269  int suborder_xyz_n_04[4*3] = {
1270  0, 0, 0, 3,
1271  0, 0, 1, 2,
1272  1, 1, 1, 0
1273  };
1274  int suborder_xyz_d_04 = 3;
1275  int suborder_w_n_04[3] = { 1, 0, 9 };
1276  int suborder_w_d_04 = 40;
1277 
1278  for ( s = 0; s < suborder_num; s++ )
1279  {
1280  for ( i = 0; i < 4; i++ )
1281  {
1282  suborder_xyz_n[i+s*4] = suborder_xyz_n_04[i+s*4];
1283  }
1284  }
1285  *suborder_xyz_d = suborder_xyz_d_04;
1286 
1287  for ( s = 0; s < suborder_num; s++ )
1288  {
1289  suborder_w_n[s] = suborder_w_n_04[s];
1290  }
1291  *suborder_w_d = suborder_w_d_04;
1292 
1293  return;
1294 }
1295 /******************************************************************************/
1296 
1297 void tetrahedron_ncc_subrule_05 ( int suborder_num, int suborder_xyz_n[],
1298  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1299 
1300 /******************************************************************************/
1301 /*
1302  Purpose:
1303 
1304  TETRAHEDRON_NCC_SUBRULE_05 returns a compressed NCC rule 5.
1305 
1306  Licensing:
1307 
1308  This code is distributed under the GNU LGPL license.
1309 
1310  Modified:
1311 
1312  31 January 2007
1313 
1314  Author:
1315 
1316  John Burkardt
1317 
1318  Reference:
1319 
1320  Peter Silvester,
1321  Symmetric Quadrature Formulae for Simplexes,
1322  Mathematics of Computation,
1323  Volume 24, Number 109, January 1970, pages 95-100.
1324 
1325  Parameters:
1326 
1327  Input, int SUBORDER_NUM, the number of suborders of the rule.
1328 
1329  Output, int SUBORDER_XYZ_N[4*SUBORDER_NUM],
1330  the numerators of the barycentric coordinates of the abscissas.
1331 
1332  Output, int *SUBORDER_XYZ_D,
1333  the denominator of the barycentric coordinates of the abscissas.
1334 
1335  Output, int SUBORDER_W_N[SUBORDER_NUM],
1336  the numerator of the suborder weights.
1337 
1338  Output, int SUBORDER_W_D,
1339  the denominator of the suborder weights.
1340 */
1341 {
1342  int i;
1343  int s;
1344  int suborder_xyz_n_05[4*5] = {
1345  0, 0, 0, 4,
1346  0, 0, 3, 1,
1347  2, 2, 0, 0,
1348  1, 1, 0, 2,
1349  1, 1, 1, 1
1350  };
1351  int suborder_xyz_d_05 = 4;
1352  int suborder_w_n_05[5] = { -5, 16, -12, 16, 128 };
1353  int suborder_w_d_05 = 420;
1354 
1355  for ( s = 0; s < suborder_num; s++ )
1356  {
1357  for ( i = 0; i < 4; i++ )
1358  {
1359  suborder_xyz_n[i+s*4] = suborder_xyz_n_05[i+s*4];
1360  }
1361  }
1362  *suborder_xyz_d = suborder_xyz_d_05;
1363 
1364  for ( s = 0; s < suborder_num; s++ )
1365  {
1366  suborder_w_n[s] = suborder_w_n_05[s];
1367  }
1368  *suborder_w_d = suborder_w_d_05;
1369 
1370  return;
1371 }
1372 /******************************************************************************/
1373 
1374 void tetrahedron_ncc_subrule_06 ( int suborder_num, int suborder_xyz_n[],
1375  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1376 
1377 /******************************************************************************/
1378 /*
1379  Purpose:
1380 
1381  TETRAHEDRON_NCC_SUBRULE_06 returns a compressed NCC rule 6.
1382 
1383  Licensing:
1384 
1385  This code is distributed under the GNU LGPL license.
1386 
1387  Modified:
1388 
1389  31 January 2007
1390 
1391  Author:
1392 
1393  John Burkardt
1394 
1395  Reference:
1396 
1397  Peter Silvester,
1398  Symmetric Quadrature Formulae for Simplexes,
1399  Mathematics of Computation,
1400  Volume 24, Number 109, January 1970, pages 95-100.
1401 
1402  Parameters:
1403 
1404  Input, int SUBORDER_NUM, the number of suborders of the rule.
1405 
1406  Output, int SUBORDER_XYZ_N[4*SUBORDER_NUM],
1407  the numerators of the barycentric coordinates of the abscissas.
1408 
1409  Output, int *SUBORDER_XYZ_D,
1410  the denominator of the barycentric coordinates of the abscissas.
1411 
1412  Output, int SUBORDER_W_N[SUBORDER_NUM],
1413  the numerator of the suborder weights.
1414 
1415  Output, int SUBORDER_W_D,
1416  the denominator of the suborder weights.
1417 */
1418 {
1419  int i;
1420  int s;
1421  int suborder_xyz_n_06[4*6] = {
1422  0, 0, 0, 5,
1423  0, 0, 4, 1,
1424  0, 0, 3, 2,
1425  1, 1, 0, 3,
1426  2, 2, 1, 0,
1427  1, 1, 1, 2
1428  };
1429  int suborder_xyz_d_06 = 5;
1430  int suborder_w_n_06[6] = { 33, -35, 35, 275, -75, 375 };
1431  int suborder_w_d_06 = 4032;
1432 
1433  for ( s = 0; s < suborder_num; s++ )
1434  {
1435  for ( i = 0; i < 4; i++ )
1436  {
1437  suborder_xyz_n[i+s*4] = suborder_xyz_n_06[i+s*4];
1438  }
1439  }
1440  *suborder_xyz_d = suborder_xyz_d_06;
1441 
1442  for ( s = 0; s < suborder_num; s++ )
1443  {
1444  suborder_w_n[s] = suborder_w_n_06[s];
1445  }
1446  *suborder_w_d = suborder_w_d_06;
1447 
1448  return;
1449 }
1450 /******************************************************************************/
1451 
1452 void tetrahedron_ncc_subrule_07 ( int suborder_num, int suborder_xyz_n[],
1453  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1454 
1455 /******************************************************************************/
1456 /*
1457  Purpose:
1458 
1459  TETRAHEDRON_NCC_SUBRULE_07 returns a compressed NCC rule 7.
1460 
1461  Licensing:
1462 
1463  This code is distributed under the GNU LGPL license.
1464 
1465  Modified:
1466 
1467  31 January 2007
1468 
1469  Author:
1470 
1471  John Burkardt
1472 
1473  Reference:
1474 
1475  Peter Silvester,
1476  Symmetric Quadrature Formulae for Simplexes,
1477  Mathematics of Computation,
1478  Volume 24, Number 109, January 1970, pages 95-100.
1479 
1480  Parameters:
1481 
1482  Input, int SUBORDER_NUM, the number of suborders of the rule.
1483 
1484  Output, int SUBORDER_XYZ_N[4*SUBORDER_NUM],
1485  the numerators of the barycentric coordinates of the abscissas.
1486 
1487  Output, int *SUBORDER_XYZ_D,
1488  the denominator of the barycentric coordinates of the abscissas.
1489 
1490  Output, int SUBORDER_W_N[SUBORDER_NUM],
1491  the numerator of the suborder weights.
1492 
1493  Output, int SUBORDER_W_D,
1494  the denominator of the suborder weights.
1495 */
1496 {
1497  int i;
1498  int s;
1499  int suborder_xyz_n_07[4*9] = {
1500  0, 0, 0, 6,
1501  0, 0, 5, 1,
1502  0, 0, 4, 2,
1503  1, 1, 0, 4,
1504  3, 3, 0, 0,
1505  3, 2, 1, 0,
1506  1, 1, 1, 3,
1507  2, 2, 2, 0,
1508  2, 2, 1, 1
1509  };
1510  int suborder_xyz_d_07 = 6;
1511  int suborder_w_n_07[9] = { -7, 24, -30, 0, 40, 30, 180, -45, 0 };
1512  int suborder_w_d_07 = 1400;
1513 
1514  for ( s = 0; s < suborder_num; s++ )
1515  {
1516  for ( i = 0; i < 4; i++ )
1517  {
1518  suborder_xyz_n[i+s*4] = suborder_xyz_n_07[i+s*4];
1519  }
1520  }
1521  *suborder_xyz_d = suborder_xyz_d_07;
1522 
1523  for ( s = 0; s < suborder_num; s++ )
1524  {
1525  suborder_w_n[s] = suborder_w_n_07[s];
1526  }
1527  *suborder_w_d = suborder_w_d_07;
1528 
1529  return;
1530 }
1531 /******************************************************************************/
1532 
1533 double tetrahedron_volume ( double tetra[3*4] )
1534 
1535 /******************************************************************************/
1536 /*
1537  Purpose:
1538 
1539  TETRAHEDRON_VOLUME computes the volume of a tetrahedron.
1540 
1541  Licensing:
1542 
1543  This code is distributed under the GNU LGPL license.
1544 
1545  Modified:
1546 
1547  06 August 2005
1548 
1549  Author:
1550 
1551  John Burkardt
1552 
1553  Parameters:
1554 
1555  Input, double TETRA[3*4], the vertices of the tetrahedron.
1556 
1557  Output, double TETRAHEDRON_VOLUME, the volume of the tetrahedron.
1558 */
1559 {
1560  double a[4*4];
1561  int i;
1562  int j;
1563  double volume;
1564 
1565  for ( i = 0; i < 3; i++ )
1566  {
1567  for ( j = 0; j < 4; j++ )
1568  {
1569  a[i+j*4] = tetra[i+j*3];
1570  }
1571  }
1572 
1573  i = 3;
1574  for ( j = 0; j < 4; j++ )
1575  {
1576  a[i+j*4] = 1.0;
1577  }
1578 
1579  volume = fabs ( r8mat_det_4d ( a ) ) / 6.0;
1580 
1581  return volume;
1582 }
1583 /******************************************************************************/
1584 
1585 // void timestamp ( )
1586 
1587 // /******************************************************************************/
1588 // /*
1589 // Purpose:
1590 
1591 // TIMESTAMP prints the current YMDHMS date as a time stamp.
1592 
1593 // Example:
1594 
1595 // 31 May 2001 09:45:54 AM
1596 
1597 // Licensing:
1598 
1599 // This code is distributed under the GNU LGPL license.
1600 
1601 // Modified:
1602 
1603 // 24 September 2003
1604 
1605 // Author:
1606 
1607 // John Burkardt
1608 
1609 // Parameters:
1610 
1611 // None
1612 // */
1613 // {
1614 // # define TIME_SIZE 40
1615 
1616 // static char time_buffer[TIME_SIZE];
1617 // const struct tm *tm;
1618 // time_t now;
1619 
1620 // now = time ( NULL );
1621 // tm = localtime ( &now );
1622 
1623 // strftime ( time_buffer, TIME_SIZE, "%d %B %Y %I:%M:%S %p", tm );
1624 
1625 // fprintf ( stdout, "%s\n", time_buffer );
1626 
1627 // return;
1628 // # undef TIME_SIZE
1629 // }
1630 
reference_to_physical_t4
void reference_to_physical_t4(double t[], int n, double ref[], double phy[])
Definition: tetrahedron_ncc_rule.c:61
tetrahedron_ncc_suborder
int * tetrahedron_ncc_suborder(int rule, int suborder_num)
Definition: tetrahedron_ncc_rule.c:709
PlasticOps::w
double w(double eqiv, double dot_tau, double f, double sigma_y, double sigma_Y)
Definition: PlasticOpsGeneric.hpp:276
tetrahedron_ncc_subrule_01
void tetrahedron_ncc_subrule_01(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
Definition: tetrahedron_ncc_rule.c:1002
tetrahedron_ncc_subrule_07
void tetrahedron_ncc_subrule_07(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
Definition: tetrahedron_ncc_rule.c:1452
tetrahedron_ncc_rule_num
int tetrahedron_ncc_rule_num()
Definition: tetrahedron_ncc_rule.c:669
order
constexpr int order
Definition: dg_projection.cpp:18
tetrahedron_ncc_subrule_03
void tetrahedron_ncc_subrule_03(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
Definition: tetrahedron_ncc_rule.c:1148
tetrahedron_ncc_rule.h
r8mat_det_4d
double r8mat_det_4d(double a[])
Definition: tetrahedron_ncc_rule.c:5
tetrahedron_ncc_subrule_04
void tetrahedron_ncc_subrule_04(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
Definition: tetrahedron_ncc_rule.c:1222
tetrahedron_volume
double tetrahedron_volume(double tetra[3 *4])
Definition: tetrahedron_ncc_rule.c:1533
tetrahedron_ncc_suborder_num
int tetrahedron_ncc_suborder_num(int rule)
Definition: tetrahedron_ncc_rule.c:810
tetrahedron_ncc_subrule
void tetrahedron_ncc_subrule(int rule, int suborder_num, double suborder_xyz[], double suborder_w[])
Definition: tetrahedron_ncc_rule.c:888
a
constexpr double a
Definition: approx_sphere.cpp:30
double
tetrahedron_ncc_degree
int tetrahedron_ncc_degree(int rule)
Definition: tetrahedron_ncc_rule.c:126
t
constexpr double t
plate stiffness
Definition: plate.cpp:59
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
convert.n
n
Definition: convert.py:82
tetrahedron_ncc_subrule_05
void tetrahedron_ncc_subrule_05(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
Definition: tetrahedron_ncc_rule.c:1297
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
tetrahedron_ncc_subrule_02
void tetrahedron_ncc_subrule_02(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
Definition: tetrahedron_ncc_rule.c:1075
tetrahedron_ncc_subrule_06
void tetrahedron_ncc_subrule_06(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
Definition: tetrahedron_ncc_rule.c:1374
tetrahedron_ncc_order_num
int tetrahedron_ncc_order_num(int rule)
Definition: tetrahedron_ncc_rule.c:180
tetrahedron_ncc_rule
void tetrahedron_ncc_rule(int rule, int order_num, double xyz[], double w[])
Definition: tetrahedron_ncc_rule.c:236
Generated by Doxygen 1.8.17 and hosted at University of Glasgow