RJF asks me to forward this, since it bounced for him, evidently.
---------- Forwarded message ---------- From: Richard Fateman <fate...@cs.berkeley.edu> Date: Mon, Apr 26, 2010 at 10:15 PM Subject: Re: [sage-devel] numerically stable fast univariate polynomial multiplication over RR[x] To: William Stein <wst...@gmail.com> Cc: sage-devel@googlegroups.com, rjf <fate...@gmail.com> William Stein wrote: > > On Mon, Apr 26, 2010 at 8:00 PM, Jason Grout > <jason-s...@creativetrax.com> wrote: >> >> When multiplying two univariate polynomials over RR[x], Sage uses the naive >> O(n^2) algorithm. The docstring for _mul_ cites numerical stability as the >> reason to not use asymptotically faster algorithms: >> >> "Here we use the naive `O(n^2)` algorithm, as asymptotically faster >> algorithms such as Karatsuba can have very inaccurate results due to >> intermediate rounding errors. " If you are doing arithmetic over the reals (RR?) then there is no rounding error. But you apparently are doing arithmetic in some approximation to the reals-- presumably doing n-bit floating-point arithmetic for some particular n. If you want to get the right answer, then convert all the numbers to exact rationals in QQ. Multiply and then convert back. Or multiply all the coefficients by 2^k so they are all integers. Do the arithmetic in ZZ. You can also multiply (approximately) by a fast Fourier transform convolution, which is reputedly rather stable. The stability may in fact directly address the notion of accuracy, at least if you are talking about dense polynomials. You can, I think, call some routine in MPFR. >> >> (this is followed by an example where Karatsuba multiplication is shown to >> mess up compared to naive multiplication) How are you judging whether one polynomial is "messed up" relative to another? I am unable to find the example or the documentation by searching from the Sage home page for f.mul_fateman or f.mul or karatsuba. None of which come up with anything. So I don't know what examples you are talking about. >> >> Can we do any better? Does anyone know of any algorithms for numerically >> stable fast univariate polynomial multiplication (better than the naive >> algorithm)? > > 1) I don't know what I'm talking about, but the first thing that pops > into my mind (in 1 second, literally) is to simply bump up the > precision some heuristic amount, then do the arithmetic using some > asymptotically fast algorithm using interval arithmetic, and see what > the precision of the result is. How do you define the precision of a polynomial? The minimum precision of any coefficient? Consider the polynomial. 10^200*x+ [-1,1] ? Is this unsatisfactory? For |x|>1, you have 200+ decimal digits. Pretty accurate. For |x|<10^(-200), you have no decimal digits. Pretty inaccurate. Using intervals you would get some probably horrendous over-estimate, on a coefficient by coefficient basis. > If it is sufficient, return the > non-interval answer. If it is not sufficient, increase the precision > and try again. By how much? Say k bits more.. > I think the complexity of arithmetic with intervals is > only a constant factor (2? 4?) compared to not using intervals, so the > complexity of what I just suggested should usually have nearly the > same O complexity as some asymptotically fast algorithm. you are assuming that arithmetic with n-bit floating-point numbers takes the same time as arithmetic with n+k bits, which is clearly wrong. and if the intervals are exponentially wider than necessary (frequently the case, 2^p for p operations) you are going to lose big. That's off the top of my head. > > 2) There is something called Fateman multiplication, which might be > relevant. I don't know what you are calling Fateman multiplication, so I can't comment. Maybe you could email me whatever documentation you are referring to. I do NOT have Sage installed. > > Didier Deshommes put it in Sage for some reason like 4 > years ago. Type f._mul_fateman? for more. Also, I've cc'd Richard > Fateman himself, in hopes you can say something about your question. > > -- William I hope this has given you something to think about. RJF -- William Stein Professor of Mathematics University of Washington http://wstein.org -- To post to this group, send an email to sage-devel@googlegroups.com To unsubscribe from this group, send an email to sage-devel+unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/sage-devel URL: http://www.sagemath.org
