1 /*
2 Copyright (C) 2015 Fredrik Johansson
3
4 This file is part of Arb.
5
6 Arb is free software: you can redistribute it and/or modify it under
7 the terms of the GNU Lesser General Public License (LGPL) as published
8 by the Free Software Foundation; either version 2.1 of the License, or
9 (at your option) any later version. See <http://www.gnu.org/licenses/>.
10 */
11
12 #include "arb_poly.h"
13 #include "arb_hypgeom.h"
14
15 void
_arb_hypgeom_gamma_upper_series(arb_ptr g,const arb_t s,arb_srcptr h,slong hlen,int regularized,slong n,slong prec)16 _arb_hypgeom_gamma_upper_series(arb_ptr g, const arb_t s, arb_srcptr h, slong hlen, int regularized, slong n, slong prec)
17 {
18 arb_t c;
19 arb_init(c);
20
21 arb_hypgeom_gamma_upper(c, s, h, regularized, prec);
22
23 hlen = FLINT_MIN(hlen, n);
24
25 if (hlen == 1)
26 {
27 _arb_vec_zero(g + 1, n - 1);
28 }
29 else
30 {
31 arb_ptr t, u, v;
32 arb_ptr w = NULL;
33
34 t = _arb_vec_init(n);
35 u = _arb_vec_init(n);
36 v = _arb_vec_init(n);
37
38 if (regularized == 2)
39 {
40 w = _arb_vec_init(n);
41 arb_neg(t, s);
42 _arb_poly_pow_arb_series(w, h, hlen, t, n, prec);
43 }
44
45 /* Gamma(s, h(x)) = -integral(h'(x) h(x)^(s-1) exp(-h(x)) */
46 arb_sub_ui(u, s, 1, prec);
47 _arb_poly_pow_arb_series(t, h, hlen, u, n, prec);
48 _arb_poly_derivative(u, h, hlen, prec);
49 _arb_poly_mullow(v, t, n, u, hlen - 1, n, prec);
50 _arb_vec_neg(t, h, hlen);
51 _arb_poly_exp_series(t, t, hlen, n, prec);
52 _arb_poly_mullow(g, v, n, t, n, n, prec);
53 _arb_poly_integral(g, g, n, prec);
54 _arb_vec_neg(g, g, n);
55
56 if (regularized == 1)
57 {
58 arb_rgamma(t, s, prec);
59 _arb_vec_scalar_mul(g, g, n, t, prec);
60 }
61 else if (regularized == 2)
62 {
63 _arb_vec_set(u, g, n);
64 _arb_poly_mullow(g, u, n, w, n, n, prec);
65 _arb_vec_clear(w, n);
66 }
67
68 _arb_vec_clear(t, n);
69 _arb_vec_clear(u, n);
70 _arb_vec_clear(v, n);
71 }
72
73 arb_swap(g, c);
74 arb_clear(c);
75 }
76
77 void
arb_hypgeom_gamma_upper_series(arb_poly_t g,const arb_t s,const arb_poly_t h,int regularized,slong n,slong prec)78 arb_hypgeom_gamma_upper_series(arb_poly_t g, const arb_t s,
79 const arb_poly_t h, int regularized, slong n, slong prec)
80 {
81 slong hlen = h->length;
82
83 if (n == 0)
84 {
85 arb_poly_zero(g);
86 return;
87 }
88
89 arb_poly_fit_length(g, n);
90
91 if (hlen == 0)
92 {
93 arb_t t;
94 arb_init(t);
95 _arb_hypgeom_gamma_upper_series(g->coeffs, s, t, 1, regularized, n, prec);
96 arb_clear(t);
97 }
98 else
99 {
100 _arb_hypgeom_gamma_upper_series(
101 g->coeffs, s, h->coeffs, hlen, regularized, n, prec);
102 }
103
104 _arb_poly_set_length(g, n);
105 _arb_poly_normalise(g);
106 }
107