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

use arb_dot in arb_poly division

Browse source
This commit is contained in:
fredrik 2018-08-20 22:43:16 +02:00
parent 10ebab1b70
commit 1764c21c03
2 changed files with 6 additions and 19 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

14
arb_poly/div_series.c
View file

@ -44,7 +44,7 @@ _arb_poly_div_series(arb_ptr Q, arb_srcptr A, slong Alen,
{ {
/* The basecase algorithm is faster for much larger Blen and n than /* The basecase algorithm is faster for much larger Blen and n than
this, but unfortunately has worse numerical stability. */ this, but unfortunately has worse numerical stability. */
slong i, j; slong i;
arb_t q; arb_t q;
arb_init(q); arb_init(q);
@ -54,16 +54,8 @@ _arb_poly_div_series(arb_ptr Q, arb_srcptr A, slong Alen,
for (i = 1; i < n; i++) for (i = 1; i < n; i++)
{ {
arb_mul(Q + i, B + 1, Q + i - 1, prec); arb_dot(Q + i, (i < Alen) ? A + i : NULL, 1,
B + 1, 1, Q + i - 1, -1, FLINT_MIN(i, Blen - 1), prec);
for (j = 2; j < FLINT_MIN(i + 1, Blen); j++)
arb_addmul(Q + i, B + j, Q + i - j, prec);
if (i < Alen)
arb_sub(Q + i, A + i, Q + i, prec);
else
arb_neg(Q + i, Q + i);
if (!arb_is_one(q)) if (!arb_is_one(q))
arb_mul(Q + i, Q + i, q, prec); arb_mul(Q + i, Q + i, q, prec);
} }

11
arb_poly/inv_series.c
View file

@ -37,7 +37,7 @@ _arb_poly_inv_series(arb_ptr Qinv,
} }
else else
{ {
slong i, j, blen; slong i, blen;
/* The basecase algorithm is faster for much larger Qlen or len than /* The basecase algorithm is faster for much larger Qlen or len than
this, but unfortunately also much less numerically stable. */ this, but unfortunately also much less numerically stable. */
@ -48,15 +48,10 @@ _arb_poly_inv_series(arb_ptr Qinv,
for (i = 1; i < blen; i++) for (i = 1; i < blen; i++)
{ {
arb_mul(Qinv + i, Q + 1, Qinv + i - 1, prec); arb_dot(Qinv + i, NULL, 1,
Q + 1, 1, Qinv + i - 1, -1, FLINT_MIN(i, Qlen - 1), prec);
for (j = 2; j < FLINT_MIN(i + 1, Qlen); j++)
arb_addmul(Qinv + i, Q + j, Qinv + i - j, prec);
if (!arb_is_one(Qinv)) if (!arb_is_one(Qinv))
arb_mul(Qinv + i, Qinv + i, Qinv, prec); arb_mul(Qinv + i, Qinv + i, Qinv, prec);
arb_neg(Qinv + i, Qinv + i);
} }
if (len > blen) if (len > blen)