1 #include "rheolef/compiler.h"
22 #include "rheolef/gamma.h"
23 #include "rheolef/jacobi.h"
24 #include "rheolef/jacobi_roots.h"
27 template <
class Size,
class OutputIterator1,
class OutputIterator2>
29 typename std::iterator_traits<OutputIterator1>::value_type
alpha,
30 typename std::iterator_traits<OutputIterator1>::value_type
beta,
31 OutputIterator1 zeta, OutputIterator2 omega)
33 typedef typename std::iterator_traits<OutputIterator1>::value_type
T;
34 check_macro (R >= 1,
"gauss_radau_jacobi: node number " << R <<
" may be >= 1");
39 w0 *= 1/(
T(1.*R)*(
alpha+
T(1.*R)));
40 for (Size k = 1; k <= size_t(static_cast<int>(
beta)); k++) {
41 num *=
T(1.*R+k)/(
alpha+
T(1.*R+k));
42 w0 *=
T(sqr(
T(
int(k))))/(
T(1.*R+k)*(
alpha+R+k));
55 for (Size r = 1; r < R; ++r)
56 omega[r] = num*(1-zeta[r])/sqr(P1(zeta[r]));
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
Float alpha[pmax+1][pmax+1]
This file is part of Rheolef.
void gauss_radau_jacobi(Size R, typename std::iterator_traits< OutputIterator1 >::value_type alpha, typename std::iterator_traits< OutputIterator1 >::value_type beta, OutputIterator1 zeta, OutputIterator2 omega)
space_mult_list< T, M > pow(const space_basic< T, M > &X, size_t n)
void jacobi_roots(Size R, T alpha, T beta, OutputIterator zeta)