v0.9.0
triangle_nco_rule.c
Go to the documentation of this file.
1 # include <stdlib.h>
2 # include <stdio.h>
3 # include <math.h>
4 # include <time.h>
5 
6 # include "triangle_nco_rule.h"
7 
8 /******************************************************************************/
9 
10 //int i4_max ( int i1, int i2 ) //NOTE: it is defined in gm_rule.c
11 
12 /******************************************************************************/
13 /*
14  Purpose:
15 
16  I4_MAX returns the maximum of two I4's.
17 
18  Licensing:
19 
20  This code is distributed under the GNU LGPL license.
21 
22  Modified:
23 
24  29 August 2006
25 
26  Author:
27 
28  John Burkardt
29 
30  Parameters:
31 
32  Input, int I1, I2, are two integers to be compared.
33 
34  Output, int I4_MAX, the larger of I1 and I2.
35 */
36 /*{
37  int value;
38 
39  if ( i2 < i1 )
40  {
41  value = i1;
42  }
43  else
44  {
45  value = i2;
46  }
47  return value;
48 }*/
49 /******************************************************************************/
50 
51 //int i4_min ( int i1, int i2 ) //NOTE: it is defined in gm_rule.c
52 
53 /******************************************************************************/
54 /*
55  Purpose:
56 
57  I4_MIN returns the smaller of two I4's.
58 
59  Licensing:
60 
61  This code is distributed under the GNU LGPL license.
62 
63  Modified:
64 
65  29 August 2006
66 
67  Author:
68 
69  John Burkardt
70 
71  Parameters:
72 
73  Input, int I1, I2, two integers to be compared.
74 
75  Output, int I4_MIN, the smaller of I1 and I2.
76 */
77 /*{
78  int value;
79 
80  if ( i1 < i2 )
81  {
82  value = i1;
83  }
84  else
85  {
86  value = i2;
87  }
88  return value;
89 }*/
90 /*******************************************************************************/
91 
92 int i4_modp ( int i, int j )
93 
94 /******************************************************************************/
95 /*
96  Purpose:
97 
98  I4_MODP returns the nonnegative remainder of I4 division.
99 
100  Discussion:
101 
102  If
103  NREM = I4_MODP ( I, J )
104  NMULT = ( I - NREM ) / J
105  then
106  I = J * NMULT + NREM
107  where NREM is always nonnegative.
108 
109  The MOD function computes a result with the same sign as the
110  quantity being divided. Thus, suppose you had an angle A,
111  and you wanted to ensure that it was between 0 and 360.
112  Then mod(A,360) would do, if A was positive, but if A
113  was negative, your result would be between -360 and 0.
114 
115  On the other hand, I4_MODP(A,360) is between 0 and 360, always.
116 
117  Example:
118 
119  I J MOD I4_MODP I4_MODP Factorization
120 
121  107 50 7 7 107 = 2 * 50 + 7
122  107 -50 7 7 107 = -2 * -50 + 7
123  -107 50 -7 43 -107 = -3 * 50 + 43
124  -107 -50 -7 43 -107 = 3 * -50 + 43
125 
126  Licensing:
127 
128  This code is distributed under the GNU LGPL license.
129 
130  Modified:
131 
132  12 January 2007
133 
134  Author:
135 
136  John Burkardt
137 
138  Parameters:
139 
140  Input, int I, the number to be divided.
141 
142  Input, int J, the number that divides I.
143 
144  Output, int I4_MODP, the nonnegative remainder when I is
145  divided by J.
146 */
147 {
148  int value;
149 
150  if ( j == 0 )
151  {
152  fprintf ( stderr, "\n" );
153  fprintf ( stderr, "I4_MODP - Fatal error!\n" );
154  fprintf ( stderr, " I4_MODP ( I, J ) called with J = %d\n", j );
155  exit ( 1 );
156  }
157 
158  value = i % j;
159 
160  if ( value < 0 )
161  {
162  value = value + abs ( j );
163  }
164 
165  return value;
166 }
167 /******************************************************************************/
168 
169 int i4_wrap ( int ival, int ilo, int ihi )
170 
171 /******************************************************************************/
172 /*
173  Purpose:
174 
175  I4_WRAP forces an I4 to lie between given limits by wrapping.
176 
177  Example:
178 
179  ILO = 4, IHI = 8
180 
181  I Value
182 
183  -2 8
184  -1 4
185  0 5
186  1 6
187  2 7
188  3 8
189  4 4
190  5 5
191  6 6
192  7 7
193  8 8
194  9 4
195  10 5
196  11 6
197  12 7
198  13 8
199  14 4
200 
201  Licensing:
202 
203  This code is distributed under the GNU LGPL license.
204 
205  Modified:
206 
207  26 December 2012
208 
209  Author:
210 
211  John Burkardt
212 
213  Parameters:
214 
215  Input, int IVAL, an integer value.
216 
217  Input, int ILO, IHI, the desired bounds for the integer value.
218 
219  Output, int I4_WRAP, a "wrapped" version of IVAL.
220 */
221 {
222  int jhi;
223  int jlo;
224  int value;
225  int wide;
226 
227  if ( ilo < ihi )
228  {
229  jlo = ilo;
230  jhi = ihi;
231  }
232  else
233  {
234  jlo = ihi;
235  jhi = ilo;
236  }
237 
238  wide = jhi + 1 - jlo;
239 
240  if ( wide == 1 )
241  {
242  value = jlo;
243  }
244  else
245  {
246  value = jlo + i4_modp ( ival - jlo, wide );
247  }
248 
249  return value;
250 }
251 /******************************************************************************/
252 
253 double r8_huge ( )
254 
255 /******************************************************************************/
256 /*
257  Purpose:
258 
259  R8_HUGE returns a "huge" R8.
260 
261  Discussion:
262 
263  The value returned by this function is NOT required to be the
264  maximum representable R8. This value varies from machine to machine,
265  from compiler to compiler, and may cause problems when being printed.
266  We simply want a "very large" but non-infinite number.
267 
268  Licensing:
269 
270  This code is distributed under the GNU LGPL license.
271 
272  Modified:
273 
274  06 October 2007
275 
276  Author:
277 
278  John Burkardt
279 
280  Parameters:
281 
282  Output, double R8_HUGE, a "huge" R8 value.
283 */
284 {
285  double value;
286 
287  value = 1.0E+30;
288 
289  return value;
290 }
291 /******************************************************************************/
292 
293 int r8_nint ( double x )
294 
295 /******************************************************************************/
296 /*
297  Purpose:
298 
299  R8_NINT returns the nearest integer to an R8.
300 
301  Example:
302 
303  X R8_NINT
304 
305  1.3 1
306  1.4 1
307  1.5 1 or 2
308  1.6 2
309  0.0 0
310  -0.7 -1
311  -1.1 -1
312  -1.6 -2
313 
314  Licensing:
315 
316  This code is distributed under the GNU LGPL license.
317 
318  Modified:
319 
320  05 May 2006
321 
322  Author:
323 
324  John Burkardt
325 
326  Parameters:
327 
328  Input, double X, the value.
329 
330  Output, int R8_NINT, the nearest integer to X.
331 */
332 {
333  int s;
334  int value;
335 
336  if ( x < 0.0 )
337  {
338  s = - 1;
339  }
340  else
341  {
342  s = + 1;
343  }
344  value = s * ( int ) ( fabs ( x ) + 0.5 );
345 
346  return value;
347 }
348 /******************************************************************************/
349 
350 void reference_to_physical_t3 ( double t[], int n, double ref[], double phy[] )
351 
352 /******************************************************************************/
353 /*
354  Purpose:
355 
356  REFERENCE_TO_PHYSICAL_T3 maps T3 reference points to physical points.
357 
358  Discussion:
359 
360  Given the vertices of an order 3 physical triangle and a point
361  (XSI,ETA) in the reference triangle, the routine computes the value
362  of the corresponding image point (X,Y) in physical space.
363 
364  Note that this routine may also be appropriate for an order 6
365  triangle, if the mapping between reference and physical space
366  is linear. This implies, in particular, that the sides of the
367  image triangle are straight and that the "midside" nodes in the
368  physical triangle are literally halfway along the sides of
369  the physical triangle.
370 
371  Reference Element T3:
372 
373  |
374  1 3
375  | |\
376  | | \
377  S | \
378  | | \
379  | | \
380  0 1-----2
381  |
382  +--0--R--1-->
383 
384  Licensing:
385 
386  This code is distributed under the GNU LGPL license.
387 
388  Modified:
389 
390  24 June 2005
391 
392  Author:
393 
394  John Burkardt
395 
396  Parameters:
397 
398  Input, double T[2*3], the coordinates of the vertices.
399  The vertices are assumed to be the images of (0,0), (1,0) and
400  (0,1) respectively.
401 
402  Input, int N, the number of objects to transform.
403 
404  Input, double REF[2*N], points in the reference triangle.
405 
406  Output, double PHY[2*N], corresponding points in the
407  physical triangle.
408 */
409 {
410  int i;
411  int j;
412 
413  for ( i = 0; i < 2; i++ )
414  {
415  for ( j = 0; j < n; j++ )
416  {
417  phy[i+j*2] = t[i+0*2] * ( 1.0 - ref[0+j*2] - ref[1+j*2] )
418  + t[i+1*2] * + ref[0+j*2]
419  + t[i+2*2] * + ref[1+j*2];
420  }
421  }
422 
423  return;
424 }
425 /******************************************************************************/
426 
427 //void timestamp ( )
428 
429 /******************************************************************************/
430 /*
431  Purpose:
432 
433  TIMESTAMP prints the current YMDHMS date as a time stamp.
434 
435  Example:
436 
437  31 May 2001 09:45:54 AM
438 
439  Licensing:
440 
441  This code is distributed under the GNU LGPL license.
442 
443  Modified:
444 
445  24 September 2003
446 
447  Author:
448 
449  John Burkardt
450 
451  Parameters:
452 
453  None
454 */
455 /*{
456 # define TIME_SIZE 40
457 
458  static char time_buffer[TIME_SIZE];
459  const struct tm *tm;
460  size_t len;
461  time_t now;
462 
463  now = time ( NULL );
464  tm = localtime ( &now );
465 
466  len = strftime ( time_buffer, TIME_SIZE, "%d %B %Y %I:%M:%S %p", tm );
467 
468  fprintf ( stdout, "%s\n", time_buffer );
469 
470  return;
471 # undef TIME_SIZE
472 }*/
473 /******************************************************************************/
474 
475 double triangle_area ( double t[2*3] )
476 
477 /*****************************************************************************80
478 
479  Purpose:
480 
481  TRIANGLE_AREA computes the area of a triangle.
482 
483  Discussion:
484 
485  If the triangle's vertices are given in counter clockwise order,
486  the area will be positive. If the triangle's vertices are given
487  in clockwise order, the area will be negative!
488 
489  An earlier version of this routine always returned the absolute
490  value of the computed area. I am convinced now that that is
491  a less useful result! For instance, by returning the signed
492  area of a triangle, it is possible to easily compute the area
493  of a nonconvex polygon as the sum of the (possibly negative)
494  areas of triangles formed by node 1 and successive pairs of vertices.
495 
496  Licensing:
497 
498  This code is distributed under the GNU LGPL license.
499 
500  Modified:
501 
502  17 October 2005
503 
504  Author:
505 
506  John Burkardt
507 
508  Parameters:
509 
510  Input, double T[2*3], the vertices of the triangle.
511 
512  Output, double TRIANGLE_AREA, the area of the triangle.
513 */
514 {
515  double area;
516 
517  area = 0.5 * (
518  t[0+0*2] * ( t[1+1*2] - t[1+2*2] ) +
519  t[0+1*2] * ( t[1+2*2] - t[1+0*2] ) +
520  t[0+2*2] * ( t[1+0*2] - t[1+1*2] ) );
521 
522  return area;
523 }
524 /******************************************************************************/
525 
526 int triangle_nco_degree ( int rule )
527 
528 /******************************************************************************/
529 /*
530  Purpose:
531 
532  TRIANGLE_NCO_DEGREE returns the degree of an NCO rule for the triangle.
533 
534  Licensing:
535 
536  This code is distributed under the GNU LGPL license.
537 
538  Modified:
539 
540  30 January 2007
541 
542  Author:
543 
544  John Burkardt
545 
546  Reference:
547 
548  Peter Silvester,
549  Symmetric Quadrature Formulae for Simplexes,
550  Mathematics of Computation,
551  Volume 24, Number 109, January 1970, pages 95-100.
552 
553  Parameters:
554 
555  Input, int RULE, the index of the rule.
556 
557  Output, int TRIANGLE_NCO_DEGREE, the polynomial degree of exactness of
558  the rule.
559 */
560 {
561  int degree;
562 
563  if ( 1 <= rule && rule <= 9 )
564  {
565  degree = rule - 1;
566  }
567  else
568  {
569  degree = -1;
570  fprintf ( stderr, "\n" );
571  fprintf ( stderr, "TRIANGLE_NCO_DEGREE - Fatal error!\n" );
572  fprintf ( stderr, " Illegal RULE = %d\n", rule );
573  exit ( 1 );
574  }
575 
576  return degree;
577 }
578 /******************************************************************************/
579 
580 int triangle_nco_order_num ( int rule )
581 
582 /******************************************************************************/
583 /*
584  Purpose:
585 
586  TRIANGLE_NCO_ORDER_NUM returns the order of an NCO rule for the triangle.
587 
588  Licensing:
589 
590  This code is distributed under the GNU LGPL license.
591 
592  Modified:
593 
594  30 January 2007
595 
596  Author:
597 
598  John Burkardt
599 
600  Reference:
601 
602  Peter Silvester,
603  Symmetric Quadrature Formulae for Simplexes,
604  Mathematics of Computation,
605  Volume 24, Number 109, January 1970, pages 95-100.
606 
607  Parameters:
608 
609  Input, int RULE, the index of the rule.
610 
611  Output, int TRIANGLE_NCO_ORDER_NUM, the order (number of points)
612  of the rule.
613 */
614 {
615  int order;
616  int order_num;
617  int *suborder;
618  int suborder_num;
619 
620  suborder_num = triangle_nco_suborder_num ( rule );
621 
622  suborder = triangle_nco_suborder ( rule, suborder_num );
623 
624  order_num = 0;
625  for ( order = 0; order < suborder_num; order++ )
626  {
627  order_num = order_num + suborder[order];
628  }
629 
630  free ( suborder );
631 
632  return order_num;
633 }
634 /******************************************************************************/
635 
636 void triangle_nco_rule ( int rule, int order_num, double xy[], double w[] )
637 
638 /******************************************************************************/
639 /*
640  Purpose:
641 
642  TRIANGLE_NCO_RULE returns the points and weights of an NCO rule.
643 
644  Licensing:
645 
646  This code is distributed under the GNU LGPL license.
647 
648  Modified:
649 
650  30 January 2007
651 
652  Author:
653 
654  John Burkardt
655 
656  Reference:
657 
658  Peter Silvester,
659  Symmetric Quadrature Formulae for Simplexes,
660  Mathematics of Computation,
661  Volume 24, Number 109, January 1970, pages 95-100.
662 
663  Parameters:
664 
665  Input, int RULE, the index of the rule.
666 
667  Input, int ORDER_NUM, the order (number of points) of the rule.
668 
669  Output, double XY[2*ORDER_NUM], the points of the rule.
670 
671  Output, double W[ORDER_NUM], the weights of the rule.
672 */
673 {
674  int k;
675  int o;
676  int s;
677  int *suborder;
678  int suborder_num;
679  double *suborder_w;
680  double *suborder_xyz;
681 /*
682  Get the suborder information.
683 */
684  suborder_num = triangle_nco_suborder_num ( rule );
685 
686  suborder_xyz = ( double * ) malloc ( 3 * suborder_num * sizeof ( double ) );
687  suborder_w = ( double * ) malloc ( suborder_num * sizeof ( double ) );
688 
689  suborder = triangle_nco_suborder ( rule, suborder_num );
690 
691  triangle_nco_subrule ( rule, suborder_num, suborder_xyz, suborder_w );
692 /*
693  Expand the suborder information to a full order rule.
694 */
695  o = 0;
696 
697  for ( s = 0; s < suborder_num; s++ )
698  {
699  if ( suborder[s] == 1 )
700  {
701  xy[0+o*2] = suborder_xyz[0+s*3];
702  xy[1+o*2] = suborder_xyz[1+s*3];
703  w[o] = suborder_w[s];
704  o = o + 1;
705  }
706  else if ( suborder[s] == 3 )
707  {
708  for ( k = 0; k < 3; k++ )
709  {
710  xy[0+o*2] = suborder_xyz [ i4_wrap(k, 0,2) + s*3 ];
711  xy[1+o*2] = suborder_xyz [ i4_wrap(k+1,0,2) + s*3 ];
712  w[o] = suborder_w[s];
713  o = o + 1;
714  }
715  }
716  else if ( suborder[s] == 6 )
717  {
718  for ( k = 0; k < 3; k++ )
719  {
720  xy[0+o*2] = suborder_xyz [ i4_wrap(k, 0,2) + s*3 ];
721  xy[1+o*2] = suborder_xyz [ i4_wrap(k+1,0,2) + s*3 ];
722  w[o] = suborder_w[s];
723  o = o + 1;
724  }
725 
726  for ( k = 0; k < 3; k++ )
727  {
728  xy[0+o*2] = suborder_xyz [ i4_wrap(k+1,0,2) + s*3 ];
729  xy[1+o*2] = suborder_xyz [ i4_wrap(k, 0,2) + s*3 ];
730  w[o] = suborder_w[s];
731  o = o + 1;
732  }
733  }
734  else
735  {
736  fprintf ( stderr, "\n" );
737  fprintf ( stderr, "TRIANGLE_NCO_RULE - Fatal error!\n" );
738  fprintf ( stderr, " Illegal SUBORDER(%d) = %d\n", s, suborder[s] );
739  exit ( 1 );
740  }
741  }
742 
743  free ( suborder );
744  free ( suborder_xyz );
745  free ( suborder_w );
746 
747  return;
748 }
749 /******************************************************************************/
750 
752 
753 /******************************************************************************/
754 /*
755  Purpose:
756 
757  TRIANGLE_NCO_RULE_NUM returns the number of NCO rules available.
758 
759  Licensing:
760 
761  This code is distributed under the GNU LGPL license.
762 
763  Modified:
764 
765  30 January 2007
766 
767  Author:
768 
769  John Burkardt
770 
771  Reference:
772 
773  Peter Silvester,
774  Symmetric Quadrature Formulae for Simplexes,
775  Mathematics of Computation,
776  Volume 24, Number 109, January 1970, pages 95-100.
777 
778  Parameters:
779 
780  Output, int TRIANGLE_NCO_RULE_NUM, the number of rules available.
781 */
782 {
783  int rule_num;
784 
785  rule_num = 9;
786 
787  return rule_num;
788 }
789 /******************************************************************************/
790 
791 int *triangle_nco_suborder ( int rule, int suborder_num )
792 
793 /******************************************************************************/
794 /*
795  Purpose:
796 
797  TRIANGLE_NCO_SUBORDER returns the suborders for an NCO rule.
798 
799  Licensing:
800 
801  This code is distributed under the GNU LGPL license.
802 
803  Modified:
804 
805  30 January 2007
806 
807  Author:
808 
809  John Burkardt
810 
811  Reference:
812 
813  Peter Silvester,
814  Symmetric Quadrature Formulae for Simplexes,
815  Mathematics of Computation,
816  Volume 24, Number 109, January 1970, pages 95-100.
817 
818  Parameters:
819 
820  Input, int RULE, the index of the rule.
821 
822  Input, int SUBORDER_NUM, the number of suborders of the rule.
823 
824  Output, int TRIANGLE_NCO_SUBORDER[SUBORDER_NUM],
825  the suborders of the rule.
826 */
827 {
828  int *suborder;
829 
830  suborder = ( int * ) malloc ( suborder_num * sizeof ( int ) );
831 
832  if ( rule == 1 )
833  {
834  suborder[0] = 1;
835  }
836  else if ( rule == 2 )
837  {
838  suborder[0] = 3;
839  }
840  else if ( rule == 3 )
841  {
842  suborder[0] = 3;
843  suborder[1] = 3;
844  }
845  else if ( rule == 4 )
846  {
847  suborder[0] = 3;
848  suborder[1] = 6;
849  suborder[2] = 1;
850  }
851  else if ( rule == 5 )
852  {
853  suborder[0] = 3;
854  suborder[1] = 6;
855  suborder[2] = 3;
856  suborder[3] = 3;
857  }
858  else if ( rule == 6 )
859  {
860  suborder[0] = 3;
861  suborder[1] = 6;
862  suborder[2] = 6;
863  suborder[3] = 3;
864  suborder[4] = 3;
865  }
866  else if ( rule == 7 )
867  {
868  suborder[0] = 3;
869  suborder[1] = 6;
870  suborder[2] = 6;
871  suborder[3] = 3;
872  suborder[4] = 3;
873  suborder[5] = 6;
874  suborder[6] = 1;
875  }
876  else if ( rule == 8 )
877  {
878  suborder[0] = 3;
879  suborder[1] = 6;
880  suborder[2] = 6;
881  suborder[3] = 3;
882  suborder[4] = 6;
883  suborder[5] = 6;
884  suborder[6] = 3;
885  suborder[7] = 3;
886  }
887  else if ( rule == 9 )
888  {
889  suborder[0] = 3;
890  suborder[1] = 6;
891  suborder[2] = 6;
892  suborder[3] = 3;
893  suborder[4] = 6;
894  suborder[5] = 6;
895  suborder[6] = 3;
896  suborder[7] = 6;
897  suborder[8] = 3;
898  suborder[9] = 3;
899  }
900  else
901  {
902  fprintf ( stderr, "\n" );
903  fprintf ( stderr, "TRIANGLE_NCO_SUBORDER - Fatal error!\n" );
904  fprintf ( stderr, " Illegal RULE = %d\n", rule );
905  exit ( 1 );
906  }
907 
908  return suborder;
909 }
910 /******************************************************************************/
911 
913 
914 /******************************************************************************/
915 /*
916  Purpose:
917 
918  TRIANGLE_NCO_SUBORDER_NUM returns the number of suborders for an NCO rule.
919 
920  Licensing:
921 
922  This code is distributed under the GNU LGPL license.
923 
924  Modified:
925 
926  30 January 2007
927 
928  Author:
929 
930  John Burkardt
931 
932  Reference:
933 
934  Peter Silvester,
935  Symmetric Quadrature Formulae for Simplexes,
936  Mathematics of Computation,
937  Volume 24, Number 109, January 1970, pages 95-100.
938 
939  Parameters:
940 
941  Input, int RULE, the index of the rule.
942 
943  Output, int TRIANGLE_NCO_SUBORDER_NUM, the number of suborders
944  of the rule.
945 */
946 {
947  int suborder_num;
948 
949  if ( rule == 1 )
950  {
951  suborder_num = 1;
952  }
953  else if ( rule == 2 )
954  {
955  suborder_num = 1;
956  }
957  else if ( rule == 3 )
958  {
959  suborder_num = 2;
960  }
961  else if ( rule == 4 )
962  {
963  suborder_num = 3;
964  }
965  else if ( rule == 5 )
966  {
967  suborder_num = 4;
968  }
969  else if ( rule == 6 )
970  {
971  suborder_num = 5;
972  }
973  else if ( rule == 7 )
974  {
975  suborder_num = 7;
976  }
977  else if ( rule == 8 )
978  {
979  suborder_num = 8;
980  }
981  else if ( rule == 9 )
982  {
983  suborder_num = 10;
984  }
985  else
986  {
987  suborder_num = -1;
988  fprintf ( stderr, "\n" );
989  fprintf ( stderr, "TRIANGLE_NCO_SUBORDER_NUM - Fatal error!\n" );
990  fprintf ( stderr," Illegal RULE = %d\n", rule );
991  exit ( 1 );
992  }
993 
994  return suborder_num;
995 }
996 /******************************************************************************/
997 
998 void triangle_nco_subrule ( int rule, int suborder_num, double suborder_xyz[],
999  double suborder_w[] )
1000 
1001 /******************************************************************************/
1002 /*
1003  Purpose:
1004 
1005  TRIANGLE_NCO_SUBRULE returns a compressed NCO rule.
1006 
1007  Licensing:
1008 
1009  This code is distributed under the GNU LGPL license.
1010 
1011  Modified:
1012 
1013  30 January 2007
1014 
1015  Author:
1016 
1017  John Burkardt
1018 
1019  Reference:
1020 
1021  Peter Silvester,
1022  Symmetric Quadrature Formulae for Simplexes,
1023  Mathematics of Computation,
1024  Volume 24, Number 109, January 1970, pages 95-100.
1025 
1026  Parameters:
1027 
1028  Input, int RULE, the index of the rule.
1029 
1030  Input, int SUBORDER_NUM, the number of suborders of the rule.
1031 
1032  Output, double SUBORDER_XYZ[3*SUBORDER_NUM],
1033  the barycentric coordinates of the abscissas.
1034 
1035  Output, double SUBORDER_W[SUBORDER_NUM], the suborder weights.
1036 */
1037 {
1038  int i;
1039  int s;
1040  int suborder_w_d;
1041  int *suborder_w_n;
1042  int suborder_xyz_d;
1043  int *suborder_xyz_n;
1044 
1045  suborder_xyz_n = ( int * ) malloc ( 3 * suborder_num * sizeof ( int ) );
1046  suborder_w_n = ( int * ) malloc ( suborder_num * sizeof ( int ) );
1047 
1048  if ( rule == 1 )
1049  {
1050  triangle_nco_subrule_01 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
1051  suborder_w_n, &suborder_w_d );
1052  }
1053  else if ( rule == 2 )
1054  {
1055  triangle_nco_subrule_02 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
1056  suborder_w_n, &suborder_w_d );
1057  }
1058  else if ( rule == 3 )
1059  {
1060  triangle_nco_subrule_03 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
1061  suborder_w_n, &suborder_w_d );
1062  }
1063  else if ( rule == 4 )
1064  {
1065  triangle_nco_subrule_04 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
1066  suborder_w_n, &suborder_w_d );
1067  }
1068  else if ( rule == 5 )
1069  {
1070  triangle_nco_subrule_05 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
1071  suborder_w_n, &suborder_w_d );
1072  }
1073  else if ( rule == 6 )
1074  {
1075  triangle_nco_subrule_06 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
1076  suborder_w_n, &suborder_w_d );
1077  }
1078  else if ( rule == 7 )
1079  {
1080  triangle_nco_subrule_07 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
1081  suborder_w_n, &suborder_w_d );
1082  }
1083  else if ( rule == 8 )
1084  {
1085  triangle_nco_subrule_08 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
1086  suborder_w_n, &suborder_w_d );
1087  }
1088  else if ( rule == 9 )
1089  {
1090  triangle_nco_subrule_09 ( suborder_num, suborder_xyz_n, &suborder_xyz_d,
1091  suborder_w_n, &suborder_w_d );
1092  }
1093  else
1094  {
1095  fprintf ( stderr, "\n" );
1096  fprintf ( stderr, "TRIANGLE_NCO_SUBRULE - Fatal error!\n" );
1097  fprintf ( stderr, " Illegal RULE = %d\n", rule );
1098  exit ( 1 );
1099  }
1100 
1101  for ( s = 0; s < suborder_num; s++ )
1102  {
1103  for ( i = 0; i < 3; i++ )
1104  {
1105  suborder_xyz[i+s*3] =
1106  ( double ) ( 1 + suborder_xyz_n[i+s*3] )
1107  / ( double ) ( 3 + suborder_xyz_d );
1108  }
1109  }
1110  for ( s = 0; s < suborder_num; s++ )
1111  {
1112  suborder_w[s] = ( double ) suborder_w_n[s] / ( double ) suborder_w_d;
1113  }
1114 
1115  free ( suborder_w_n );
1116  free ( suborder_xyz_n );
1117 
1118  return;
1119 }
1120 /******************************************************************************/
1121 
1122 void triangle_nco_subrule_01 ( int suborder_num, int suborder_xyz_n[],
1123  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1124 
1125 /******************************************************************************/
1126 /*
1127  Purpose:
1128 
1129  TRIANGLE_NCO_SUBRULE_01 returns a compressed NCO rule 1.
1130 
1131  Licensing:
1132 
1133  This code is distributed under the GNU LGPL license.
1134 
1135  Modified:
1136 
1137  30 January 2007
1138 
1139  Author:
1140 
1141  John Burkardt
1142 
1143  Reference:
1144 
1145  Peter Silvester,
1146  Symmetric Quadrature Formulae for Simplexes,
1147  Mathematics of Computation,
1148  Volume 24, Number 109, January 1970, pages 95-100.
1149 
1150  Parameters:
1151 
1152  Input, int SUBORDER_NUM, the number of suborders of the rule.
1153 
1154  Output, int SUBORDER_XYZ_N[3*SUBORDER_NUM],
1155  the numerators of the barycentric coordinates of the abscissas.
1156 
1157  Output, int *SUBORDER_XYZ_D,
1158  the denominator of the barycentric coordinates of the abscissas.
1159 
1160  Output, int SUBORDER_W_N[SUBORDER_NUM],
1161  the numerator of the suborder weights.
1162 
1163  Output, int SUBORDER_W_D,
1164  the denominator of the suborder weights.
1165 */
1166 {
1167  int i;
1168  int s;
1169  int suborder_xyz_n_01[3*1] = {
1170  0, 0, 0
1171  };
1172  int suborder_xyz_d_01 = 0;
1173  int suborder_w_n_01[1] = { 1 };
1174  int suborder_w_d_01 = 1;
1175 
1176  for ( s = 0; s < suborder_num; s++ )
1177  {
1178  for ( i = 0; i < 3; i++ )
1179  {
1180  suborder_xyz_n[i+s*3] = suborder_xyz_n_01[i+s*3];
1181  }
1182  }
1183  *suborder_xyz_d = suborder_xyz_d_01;
1184 
1185  for ( s = 0; s < suborder_num; s++ )
1186  {
1187  suborder_w_n[s] = suborder_w_n_01[s];
1188  }
1189  *suborder_w_d = suborder_w_d_01;
1190 
1191  return;
1192 }
1193 /******************************************************************************/
1194 
1195 void triangle_nco_subrule_02 ( int suborder_num, int suborder_xyz_n[],
1196  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1197 
1198 /******************************************************************************/
1199 /*
1200  Purpose:
1201 
1202  TRIANGLE_NCO_SUBRULE_02 returns a compressed NCO rule 2.
1203 
1204  Licensing:
1205 
1206  This code is distributed under the GNU LGPL license.
1207 
1208  Modified:
1209 
1210  30 January 2007
1211 
1212  Author:
1213 
1214  John Burkardt
1215 
1216  Reference:
1217 
1218  Peter Silvester,
1219  Symmetric Quadrature Formulae for Simplexes,
1220  Mathematics of Computation,
1221  Volume 24, Number 109, January 1970, pages 95-100.
1222 
1223  Parameters:
1224 
1225  Input, int SUBORDER_NUM, the number of suborders of the rule.
1226 
1227  Output, int SUBORDER_XYZ_N[3*SUBORDER_NUM],
1228  the numerators of the barycentric coordinates of the abscissas.
1229 
1230  Output, int *SUBORDER_XYZ_D,
1231  the denominator of the barycentric coordinates of the abscissas.
1232 
1233  Output, int SUBORDER_W_N[SUBORDER_NUM],
1234  the numerator of the suborder weights.
1235 
1236  Output, int SUBORDER_W_D,
1237  the denominator of the suborder weights.
1238 */
1239 {
1240  int i;
1241  int s;
1242  int suborder_xyz_n_02[3*1] = {
1243  1, 0, 0
1244  };
1245  int suborder_xyz_d_02 = 1;
1246  int suborder_w_n_02[1] = { 1 };
1247  int suborder_w_d_02 = 3;
1248 
1249  for ( s = 0; s < suborder_num; s++ )
1250  {
1251  for ( i = 0; i < 3; i++ )
1252  {
1253  suborder_xyz_n[i+s*3] = suborder_xyz_n_02[i+s*3];
1254  }
1255  }
1256  *suborder_xyz_d = suborder_xyz_d_02;
1257 
1258  for ( s = 0; s < suborder_num; s++ )
1259  {
1260  suborder_w_n[s] = suborder_w_n_02[s];
1261  }
1262  *suborder_w_d = suborder_w_d_02;
1263 
1264  return;
1265 }
1266 /******************************************************************************/
1267 
1268 void triangle_nco_subrule_03 ( int suborder_num, int suborder_xyz_n[],
1269  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1270 
1271 /******************************************************************************/
1272 /*
1273  Purpose:
1274 
1275  TRIANGLE_NCO_SUBRULE_03 returns a compressed NCO rule 3.
1276 
1277  Licensing:
1278 
1279  This code is distributed under the GNU LGPL license.
1280 
1281  Modified:
1282 
1283  30 January 2007
1284 
1285  Author:
1286 
1287  John Burkardt
1288 
1289  Reference:
1290 
1291  Peter Silvester,
1292  Symmetric Quadrature Formulae for Simplexes,
1293  Mathematics of Computation,
1294  Volume 24, Number 109, January 1970, pages 95-100.
1295 
1296  Parameters:
1297 
1298  Input, int SUBORDER_NUM, the number of suborders of the rule.
1299 
1300  Output, int SUBORDER_XYZ_N[3*SUBORDER_NUM],
1301  the numerators of the barycentric coordinates of the abscissas.
1302 
1303  Output, int *SUBORDER_XYZ_D,
1304  the denominator of the barycentric coordinates of the abscissas.
1305 
1306  Output, int SUBORDER_W_N[SUBORDER_NUM],
1307  the numerator of the suborder weights.
1308 
1309  Output, int SUBORDER_W_D,
1310  the denominator of the suborder weights.
1311 */
1312 {
1313  int i;
1314  int s;
1315  int suborder_xyz_n_03[3*2] = {
1316  2, 0, 0,
1317  1, 1, 0
1318  };
1319  int suborder_xyz_d_03 = 2;
1320  int suborder_w_n_03[2] = { 7, -3 };
1321  int suborder_w_d_03 = 12;
1322 
1323  for ( s = 0; s < suborder_num; s++ )
1324  {
1325  for ( i = 0; i < 3; i++ )
1326  {
1327  suborder_xyz_n[i+s*3] = suborder_xyz_n_03[i+s*3];
1328  }
1329  }
1330  *suborder_xyz_d = suborder_xyz_d_03;
1331 
1332  for ( s = 0; s < suborder_num; s++ )
1333  {
1334  suborder_w_n[s] = suborder_w_n_03[s];
1335  }
1336  *suborder_w_d = suborder_w_d_03;
1337 
1338  return;
1339 }
1340 /******************************************************************************/
1341 
1342 void triangle_nco_subrule_04 ( int suborder_num, int suborder_xyz_n[],
1343  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1344 
1345 /******************************************************************************/
1346 /*
1347  Purpose:
1348 
1349  TRIANGLE_NCO_SUBRULE_04 returns a compressed NCO rule 4.
1350 
1351  Licensing:
1352 
1353  This code is distributed under the GNU LGPL license.
1354 
1355  Modified:
1356 
1357  30 January 2007
1358 
1359  Author:
1360 
1361  John Burkardt
1362 
1363  Reference:
1364 
1365  Peter Silvester,
1366  Symmetric Quadrature Formulae for Simplexes,
1367  Mathematics of Computation,
1368  Volume 24, Number 109, January 1970, pages 95-100.
1369 
1370  Parameters:
1371 
1372  Input, int SUBORDER_NUM, the number of suborders of the rule.
1373 
1374  Output, int SUBORDER_XYZ_N[3*SUBORDER_NUM],
1375  the numerators of the barycentric coordinates of the abscissas.
1376 
1377  Output, int *SUBORDER_XYZ_D,
1378  the denominator of the barycentric coordinates of the abscissas.
1379 
1380  Output, int SUBORDER_W_N[SUBORDER_NUM],
1381  the numerator of the suborder weights.
1382 
1383  Output, int SUBORDER_W_D,
1384  the denominator of the suborder weights.
1385 */
1386 {
1387  int i;
1388  int s;
1389  int suborder_xyz_n_04[3*3] = {
1390  3, 0, 0,
1391  2, 1, 0,
1392  1, 1, 1
1393  };
1394  int suborder_xyz_d_04 = 3;
1395  int suborder_w_n_04[3] = { 8, 3, -12 };
1396  int suborder_w_d_04 = 30;
1397 
1398  for ( s = 0; s < suborder_num; s++ )
1399  {
1400  for ( i = 0; i < 3; i++ )
1401  {
1402  suborder_xyz_n[i+s*3] = suborder_xyz_n_04[i+s*3];
1403  }
1404  }
1405  *suborder_xyz_d = suborder_xyz_d_04;
1406 
1407  for ( s = 0; s < suborder_num; s++ )
1408  {
1409  suborder_w_n[s] = suborder_w_n_04[s];
1410  }
1411  *suborder_w_d = suborder_w_d_04;
1412 
1413  return;
1414 }
1415 /******************************************************************************/
1416 
1417 void triangle_nco_subrule_05 ( int suborder_num, int suborder_xyz_n[],
1418  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1419 
1420 /******************************************************************************/
1421 /*
1422  Purpose:
1423 
1424  TRIANGLE_NCO_SUBRULE_05 returns a compressed NCO rule 5.
1425 
1426  Licensing:
1427 
1428  This code is distributed under the GNU LGPL license.
1429 
1430  Modified:
1431 
1432  30 January 2007
1433 
1434  Author:
1435 
1436  John Burkardt
1437 
1438  Reference:
1439 
1440  Peter Silvester,
1441  Symmetric Quadrature Formulae for Simplexes,
1442  Mathematics of Computation,
1443  Volume 24, Number 109, January 1970, pages 95-100.
1444 
1445  Parameters:
1446 
1447  Input, int SUBORDER_NUM, the number of suborders of the rule.
1448 
1449  Output, int SUBORDER_XYZ_N[3*SUBORDER_NUM],
1450  the numerators of the barycentric coordinates of the abscissas.
1451 
1452  Output, int *SUBORDER_XYZ_D,
1453  the denominator of the barycentric coordinates of the abscissas.
1454 
1455  Output, int SUBORDER_W_N[SUBORDER_NUM],
1456  the numerator of the suborder weights.
1457 
1458  Output, int SUBORDER_W_D,
1459  the denominator of the suborder weights.
1460 */
1461 {
1462  int i;
1463  int s;
1464  int suborder_xyz_n_05[3*4] = {
1465  4, 0, 0,
1466  3, 1, 0,
1467  2, 2, 0,
1468  2, 1, 1
1469  };
1470  int suborder_xyz_d_05 = 4;
1471  int suborder_w_n_05[4] = { 307, -316, 629, -64 };
1472  int suborder_w_d_05 = 720;
1473 
1474  for ( s = 0; s < suborder_num; s++ )
1475  {
1476  for ( i = 0; i < 3; i++ )
1477  {
1478  suborder_xyz_n[i+s*3] = suborder_xyz_n_05[i+s*3];
1479  }
1480  }
1481  *suborder_xyz_d = suborder_xyz_d_05;
1482 
1483  for ( s = 0; s < suborder_num; s++ )
1484  {
1485  suborder_w_n[s] = suborder_w_n_05[s];
1486  }
1487  *suborder_w_d = suborder_w_d_05;
1488 
1489  return;
1490 }
1491 /******************************************************************************/
1492 
1493 void triangle_nco_subrule_06 ( int suborder_num, int suborder_xyz_n[],
1494  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1495 
1496 /******************************************************************************/
1497 /*
1498  Purpose:
1499 
1500  TRIANGLE_NCO_SUBRULE_06 returns a compressed NCO rule 6.
1501 
1502  Licensing:
1503 
1504  This code is distributed under the GNU LGPL license.
1505 
1506  Modified:
1507 
1508  30 January 2007
1509 
1510  Author:
1511 
1512  John Burkardt
1513 
1514  Reference:
1515 
1516  Peter Silvester,
1517  Symmetric Quadrature Formulae for Simplexes,
1518  Mathematics of Computation,
1519  Volume 24, Number 109, January 1970, pages 95-100.
1520 
1521  Parameters:
1522 
1523  Input, int SUBORDER_NUM, the number of suborders of the rule.
1524 
1525  Output, int SUBORDER_XYZ_N[3*SUBORDER_NUM],
1526  the numerators of the barycentric coordinates of the abscissas.
1527 
1528  Output, int *SUBORDER_XYZ_D,
1529  the denominator of the barycentric coordinates of the abscissas.
1530 
1531  Output, int SUBORDER_W_N[SUBORDER_NUM],
1532  the numerator of the suborder weights.
1533 
1534  Output, int SUBORDER_W_D,
1535  the denominator of the suborder weights.
1536 */
1537 {
1538  int i;
1539  int s;
1540  int suborder_xyz_n_06[3*5] = {
1541  5, 0, 0,
1542  4, 1, 0,
1543  3, 2, 0,
1544  3, 1, 1,
1545  2, 2, 1
1546  };
1547  int suborder_xyz_d_06 = 5;
1548  int suborder_w_n_06[5] = { 71, -13, 57, -167, 113 };
1549  int suborder_w_d_06 = 315;
1550 
1551  for ( s = 0; s < suborder_num; s++ )
1552  {
1553  for ( i = 0; i < 3; i++ )
1554  {
1555  suborder_xyz_n[i+s*3] = suborder_xyz_n_06[i+s*3];
1556  }
1557  }
1558  *suborder_xyz_d = suborder_xyz_d_06;
1559 
1560  for ( s = 0; s < suborder_num; s++ )
1561  {
1562  suborder_w_n[s] = suborder_w_n_06[s];
1563  }
1564  *suborder_w_d = suborder_w_d_06;
1565 
1566  return;
1567 }
1568 /******************************************************************************/
1569 
1570 void triangle_nco_subrule_07 ( int suborder_num, int suborder_xyz_n[],
1571  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1572 
1573 /******************************************************************************/
1574 /*
1575  Purpose:
1576 
1577  TRIANGLE_NCO_SUBRULE_07 returns a compressed NCO rule 7.
1578 
1579  Licensing:
1580 
1581  This code is distributed under the GNU LGPL license.
1582 
1583  Modified:
1584 
1585  30 January 2007
1586 
1587  Author:
1588 
1589  John Burkardt
1590 
1591  Reference:
1592 
1593  Peter Silvester,
1594  Symmetric Quadrature Formulae for Simplexes,
1595  Mathematics of Computation,
1596  Volume 24, Number 109, January 1970, pages 95-100.
1597 
1598  Parameters:
1599 
1600  Input, int SUBORDER_NUM, the number of suborders of the rule.
1601 
1602  Output, int SUBORDER_XYZ_N[3*SUBORDER_NUM],
1603  the numerators of the barycentric coordinates of the abscissas.
1604 
1605  Output, int *SUBORDER_XYZ_D,
1606  the denominator of the barycentric coordinates of the abscissas.
1607 
1608  Output, int SUBORDER_W_N[SUBORDER_NUM],
1609  the numerator of the suborder weights.
1610 
1611  Output, int SUBORDER_W_D,
1612  the denominator of the suborder weights.
1613 */
1614 {
1615  int i;
1616  int s;
1617  int suborder_xyz_n_07[3*7] = {
1618  6, 0, 0,
1619  5, 1, 0,
1620  4, 2, 0,
1621  4, 1, 1,
1622  3, 3, 0,
1623  3, 2, 1,
1624  2, 2, 2
1625  };
1626  int suborder_xyz_d_07 = 6;
1627  int suborder_w_n_07[7] = { 767, -1257, 2901, 387, -3035, -915, 3509 };
1628  int suborder_w_d_07 = 2240;
1629 
1630  for ( s = 0; s < suborder_num; s++ )
1631  {
1632  for ( i = 0; i < 3; i++ )
1633  {
1634  suborder_xyz_n[i+s*3] = suborder_xyz_n_07[i+s*3];
1635  }
1636  }
1637  *suborder_xyz_d = suborder_xyz_d_07;
1638 
1639  for ( s = 0; s < suborder_num; s++ )
1640  {
1641  suborder_w_n[s] = suborder_w_n_07[s];
1642  }
1643  *suborder_w_d = suborder_w_d_07;
1644 
1645  return;
1646 }
1647 /******************************************************************************/
1648 
1649 void triangle_nco_subrule_08 ( int suborder_num, int suborder_xyz_n[],
1650  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1651 
1652 /******************************************************************************/
1653 /*
1654  Purpose:
1655 
1656  TRIANGLE_NCO_SUBRULE_08 returns a compressed NCO rule 8.
1657 
1658  Licensing:
1659 
1660  This code is distributed under the GNU LGPL license.
1661 
1662  Modified:
1663 
1664  30 January 2007
1665 
1666  Author:
1667 
1668  John Burkardt
1669 
1670  Reference:
1671 
1672  Peter Silvester,
1673  Symmetric Quadrature Formulae for Simplexes,
1674  Mathematics of Computation,
1675  Volume 24, Number 109, January 1970, pages 95-100.
1676 
1677  Parameters:
1678 
1679  Input, int SUBORDER_NUM, the number of suborders of the rule.
1680 
1681  Output, int SUBORDER_XYZ_N[3*SUBORDER_NUM],
1682  the numerators of the barycentric coordinates of the abscissas.
1683 
1684  Output, int *SUBORDER_XYZ_D,
1685  the denominator of the barycentric coordinates of the abscissas.
1686 
1687  Output, int SUBORDER_W_N[SUBORDER_NUM],
1688  the numerator of the suborder weights.
1689 
1690  Output, int SUBORDER_W_D,
1691  the denominator of the suborder weights.
1692 */
1693 {
1694  int i;
1695  int s;
1696  int suborder_xyz_n_08[3*8] = {
1697  7, 0, 0,
1698  6, 1, 0,
1699  5, 2, 0,
1700  5, 1, 1,
1701  4, 3, 0,
1702  4, 2, 1,
1703  3, 3, 1,
1704  3, 2, 2
1705  };
1706  int suborder_xyz_d_08 = 7;
1707  int suborder_w_n_08[8] = { 898, -662, 1573, -2522, -191, 2989, -5726, 1444 };
1708  int suborder_w_d_08 = 4536;
1709 
1710  for ( s = 0; s < suborder_num; s++ )
1711  {
1712  for ( i = 0; i < 3; i++ )
1713  {
1714  suborder_xyz_n[i+s*3] = suborder_xyz_n_08[i+s*3];
1715  }
1716  }
1717  *suborder_xyz_d = suborder_xyz_d_08;
1718 
1719  for ( s = 0; s < suborder_num; s++ )
1720  {
1721  suborder_w_n[s] = suborder_w_n_08[s];
1722  }
1723  *suborder_w_d = suborder_w_d_08;
1724 
1725  return;
1726 }
1727 /******************************************************************************/
1728 
1729 void triangle_nco_subrule_09 ( int suborder_num, int suborder_xyz_n[],
1730  int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d )
1731 
1732 /******************************************************************************/
1733 /*
1734  Purpose:
1735 
1736  TRIANGLE_NCO_SUBRULE_09 returns a compressed NCO rule 9.
1737 
1738  Licensing:
1739 
1740  This code is distributed under the GNU LGPL license.
1741 
1742  Modified:
1743 
1744  30 January 2007
1745 
1746  Author:
1747 
1748  John Burkardt
1749 
1750  Reference:
1751 
1752  Peter Silvester,
1753  Symmetric Quadrature Formulae for Simplexes,
1754  Mathematics of Computation,
1755  Volume 24, Number 109, January 1970, pages 95-100.
1756 
1757  Parameters:
1758 
1759  Input, int SUBORDER_NUM, the number of suborders of the rule.
1760 
1761  Output, int SUBORDER_XYZ_N[3*SUBORDER_NUM],
1762  the numerators of the barycentric coordinates of the abscissas.
1763 
1764  Output, int *SUBORDER_XYZ_D,
1765  the denominator of the barycentric coordinates of the abscissas.
1766 
1767  Output, int SUBORDER_W_N[SUBORDER_NUM],
1768  the numerator of the suborder weights.
1769 
1770  Output, int SUBORDER_W_D,
1771  the denominator of the suborder weights.
1772 */
1773 {
1774  int i;
1775  int s;
1776  int suborder_xyz_n_09[3*10] = {
1777  8, 0, 0,
1778  7, 1, 0,
1779  6, 2, 0,
1780  6, 1, 1,
1781  5, 3, 0,
1782  5, 2, 1,
1783  4, 4, 0,
1784  4, 3, 1,
1785  4, 2, 2,
1786  3, 3, 2
1787  };
1788  int suborder_xyz_d_09 = 8;
1789  int suborder_w_n_09[10] = {
1790  1051445, -2366706, 6493915, 1818134, -9986439,-3757007, 12368047,
1791  478257, 10685542, -6437608 };
1792  int suborder_w_d_09 = 3628800;
1793 
1794  for ( s = 0; s < suborder_num; s++ )
1795  {
1796  for ( i = 0; i < 3; i++ )
1797  {
1798  suborder_xyz_n[i+s*3] = suborder_xyz_n_09[i+s*3];
1799  }
1800  }
1801  *suborder_xyz_d = suborder_xyz_d_09;
1802 
1803  for ( s = 0; s < suborder_num; s++ )
1804  {
1805  suborder_w_n[s] = suborder_w_n_09[s];
1806  }
1807  *suborder_w_d = suborder_w_d_09;
1808 
1809  return;
1810 }
int i4_modp(int i, int j)
void triangle_nco_subrule_04(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
int triangle_nco_suborder_num(int rule)
void triangle_nco_subrule_01(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
int r8_nint(double x)
int triangle_nco_rule_num()
double r8_huge()
int * triangle_nco_suborder(int rule, int suborder_num)
void triangle_nco_subrule_02(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
double triangle_area(double t[2 *3])
void triangle_nco_subrule_05(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
void triangle_nco_subrule_09(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
int triangle_nco_order_num(int rule)
void triangle_nco_subrule_03(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
int i4_wrap(int ival, int ilo, int ihi)
void reference_to_physical_t3(double t[], int n, double ref[], double phy[])
void triangle_nco_subrule(int rule, int suborder_num, double suborder_xyz[], double suborder_w[])
void triangle_nco_rule(int rule, int order_num, double xy[], double w[])
int triangle_nco_degree(int rule)
void triangle_nco_subrule_08(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
void triangle_nco_subrule_07(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)
void triangle_nco_subrule_06(int suborder_num, int suborder_xyz_n[], int *suborder_xyz_d, int suborder_w_n[], int *suborder_w_d)