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