hiro/arb
1
0
Fork 0
mirror of https://github.com/vale981/arb synced 2025-03-05 17:31:38 -05:00
Code Issues Projects Releases Packages Wiki Activity Actions

TST: lower incomplete gamma series with nonpositive integer parameter

Browse source
This commit is contained in:
alex 2016-04-26 21:45:41 -04:00
parent 70e43f9da9
commit abcaa05e61
2 changed files with 52 additions and 23 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

3
acb_hypgeom/gamma_lower_series.c
View file

@ -99,8 +99,7 @@ acb_hypgeom_gamma_lower_series(acb_poly_t g, const acb_t s,
arf_sgn(arb_midref(acb_realref(s))) <= 0 &&
arf_cmpabs_2exp_si(arb_midref(acb_realref(s)), 30) < 0)
{
/* fixme: deliberate typo, add negative sign */
slong exp = arf_get_si(arb_midref(acb_realref(s)), ARF_RND_DOWN);
slong exp = -arf_get_si(arb_midref(acb_realref(s)), ARF_RND_DOWN);
acb_poly_pow_ui_trunc_binexp(g, h, (ulong) exp, n, prec);
return;
}

72
acb_hypgeom/test/t-gamma_lower_series.c
View file

@ -64,7 +64,16 @@ int main()
acb_poly_randtest(S, state, m, bits1, 3);
acb_poly_randtest(A, state, m, bits1, 3);
acb_poly_randtest(B, state, m, bits1, 3);
acb_randtest(s, state, bits1, 3);
if (n_randint(state, 10))
{
acb_randtest(s, state, bits1, 3);
}
else
{
slong p = n_randint(state, 100);
acb_set_d(s, -p);
}
acb_hypgeom_gamma_lower_series(A, s, S, regularized, n1, bits2);
acb_hypgeom_gamma_lower_series(B, s, S, regularized, n2, bits3);
@ -112,33 +121,54 @@ int main()
}
else if (regularized == 2)
{
/* (h(x)^-s integral(f(h(x)) h'(x)) / c)' =
* h(x)^-(s+1) (h(x) f(h(x)) - s integral(f(h(x)) h'(x))) h'(x) / c */
acb_poly_t D;
acb_poly_init(D);
if (acb_is_int(s) &&
arf_sgn(arb_midref(acb_realref(s))) <= 0 &&
arf_cmpabs_2exp_si(arb_midref(acb_realref(s)), 30) < 0)
{
/* (h(x)^-s)' = -s h(x)^-(s+1) h'(x) */
acb_poly_derivative(C, S, bits2);
acb_poly_derivative(B, S, bits2);
acb_poly_mullow(D, C, B, n1, bits2);
acb_poly_integral(D, D, bits2);
_acb_vec_scalar_mul(D->coeffs, D->coeffs, D->length, s, bits2);
acb_poly_mullow(C, C, S, n1, bits2);
acb_poly_sub(D, C, D, bits2);
acb_add_ui(t, s, 1, bits2);
acb_neg(t, t);
acb_poly_pow_acb_series(B, S, t, n1, bits2);
acb_add_ui(t, s, 1, bits2);
acb_neg(t, t);
acb_poly_pow_acb_series(B, S, t, n1, bits2);
acb_poly_mullow(C, C, B, n1, bits2);
acb_poly_mullow(C, D, B, n1, bits2);
_acb_vec_scalar_mul(C->coeffs, C->coeffs, C->length, s, bits2);
acb_poly_neg(C, C);
acb_poly_derivative(B, S, bits2);
acb_poly_mullow(C, C, B, n1, bits2);
acb_poly_truncate(C, n1 - 1);
}
else
{
/* (h(x)^-s integral(f(h(x)) h'(x)) / c)' =
* h(x)^-(s+1) (h(x) f(h(x)) - s integral(f(h(x)) h'(x))) h'(x) / c */
acb_poly_t D;
acb_poly_init(D);
acb_gamma(c, s, bits2);
_acb_vec_scalar_div(C->coeffs, C->coeffs, C->length, c, bits2);
acb_poly_derivative(B, S, bits2);
acb_poly_mullow(D, C, B, n1, bits2);
acb_poly_integral(D, D, bits2);
_acb_vec_scalar_mul(D->coeffs, D->coeffs, D->length, s, bits2);
acb_poly_mullow(C, C, S, n1, bits2);
acb_poly_sub(D, C, D, bits2);
acb_poly_truncate(C, n1 - 1);
acb_add_ui(t, s, 1, bits2);
acb_neg(t, t);
acb_poly_pow_acb_series(B, S, t, n1, bits2);
acb_poly_clear(D);
acb_poly_mullow(C, D, B, n1, bits2);
acb_poly_derivative(B, S, bits2);
acb_poly_mullow(C, C, B, n1, bits2);
acb_gamma(c, s, bits2);
_acb_vec_scalar_div(C->coeffs, C->coeffs, C->length, c, bits2);
acb_poly_truncate(C, n1 - 1);
acb_poly_clear(D);
}
}
acb_poly_derivative(B, A, bits2);