ROL
ROL_SandiaRules2Def.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ************************************************************************
3 //
4 // Rapid Optimization Library (ROL) Package
5 // Copyright (2014) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact lead developers:
38 // Drew Kouri (dpkouri@sandia.gov) and
39 // Denis Ridzal (dridzal@sandia.gov)
40 //
41 // ************************************************************************
42 // @HEADER
43 
44 
45 namespace ROL {
46 
47 //****************************************************************************80
48 void SandiaRules2::ccn_points ( int n, int dim, double x[] )
49 //****************************************************************************80
50 //
51 // Purpose:
52 //
53 // CCN_POINTS computes nested Clenshaw Curtis quadrature points.
54 //
55 // Licensing:
56 //
57 // This code is distributed under the GNU LGPL license.
58 //
59 // Modified:
60 //
61 // 07 March 2011
62 //
63 // Author:
64 //
65 // John Burkardt
66 //
67 // Parameters:
68 //
69 // Input, int N, the order.
70 //
71 // Input, int DIM, the spatial dimension represented by this rule,
72 // in cases where a multidimensional product rule is being formed.
73 //
74 // Output, double X[N], the abscissas.
75 //
76 {
78  return;
79 }
80 
81 //****************************************************************************80
82 void SandiaRules2::ccn_weights ( int n, int dim, double w[] )
83 //****************************************************************************80
84 //
85 // Purpose:
86 //
87 // CCN_WEIGHTS computes nested Clenshaw Curtis quadrature weights.
88 //
89 // Licensing:
90 //
91 // This code is distributed under the GNU LGPL license.
92 //
93 // Modified:
94 //
95 // 07 March 2011
96 //
97 // Author:
98 //
99 // John Burkardt
100 //
101 // Parameters:
102 //
103 // Input, int N, the order.
104 //
105 // Input, int DIM, the spatial dimension represented by this rule,
106 // in cases where a multidimensional product rule is being formed.
107 //
108 // Output, double W[N], the weights.
109 //
110 {
112  return;
113 }
114 
115 //****************************************************************************80
116 void SandiaRules2::clenshaw_curtis_points ( int n, int dim, double x[] )
117 //****************************************************************************80
118 //
119 // Purpose:
120 //
121 // CLENSHAW_CURTIS_POINTS computes Clenshaw Curtis quadrature points.
122 //
123 // Discussion:
124 //
125 // Our convention is that the abscissas are numbered from left to right.
126 //
127 // This rule is defined on [-1,1].
128 //
129 // Licensing:
130 //
131 // This code is distributed under the GNU LGPL license.
132 //
133 // Modified:
134 //
135 // 02 April 2010
136 //
137 // Author:
138 //
139 // John Burkardt
140 //
141 // Parameters:
142 //
143 // Input, int N, the order.
144 //
145 // Input, int DIM, the spatial dimension represented by this rule,
146 // in cases where a multidimensional product rule is being formed.
147 //
148 // Output, double X[N], the abscissas.
149 //
150 {
151  int index;
152  double pi = 3.141592653589793;
153 
154  if ( n < 1 )
155  {
156  std::cerr << "\n";
157  std::cerr << "CLENSHAW_CURTIS_POINTS - Fatal error!\n";
158  std::cerr << " Order is less than 1.\n";
159  std::exit ( 1 );
160  }
161  else if ( n == 1 )
162  {
163  x[0] = 0.0;
164  }
165  else
166  {
167  for ( index = 1; index <= n; index++ )
168  {
169  x[index-1] = std::cos ( ( double ) ( n - index ) * pi
170  / ( double ) ( n - 1 ) );
171  }
172  x[0] = -1.0;
173  if ( ( n % 2 ) == 1 )
174  {
175  x[(n-1)/2] = 0.0;
176  }
177  x[n-1] = +1.0;
178  }
179  return;
180 }
181 
182 //****************************************************************************80
183 void SandiaRules2::clenshaw_curtis_weights ( int n, int dim, double w[] )
184 //****************************************************************************80
185 //
186 // Purpose:
187 //
188 // CLENSHAW_CURTIS_WEIGHTS computes Clenshaw Curtis quadrature weights.
189 //
190 // Licensing:
191 //
192 // This code is distributed under the GNU LGPL license.
193 //
194 // Modified:
195 //
196 // 02 April 2010
197 //
198 // Author:
199 //
200 // John Burkardt
201 //
202 // Parameters:
203 //
204 // Input, int N, the order.
205 //
206 // Input, int DIM, the spatial dimension represented by this rule,
207 // in cases where a multidimensional product rule is being formed.
208 //
209 // Output, double W[N], the weights.
210 //
211 {
212  double b;
213  int i;
214  int j;
215  double pi = 3.141592653589793;
216  double theta;
217 
218  if ( n < 1 )
219  {
220  std::cerr << "\n";
221  std::cerr << "CLENSHAW_CURTIS_WEIGHTS - Fatal error!\n";
222  std::cerr << " Order is less than 1.\n";
223  std::exit ( 1 );
224  }
225  else if ( n == 1 )
226  {
227  w[0] = 2.0;
228  return;
229  }
230 
231  for ( i = 1; i <= n; i++ )
232  {
233  theta = ( double ) ( i - 1 ) * pi / ( double ) ( n - 1 );
234 
235  w[i-1] = 1.0;
236 
237  for ( j = 1; j <= ( n - 1 ) / 2; j++ )
238  {
239  if ( 2 * j == ( n - 1 ) )
240  {
241  b = 1.0;
242  }
243  else
244  {
245  b = 2.0;
246  }
247 
248  w[i-1] = w[i-1] - b * std::cos ( 2.0 * ( double ) ( j ) * theta )
249  / ( double ) ( 4 * j * j - 1 );
250  }
251  }
252 
253  w[0] = w[0] / ( double ) ( n - 1 );
254  for ( i = 1; i < n - 1; i++ )
255  {
256  w[i] = 2.0 * w[i] / ( double ) ( n - 1 );
257  }
258  w[n-1] = w[n-1] / ( double ) ( n - 1 );
259 
260  return;
261 }
262 
263 //****************************************************************************80
264 void SandiaRules2::fejer2_points ( int n, int dim, double x[] )
265 //****************************************************************************80
266 //
267 // Purpose:
268 //
269 // FEJER2_POINTS computes Fejer type 2 quadrature points.
270 //
271 // Discussion:
272 //
273 // Our convention is that the abscissas are numbered from left to right.
274 //
275 // The rule is defined on [-1,1].
276 //
277 // Licensing:
278 //
279 // This code is distributed under the GNU LGPL license.
280 //
281 // Modified:
282 //
283 // 02 April 2010
284 //
285 // Author:
286 //
287 // John Burkardt
288 //
289 // Parameters:
290 //
291 // Input, int N, the order.
292 //
293 // Input, int DIM, the spatial dimension represented by this rule,
294 // in cases where a multidimensional product rule is being formed.
295 //
296 // Output, double X[N], the abscissas.
297 //
298 {
299  int index;
300  double pi = 3.141592653589793;
301 
302  if ( n < 1 )
303  {
304  std::cerr << "\n";
305  std::cerr << "FEJER2_POINTS - Fatal error!\n";
306  std::cerr << " Order is less than 1.\n";
307  std::exit ( 1 );
308  }
309  else if ( n == 1 )
310  {
311  x[0] = 0.0;
312  }
313  else
314  {
315  for ( index = 1; index <= n; index++ )
316  {
317  x[index-1] = std::cos ( ( double ) ( n + 1 - index ) * pi
318  / ( double ) ( n + 1 ) );
319  }
320  if ( ( n % 2 ) == 1 )
321  {
322  x[(n-1)/2] = 0.0;
323  }
324  }
325  return;
326 }
327 
328 //****************************************************************************80
329 void SandiaRules2::fejer2_weights ( int n, int dim, double w[] )
330 //****************************************************************************80
331 //
332 // Purpose:
333 //
334 // FEJER2_WEIGHTS computes Fejer type 2 quadrature weights.
335 //
336 // Licensing:
337 //
338 // This code is distributed under the GNU LGPL license.
339 //
340 // Modified:
341 //
342 // 02 April 2010
343 //
344 // Author:
345 //
346 // John Burkardt
347 //
348 // Parameters:
349 //
350 // Input, int N, the order.
351 //
352 // Input, int DIM, the spatial dimension represented by this rule,
353 // in cases where a multidimensional product rule is being formed.
354 //
355 // Output, double W[N], the weights.
356 //
357 {
358  int i;
359  int j;
360  double p;
361  double pi = 3.141592653589793;
362  double theta;
363 
364  if ( n < 1 )
365  {
366  std::cerr << "\n";
367  std::cerr << "FEJER2_WEIGHTS - Fatal error!\n";
368  std::cerr << " N < 1.\n";
369  std::exit ( 1 );
370  }
371  else if ( n == 1 )
372  {
373  w[0] = 2.0;
374  }
375  else if ( n == 2 )
376  {
377  w[0] = 1.0;
378  w[1] = 1.0;
379  }
380  else
381  {
382  for ( i = 1; i <= n; i++ )
383  {
384  theta = ( double ) ( n + 1 - i ) * pi
385  / ( double ) ( n + 1 );
386 
387  w[i-1] = 1.0;
388 
389  for ( j = 1; j <= ( ( n - 1 ) / 2 ); j++ )
390  {
391  w[i-1] = w[i-1] - 2.0 * std::cos ( 2.0 * ( double ) ( j ) * theta )
392  / ( double ) ( 4 * j * j - 1 );
393  }
394  p = 2.0 * ( double ) ( ( ( n + 1 ) / 2 ) ) - 1.0;
395  w[i-1] = w[i-1] - std::cos ( ( p + 1.0 ) * theta ) / p;
396  }
397  for ( i = 0; i < n; i++ )
398  {
399  w[i] = 2.0 * w[i] / ( double ) ( n + 1 );
400  }
401  }
402  return;
403 }
404 
405 //****************************************************************************80
406 void SandiaRules2::gen_hermite_points ( int n, int dim, double x[] )
407 //****************************************************************************80
408 //
409 // Purpose:
410 //
411 // GEN_HERMITE_POINTS: Generalized Hermite quadrature points.
412 //
413 // Discussion:
414 //
415 // This function assumes the existence of a function:
416 // double parameter ( int dim, int offset )
417 // which can supply the value of the ALPHA parameter.
418 //
419 // Licensing:
420 //
421 // This code is distributed under the GNU LGPL license.
422 //
423 // Modified:
424 //
425 // 02 April 2010
426 //
427 // Author:
428 //
429 // John Burkardt
430 //
431 // Parameters:
432 //
433 // Input, int N, the order.
434 //
435 // Input, int DIM, the spatial dimension represented by this rule,
436 // in cases where a multidimensional product rule is being formed.
437 //
438 // Output, double X[N], the abscissas.
439 //
440 {
441  double alpha;
442 
443  alpha = parameter ( dim, 0 );
444 
446 
447  return;
448 }
449 
450 //****************************************************************************80
451 void SandiaRules2::gen_hermite_weights ( int n, int dim, double w[] )
452 //****************************************************************************80
453 //
454 // Purpose:
455 //
456 // GEN_HERMITE_WEIGHTS: Generalized Hermite quadrature weights.
457 //
458 // Discussion:
459 //
460 // This function assumes the existence of a function:
461 // double parameter ( int dim, int offset )
462 // which can supply the value of the ALPHA parameter.
463 //
464 // Licensing:
465 //
466 // This code is distributed under the GNU LGPL license.
467 //
468 // Modified:
469 //
470 // 02 April 2010
471 //
472 // Author:
473 //
474 // John Burkardt
475 //
476 // Parameters:
477 //
478 // Input, int N, the order.
479 //
480 // Input, int DIM, the spatial dimension represented by this rule,
481 // in cases where a multidimensional product rule is being formed.
482 //
483 // Output, double W[N], the weights.
484 //
485 {
486  double alpha;
487  alpha = parameter ( dim, 0 );
489  return;
490 }
491 
492 //****************************************************************************80
493 void SandiaRules2::gen_laguerre_points ( int n, int dim, double x[] )
494 //****************************************************************************80
495 //
496 // Purpose:
497 //
498 // GEN_LAGUERRE_POINTS: Generalized Laguerre quadrature points.
499 //
500 // Discussion:
501 //
502 // This function assumes the existence of a function:
503 // double parameter ( int dim, int offset )
504 // which can supply the value of the ALPHA parameter.
505 //
506 // Licensing:
507 //
508 // This code is distributed under the GNU LGPL license.
509 //
510 // Modified:
511 //
512 // 02 April 2010
513 //
514 // Author:
515 //
516 // John Burkardt
517 //
518 // Parameters:
519 //
520 // Input, int N, the order.
521 //
522 // Input, int DIM, the spatial dimension represented by this rule,
523 // in cases where a multidimensional product rule is being formed.
524 //
525 // Output, double X[N], the abscissas.
526 //
527 {
528  double alpha;
529  alpha = parameter ( dim, 0 );
531  return;
532 }
533 
534 //****************************************************************************80
535 void SandiaRules2::gen_laguerre_weights ( int n, int dim, double w[] )
536 //****************************************************************************80
537 //
538 // Purpose:
539 //
540 // GEN_LAGUERRE_WEIGHTS: Generalized Laguerre quadrature weights.
541 //
542 // Discussion:
543 //
544 // This function assumes the existence of a function:
545 // double parameter ( int dim, int offset )
546 // which can supply the value of the ALPHA parameter.
547 //
548 // Licensing:
549 //
550 // This code is distributed under the GNU LGPL license.
551 //
552 // Modified:
553 //
554 // 02 April 2010
555 //
556 // Author:
557 //
558 // John Burkardt
559 //
560 // Parameters:
561 //
562 // Input, int N, the order.
563 //
564 // Input, int DIM, the spatial dimension represented by this rule,
565 // in cases where a multidimensional product rule is being formed.
566 //
567 // Output, double W[N], the weights.
568 //
569 {
570  double alpha;
571  alpha = parameter ( dim, 0 );
573  return;
574 }
575 
576 //****************************************************************************80
577 void SandiaRules2::hcc_points ( int n, int dim, double x[] )
578 //****************************************************************************80
579 //
580 // Purpose:
581 //
582 // HCC_POINTS computes Hermite-Cubic-Chebyshev-Spacing quadrature points.
583 //
584 // Discussion:
585 //
586 // For the HCC rule, we assume that an interval has been divided by
587 // M nodes X into Chebyshev-spaced subintervals, and that at each
588 // abscissa both function and derivative information is available.
589 // The piecewise cubic Hermite interpolant is constructed for this data.
590 // The quadrature rule uses N = 2 * M abscissas, where each node is
591 // listed twice, and the weights occur in pairs, with the first multiplying
592 // the function value and the second the derivative.
593 //
594 // Licensing:
595 //
596 // This code is distributed under the GNU LGPL license.
597 //
598 // Modified:
599 //
600 // 17 March 2011
601 //
602 // Author:
603 //
604 // John Burkardt
605 //
606 // Parameters:
607 //
608 // Input, int N, the order.
609 //
610 // Input, int DIM, the spatial dimension represented by this rule,
611 // in cases where a multidimensional product rule is being formed.
612 //
613 // Output, double X[N], the abscissas.
614 //
615 {
617  return;
618 }
619 
620 //****************************************************************************80
621 void SandiaRules2::hcc_weights ( int n, int dim, double w[] )
622 //****************************************************************************80
623 //
624 // Purpose:
625 //
626 // HCC_WEIGHTS computes Hermite-Cubic-Chebyshev-Spacing quadrature weights.
627 //
628 // Discussion:
629 //
630 // For the HCC rule, we assume that an interval has been divided by
631 // M nodes X into Chebyshev-spaced subintervals, and that at each
632 // abscissa both function and derivative information is available.
633 // The piecewise cubic Hermite interpolant is constructed for this data.
634 // The quadrature rule uses N = 2 * M abscissas, where each node is
635 // listed twice, and the weights occur in pairs, with the first multiplying
636 // the function value and the second the derivative.
637 //
638 // Licensing:
639 //
640 // This code is distributed under the GNU LGPL license.
641 //
642 // Modified:
643 //
644 // 17 March 2011
645 //
646 // Author:
647 //
648 // John Burkardt
649 //
650 // Parameters:
651 //
652 // Input, int N, the order.
653 //
654 // Input, int DIM, the spatial dimension represented by this rule,
655 // in cases where a multidimensional product rule is being formed.
656 //
657 // Output, double W[N], the weights.
658 //
659 {
661  return;
662 }
663 
664 //****************************************************************************80
665 void SandiaRules2::hce_points ( int n, int dim, double x[] )
666 //****************************************************************************80
667 //
668 // Purpose:
669 //
670 // HCE_POINTS computes Hermite-Cubic-Equal-Spacing quadrature points.
671 //
672 // Discussion:
673 //
674 // For the HCE rule, we assume that an interval has been divided by
675 // M nodes X into equally spaced subintervals, and that at each
676 // abscissa both function and derivative information is available.
677 // The piecewise cubic Hermite interpolant is constructed for this data.
678 // The quadrature rule uses N = 2 * M abscissas, where each node is
679 // listed twice, and the weights occur in pairs, with the first multiplying
680 // the function value and the second the derivative.
681 //
682 // Licensing:
683 //
684 // This code is distributed under the GNU LGPL license.
685 //
686 // Modified:
687 //
688 // 07 March 2011
689 //
690 // Author:
691 //
692 // John Burkardt
693 //
694 // Parameters:
695 //
696 // Input, int N, the order.
697 //
698 // Input, int DIM, the spatial dimension represented by this rule,
699 // in cases where a multidimensional product rule is being formed.
700 //
701 // Output, double X[N], the abscissas.
702 //
703 {
705  return;
706 }
707 
708 //****************************************************************************80
709 void SandiaRules2::hce_weights ( int n, int dim, double w[] )
710 //****************************************************************************80
711 //
712 // Purpose:
713 //
714 // HCE_WEIGHTS computes Hermite-Cubic-Equal-Spacing quadrature weights.
715 //
716 // Discussion:
717 //
718 // For the HCE rule, we assume that an interval has been divided by
719 // M nodes X into equally spaced subintervals, and that at each
720 // abscissa both function and derivative information is available.
721 // The piecewise cubic Hermite interpolant is constructed for this data.
722 // The quadrature rule uses N = 2 * M abscissas, where each node is
723 // listed twice, and the weights occur in pairs, with the first multiplying
724 // the function value and the second the derivative.
725 //
726 // Licensing:
727 //
728 // This code is distributed under the GNU LGPL license.
729 //
730 // Modified:
731 //
732 // 07 March 2011
733 //
734 // Author:
735 //
736 // John Burkardt
737 //
738 // Parameters:
739 //
740 // Input, int N, the order.
741 //
742 // Input, int DIM, the spatial dimension represented by this rule,
743 // in cases where a multidimensional product rule is being formed.
744 //
745 // Output, double W[N], the weights.
746 //
747 {
749  return;
750 }
751 
752 //****************************************************************************80
753 void SandiaRules2::hermite_genz_keister_points ( int n, int dim, double x[] )
754 //****************************************************************************80
755 //
756 // Purpose:
757 //
758 // HERMITE_GENZ_KEISTER_POINTS looks up Genz-Keister Hermite abscissas.
759 //
760 // Licensing:
761 //
762 // This code is distributed under the GNU LGPL license.
763 //
764 // Modified:
765 //
766 // 07 June 2010
767 //
768 // Author:
769 //
770 // John Burkardt
771 //
772 // Parameters:
773 //
774 // Input, int N, the order.
775 // The allowed orders are 1, 3, 9, 19 and 35.
776 //
777 // Input, int DIM, the spatial dimension represented by this rule,
778 // in cases where a multidimensional product rule is being formed.
779 //
780 // Output, double X[N], the abscissas.
781 //
782 {
784  return;
785 }
786 
787 //****************************************************************************80
788 void SandiaRules2::hermite_genz_keister_weights ( int n, int dim, double w[] )
789 //****************************************************************************80
790 //
791 // Purpose:
792 //
793 // HERMITE_GENZ_KEISTER_WEIGHTS looks up Genz-Keister Hermite weights.
794 //
795 // Licensing:
796 //
797 // This code is distributed under the GNU LGPL license.
798 //
799 // Modified:
800 //
801 // 07 June 2010
802 //
803 // Author:
804 //
805 // John Burkardt
806 //
807 // Parameters:
808 //
809 // Input, int N, the order.
810 // The allowed orders are 1, 3, 9, 19 and 35.
811 //
812 // Input, int DIM, the spatial dimension represented by this rule,
813 // in cases where a multidimensional product rule is being formed.
814 //
815 // Output, double W[N], the weights.
816 //
817 {
819  return;
820 }
821 
822 //****************************************************************************80
823 void SandiaRules2::hermite_points ( int n, int dim, double x[] )
824 //****************************************************************************80
825 //
826 // Purpose:
827 //
828 // HERMITE_POINTS computes Hermite quadrature points.
829 //
830 // Licensing:
831 //
832 // This code is distributed under the GNU LGPL license.
833 //
834 // Modified:
835 //
836 // 02 April 2010
837 //
838 // Author:
839 //
840 // John Burkardt
841 //
842 // Parameters:
843 //
844 // Input, int N, the order.
845 //
846 // Input, int DIM, the spatial dimension represented by this rule,
847 // in cases where a multidimensional product rule is being formed.
848 //
849 // Output, double X[N], the abscissas.
850 //
851 {
853  return;
854 }
855 
856 //****************************************************************************80
857 void SandiaRules2::hermite_weights ( int n, int dim, double w[] )
858 //****************************************************************************80
859 //
860 // Purpose:
861 //
862 // HERMITE_WEIGHTS computes Hermite quadrature weights.
863 //
864 // Licensing:
865 //
866 // This code is distributed under the GNU LGPL license.
867 //
868 // Modified:
869 //
870 // 02 April 2010
871 //
872 // Author:
873 //
874 // John Burkardt
875 //
876 // Parameters:
877 //
878 // Input, int N, the order.
879 //
880 // Input, int DIM, the spatial dimension represented by this rule,
881 // in cases where a multidimensional product rule is being formed.
882 //
883 // Output, double W[N], the weights.
884 //
885 {
887  return;
888 }
889 
890 //****************************************************************************80
891 void SandiaRules2::jacobi_points ( int n, int dim, double x[] )
892 //****************************************************************************80
893 //
894 // Purpose:
895 //
896 // JACOBI_POINTS computes Jacobi quadrature points.
897 //
898 // Discussion:
899 //
900 // This function assumes the existence of a function:
901 // double parameter ( int dim, int offset )
902 // which can supply the values of the ALPHA and BETA parameters.
903 //
904 // Licensing:
905 //
906 // This code is distributed under the GNU LGPL license.
907 //
908 // Modified:
909 //
910 // 02 April 2010
911 //
912 // Author:
913 //
914 // John Burkardt
915 //
916 // Parameters:
917 //
918 // Input, int N, the order.
919 //
920 // Input, int DIM, the spatial dimension represented by this rule,
921 // in cases where a multidimensional product rule is being formed.
922 //
923 // Output, double X[N], the abscissas.
924 //
925 {
926  double alpha;
927  double beta;
928 
929  alpha = parameter ( dim, 0 );
930  beta = parameter ( dim, 1 );
931 
932  SandiaRules::jacobi_compute_points ( n, alpha, beta, x );
933 
934  return;
935 }
936 
937 //****************************************************************************80
938 void SandiaRules2::jacobi_weights ( int n, int dim, double w[] )
939 //****************************************************************************80
940 //
941 // Purpose:
942 //
943 // JACOBI_WEIGHTS computes Jacobi quadrature weights.
944 //
945 // Discussion:
946 //
947 // This function assumes the existence of a function:
948 // double parameter ( int dim, int offset )
949 // which can supply the values of the ALPHA and BETA parameters.
950 //
951 // Licensing:
952 //
953 // This code is distributed under the GNU LGPL license.
954 //
955 // Modified:
956 //
957 // 02 April 2010
958 //
959 // Author:
960 //
961 // John Burkardt
962 //
963 // Parameters:
964 //
965 // Input, int N, the order.
966 //
967 // Input, int DIM, the spatial dimension represented by this rule,
968 // in cases where a multidimensional product rule is being formed.
969 //
970 // Output, double W[N], the weights.
971 //
972 {
973  double alpha;
974  double beta;
975 
976  alpha = parameter ( dim, 0 );
977  beta = parameter ( dim, 1 );
978 
979  SandiaRules::jacobi_compute_weights ( n, alpha, beta, w );
980 
981  return;
982 }
983 
984 //****************************************************************************80
985 void SandiaRules2::laguerre_points ( int n, int dim, double x[] )
986 //****************************************************************************80
987 //
988 // Purpose:
989 //
990 // LAGUERRE_POINTS computes Laguerre quadrature points.
991 //
992 // Licensing:
993 //
994 // This code is distributed under the GNU LGPL license.
995 //
996 // Modified:
997 //
998 // 02 April 2010
999 //
1000 // Author:
1001 //
1002 // John Burkardt
1003 //
1004 // Parameters:
1005 //
1006 // Input, int N, the order.
1007 //
1008 // Input, int DIM, the spatial dimension represented by this rule,
1009 // in cases where a multidimensional product rule is being formed.
1010 //
1011 // Output, double X[N], the abscissas.
1012 //
1013 {
1014  laguerre_compute_points ( n, x );
1015  return;
1016 }
1017 
1018 //****************************************************************************80
1019 void SandiaRules2::laguerre_weights ( int n, int dim, double w[] )
1020 //****************************************************************************80
1021 //
1022 // Purpose:
1023 //
1024 // LAGUERRE_WEIGHTS computes Laguerre quadrature weights.
1025 //
1026 // Licensing:
1027 //
1028 // This code is distributed under the GNU LGPL license.
1029 //
1030 // Modified:
1031 //
1032 // 02 April 2010
1033 //
1034 // Author:
1035 //
1036 // John Burkardt
1037 //
1038 // Parameters:
1039 //
1040 // Input, int N, the order.
1041 //
1042 // Input, int DIM, the spatial dimension represented by this rule,
1043 // in cases where a multidimensional product rule is being formed.
1044 //
1045 // Output, double W[N], the weights.
1046 //
1047 {
1049  return;
1050 }
1051 
1052 //****************************************************************************80
1053 void SandiaRules2::legendre_points ( int n, int dim, double x[] )
1054 //****************************************************************************80
1055 //
1056 // Purpose:
1057 //
1058 // LEGENDRE_POINTS computes Legendre quadrature points.
1059 //
1060 // Licensing:
1061 //
1062 // This code is distributed under the GNU LGPL license.
1063 //
1064 // Modified:
1065 //
1066 // 02 April 2010
1067 //
1068 // Author:
1069 //
1070 // John Burkardt
1071 //
1072 // Parameters:
1073 //
1074 // Input, int N, the order.
1075 //
1076 // Input, int DIM, the spatial dimension represented by this rule,
1077 // in cases where a multidimensional product rule is being formed.
1078 //
1079 // Output, double X[N], the abscissas.
1080 //
1081 {
1083  return;
1084 }
1085 
1086 //****************************************************************************80
1087 void SandiaRules2::legendre_weights ( int n, int dim, double w[] )
1088 //****************************************************************************80
1089 //
1090 // Purpose:
1091 //
1092 // LEGENDRE_WEIGHTS computes Legendre quadrature weights.
1093 //
1094 // Licensing:
1095 //
1096 // This code is distributed under the GNU LGPL license.
1097 //
1098 // Modified:
1099 //
1100 // 02 April 2010
1101 //
1102 // Author:
1103 //
1104 // John Burkardt
1105 //
1106 // Parameters:
1107 //
1108 // Input, int N, the order.
1109 //
1110 // Input, int DIM, the spatial dimension represented by this rule,
1111 // in cases where a multidimensional product rule is being formed.
1112 //
1113 // Output, double W[N], the weights.
1114 //
1115 {
1117  return;
1118 }
1119 
1120 //****************************************************************************80
1121 void SandiaRules2::ncc_points ( int n, int dim, double x[] )
1122 //****************************************************************************80
1123 //
1124 // Purpose:
1125 //
1126 // NCC_POINTS computes Newton Cotes Closed quadrature points.
1127 //
1128 // Licensing:
1129 //
1130 // This code is distributed under the GNU LGPL license.
1131 //
1132 // Modified:
1133 //
1134 // 03 August 2011
1135 //
1136 // Author:
1137 //
1138 // John Burkardt
1139 //
1140 // Parameters:
1141 //
1142 // Input, int N, the order.
1143 //
1144 // Input, int DIM, the spatial dimension represented by this rule,
1145 // in cases where a multidimensional product rule is being formed.
1146 //
1147 // Output, double X[N], the abscissas.
1148 //
1149 {
1151  return;
1152 }
1153 
1154 //****************************************************************************80
1155 void SandiaRules2::ncc_weights ( int n, int dim, double w[] )
1156 //****************************************************************************80
1157 //
1158 // Purpose:
1159 //
1160 // NCC_WEIGHTS computes Newton Cotes Closed quadrature weights.
1161 //
1162 // Licensing:
1163 //
1164 // This code is distributed under the GNU LGPL license.
1165 //
1166 // Modified:
1167 //
1168 // 03 August 2011
1169 //
1170 // Author:
1171 //
1172 // John Burkardt
1173 //
1174 // Parameters:
1175 //
1176 // Input, int N, the order.
1177 //
1178 // Input, int DIM, the spatial dimension represented by this rule,
1179 // in cases where a multidimensional product rule is being formed.
1180 //
1181 // Output, double W[N], the weights.
1182 //
1183 {
1185 }
1186 
1187 //****************************************************************************80
1188 void SandiaRules2::nco_points ( int n, int dim, double x[] )
1189 //****************************************************************************80
1190 //
1191 // Purpose:
1192 //
1193 // NCO_POINTS computes Newton Cotes Open quadrature points.
1194 //
1195 // Licensing:
1196 //
1197 // This code is distributed under the GNU LGPL license.
1198 //
1199 // Modified:
1200 //
1201 // 03 August 2011
1202 //
1203 // Author:
1204 //
1205 // John Burkardt
1206 //
1207 // Parameters:
1208 //
1209 // Input, int N, the order.
1210 //
1211 // Input, int DIM, the spatial dimension represented by this rule,
1212 // in cases where a multidimensional product rule is being formed.
1213 //
1214 // Output, double X[N], the abscissas.
1215 //
1216 {
1218  return;
1219 }
1220 
1221 //****************************************************************************80
1222 void SandiaRules2::nco_weights ( int n, int dim, double w[] )
1223 //****************************************************************************80
1224 //
1225 // Purpose:
1226 //
1227 // NCO_WEIGHTS computes Newton Cotes Open quadrature weights.
1228 //
1229 // Licensing:
1230 //
1231 // This code is distributed under the GNU LGPL license.
1232 //
1233 // Modified:
1234 //
1235 // 03 August 2011
1236 //
1237 // Author:
1238 //
1239 // John Burkardt
1240 //
1241 // Parameters:
1242 //
1243 // Input, int N, the order.
1244 //
1245 // Input, int DIM, the spatial dimension represented by this rule,
1246 // in cases where a multidimensional product rule is being formed.
1247 //
1248 // Output, double W[N], the weights.
1249 //
1250 {
1252  return;
1253 }
1254 
1255 //****************************************************************************80
1256 void SandiaRules2::patterson_points ( int n, int dim, double x[] )
1257 //****************************************************************************80
1258 //
1259 // Purpose:
1260 //
1261 // PATTERSON_POINTS looks up Patterson quadrature points.
1262 //
1263 // Discussion:
1264 //
1265 // Our convention is that the abscissas are numbered from left to right.
1266 //
1267 // The rule is defined on [-1,1],
1268 //
1269 // Licensing:
1270 //
1271 // This code is distributed under the GNU LGPL license.
1272 //
1273 // Modified:
1274 //
1275 // 02 April 2010
1276 //
1277 // Author:
1278 //
1279 // John Burkardt
1280 //
1281 // Parameters:
1282 //
1283 // Input, int N, the order.
1284 //
1285 // Input, int DIM, the spatial dimension represented by this rule,
1286 // in cases where a multidimensional product rule is being formed.
1287 //
1288 // Output, double X[N], the abscissas.
1289 //
1290 {
1292  return;
1293 }
1294 
1295 //****************************************************************************80
1296 void SandiaRules2::patterson_weights ( int n, int dim, double w[] )
1297 //****************************************************************************80
1298 //
1299 // Purpose:
1300 //
1301 // PATTERSON_WEIGHTS looks up Patterson quadrature weights.
1302 //
1303 // Discussion:
1304 //
1305 // The allowed orders are 1, 3, 7, 15, 31, 63 and 127.
1306 //
1307 // The weights are positive, symmetric and should sum to 2.
1308 //
1309 // The user must preallocate space for the output array W.
1310 //
1311 // Licensing:
1312 //
1313 // This code is distributed under the GNU LGPL license.
1314 //
1315 // Modified:
1316 //
1317 // 02 April 2010
1318 //
1319 // Author:
1320 //
1321 // John Burkardt
1322 //
1323 // Parameters:
1324 //
1325 // Input, int N, the order.
1326 //
1327 // Input, int DIM, the spatial dimension represented by this rule,
1328 // in cases where a multidimensional product rule is being formed.
1329 //
1330 // Output, double W[N], the weights.
1331 //
1332 {
1334  return;
1335 }
1336 
1337 //****************************************************************************80
1338 double SandiaRules2::parameter ( int dim, int offset )
1339 //****************************************************************************80
1340 //
1341 // Purpose:
1342 //
1343 // PARAMETER is a user-supplied routine to retrieve parameters.
1344 //
1345 // Discussion:
1346 //
1347 // The declaration for this function is in SANDIA_RULES.H
1348 //
1349 // Licensing:
1350 //
1351 // This code is distributed under the GNU LGPL license.
1352 //
1353 // Modified:
1354 //
1355 // 02 April 2010
1356 //
1357 // Author:
1358 //
1359 // John Burkardt
1360 //
1361 // Parameters:
1362 //
1363 // Input, int DIM, the spatial dimension.
1364 //
1365 // Input, int OFFSET, the offset of the parameter within the
1366 // spatial dimension.
1367 //
1368 // Output, double PARAMETER, the value of the OFFSET-th parameter
1369 // associated with the DIM-th dimension.
1370 //
1371 {
1372  return 0.0;
1373 }
1374 
1375 } // namespace ROL
void hce_weights(int n, int dim, double w[])
void nco_compute_weights(int n, double w[])
void patterson_lookup_weights(int n, double w[])
void fejer2_points(int n, int dim, double x[])
void laguerre_compute_points(int order, double x[])
void nco_weights(int n, int dim, double w[])
void gen_laguerre_compute_points(int order, double alpha, double x[])
void hce_points(int n, int dim, double x[])
void gen_hermite_compute_weights(int order, double alpha, double w[])
void hcc_compute_points(int n, double x[])
void patterson_points(int n, int dim, double x[])
void nco_points(int n, int dim, double x[])
void jacobi_compute_points(int order, double alpha, double beta, double x[])
void ncc_weights(int n, int dim, double w[])
void ccn_points(int n, int dim, double x[])
void laguerre_weights(int n, int dim, double w[])
void hce_compute_points(int n, double x[])
void jacobi_compute_weights(int order, double alpha, double beta, double w[])
void hermite_compute_points(int order, double x[])
void jacobi_weights(int n, int dim, double w[])
double parameter(int dim, int offset)
void hermite_points(int n, int dim, double x[])
void hcc_weights(int n, int dim, double w[])
void legendre_compute_weights(int order, double w[])
void laguerre_points(int n, int dim, double x[])
void clenshaw_curtis_weights(int n, int dim, double w[])
void ccn_weights(int n, int dim, double w[])
void patterson_lookup_points(int n, double x[])
void laguerre_compute_weights(int order, double w[])
void gen_laguerre_compute_weights(int order, double alpha, double w[])
void gen_hermite_compute_points(int order, double alpha, double x[])
void clenshaw_curtis_points(int n, int dim, double x[])
void hermite_weights(int n, int dim, double w[])
void ccn_compute_points(int n, double x[])
void hcc_points(int n, int dim, double x[])
void hce_compute_weights(int n, double w[])
void gen_laguerre_points(int n, int dim, double x[])
void ncc_points(int n, int dim, double x[])
void hcc_compute_weights(int n, double w[])
void gen_laguerre_weights(int n, int dim, double w[])
void nco_compute_points(int n, double x[])
void hermite_genz_keister_lookup_points(int n, double x[])
void ncc_compute_weights(int n, double w[])
void patterson_weights(int n, int dim, double w[])
void gen_hermite_points(int n, int dim, double x[])
void gen_hermite_weights(int n, int dim, double w[])
void hermite_genz_keister_lookup_weights(int n, double w[])
void legendre_points(int n, int dim, double x[])
void hermite_genz_keister_points(int n, int dim, double x[])
void hermite_genz_keister_weights(int n, int dim, double w[])
void legendre_compute_points(int order, double x[])
void hermite_compute_weights(int order, double w[])
void ncc_compute_points(int n, double x[])
void jacobi_points(int n, int dim, double x[])
void ccn_compute_weights(int n, double w[])
void legendre_weights(int n, int dim, double w[])
void fejer2_weights(int n, int dim, double w[])