src/polynomial/misc/cl_UP_I_hermite.cc

   1 // hermite().
   2
   3 // General includes.
   4 #include "base/cl_sysdep.h"
   5
   6 // Specification.
   7 #include "cln/univpoly_integer.h"
   8
   9
  10 // Implementation.
  11
  12 #include "cln/integer.h"
  13
  14 namespace cln {
  15
  16 const cl_UP_I hermite (sintL n)
  17 {
  18 // The Hermite polynomials H_n(x) are defined as
  19 //
  20 //                             ( d  ) n
  21 //    H_n(x) = (-1)^n exp(x^2) (----)   exp(- x^2)
  22 //                             ( dx )
  23 //
  24 // They satisfy the recurrence relation
  25 //
  26 //    H_0(x) = 1
  27 //    H_{n+1}(x) = 2x H_n(x) - 2n H_{n-1}(x) for n >= 0.
  28 //
  29 // Theorem:
  30 //    H_n(x) satisfies the differential equation
  31 //    H_n''(x) - 2x*H_n'(x) + 2n*H_n(x) = 0.
  32 //
  33 // Proof: See elsewhere.
  34 //
  35 // Corollary:
  36 //    The coefficients c_{n,k} of H_n(x) = sum(k=0..n, c_{n,k} x^k)
  37 //    satisfy:
  38 //       c_{n,n} = 2^n,
  39 //       c_{n,n-1} = 0,
  40 //       c_{n,k} = (k+1)(k+2)/2(k-n)*c_{n,k+2}
  41 //
  42 // It follows that for n>=0
  43 //
  44 //    H_n(x) = sum(j=0..floor(n/2), (-1)^j n!/j!(n-2j)! 2^(n-2j) x^(n-2j))
  45 //
  46         var cl_univpoly_integer_ring R = find_univpoly_ring(cl_I_ring);
  47         var cl_UP_I h = R->create(n);
  48         var sintL k = n;
  49         var cl_I c_k = ash(1,n);
  50         for (;;) {
  51                 h.set_coeff(k,c_k);
  52                 k = k-2;
  53                 if (k < 0)
  54                         break;
  55                 c_k = exquo((cl_I)(k+1) * (cl_I)(k+2) * c_k,
  56                             2*(cl_I)(k-n));
  57         }
  58         h.finalize();
  59         return h;
  60 }
  61
  62 }  // namespace cln