11#ifndef _PFASST__QUADRATURE__GAUSS_RADAU_HPP_
22#define _PFASST__QUADRATURE__GAUSS_RADAU_HPP_
33
4- #include < cassert>
5- #include < vector>
6-
7- #include " ../interfaces.hpp"
8- #include " polynomial.hpp"
94#include " interface.hpp"
105
11- using namespace std ;
12-
136
147namespace pfasst
158{
@@ -27,47 +20,24 @@ namespace pfasst
2720
2821 public:
2922 // ! @{
30- explicit GaussRadau (const size_t num_nodes)
31- : IQuadrature<precision>(num_nodes)
32- {
33- if (this ->num_nodes < 2 ) {
34- throw invalid_argument (" Gauss-Radau quadrature requires at least two quadrature nodes." );
35- }
36- this ->compute_nodes ();
37- this ->compute_weights ();
38- }
39-
23+ explicit GaussRadau (const size_t num_nodes);
4024 GaussRadau () = default ;
41-
4225 virtual ~GaussRadau () = default ;
4326 // ! @}
4427
4528 // ! @{
46- virtual bool left_is_node () const { return LEFT_IS_NODE; }
47-
48- virtual bool right_is_node () const { return RIGHT_IS_NODE; }
29+ virtual bool left_is_node () const override ;
30+ virtual bool right_is_node () const override ;
4931 // ! @}
5032
5133 protected:
5234 // ! @{
53- virtual void compute_nodes () override
54- {
55- this ->nodes = vector<precision>(this ->num_nodes , precision (0.0 ));
56- auto l = Polynomial<precision>::legendre (this ->num_nodes );
57- auto lm1 = Polynomial<precision>::legendre (this ->num_nodes - 1 );
58-
59- for (size_t i = 0 ; i < this ->num_nodes ; i++) {
60- l[i] += lm1[i];
61- }
62- auto roots = l.roots ();
63- for (size_t j = 1 ; j < this ->num_nodes ; j++) {
64- this ->nodes [j - 1 ] = 0.5 * (1.0 - roots[this ->num_nodes - j]);
65- }
66- this ->nodes .back () = 1.0 ;
67- }
35+ virtual void compute_nodes () override ;
6836 // ! @}
6937 };
7038 } // ::pfasst::quadrature
7139} // ::pfasst
7240
41+ #include " gauss_radau_impl.hpp"
42+
7343#endif // _PFASST__QUADRATURE__GAUSS_RADAU_HPP_
0 commit comments