On Mon, 18 Aug 2008, John Cremona wrote: > > 2008/8/18 Bill Hart <[EMAIL PROTECTED]>: >> >> Oh sorry, I get you. You mean that if the first term *is* a square >> then the power series method will *always* find a purported square >> root for you. > > Yes, that's what I meant. > >> >> So what I meant is that as soon as the power series algorithm returns >> a non-zero coefficient past the n/2-th coefficient you know your >> original polynomial was not a square. This may occur if your Newton >> iteration exceeds this many coefficients of the power series square >> root. So I think what I said makes sense. >> > > That's true.
The power series sqrt code alreay checks to see if an exact square was found (note that the result is not +O(x^20)). It may not abort as early as it could though. > >> Of course to prove it is a polynomial square root, if it is, you need >> to do an n/2*n/2 truncated squaring. I think that is less than what is >> required to complete the final step of the Newton iteration. >> > > John > >> Bill. >> >> On 18 Aug, 20:39, Bill Hart <[EMAIL PROTECTED]> wrote: >>> If a polynomial is a perfect square, is it not necessarily true that >>> the constant term is a square? Thus you already bail out if this is >>> not the case. >>> >>> I wonder when the time will come that virtually everything that can be >>> implemented, is implemented in Sage. The the algorithm will then >>> always be, "let's ask Sage.... nope, an algorithm does not exist!" >>> >>> Bill. >>> >>> On 18 Aug, 18:32, "John Cremona" <[EMAIL PROTECTED]> wrote: >>> >>> >>> >>>> The power series method would apply to any polynomial whose constant >>>> term is a square (or rather, whose lowest degree monomial has the form >>>> (square coefficient)*(even power of x). Viewed as a power series such >>>> a thing is always a square, as a simple induction constructs a square >>>> root. We already have this: >>> >>>> sage: R.<x>=PowerSeriesRing(QQ) >>>> sage: f=R(1+2*x+3*x^2) >>>> sage: f.sqrt() >>>> 1 + x + x^2 - x^3 + 1/2*x^4 + 1/2*x^5 - 3/2*x^6 + 3/2*x^7 + 3/8*x^8 - >>>> 29/8*x^9 + 43/8*x^10 - 11/8*x^11 - 155/16*x^12 + 325/16*x^13 - >>>> 215/16*x^14 - 393/16*x^15 + 9891/128*x^16 - 10269/128*x^17 - >>>> 5901/128*x^18 + 36813/128*x^19 + O(x^20) >>>> sage: g=f^2 >>>> sage: g.sqrt() >>>> 1 + 2*x + 3*x^2 >>>> sage: PolynomialRing(QQ,'X')(g.sqrt().list()) >>>> 3*X^2 + 2*X + 1 >>> >>>> John >>> >>>> 2008/8/18 Bill Hart <[EMAIL PROTECTED]>: >>> >>>>> On 18 Aug, 17:39, "Joel B. Mohler" <[EMAIL PROTECTED]> wrote: >>>>>> On Monday 18 August 2008 11:44:11 am Bill Hart wrote: >>> >>>>>>> Assuming FLINT is in fact being used for the GCD in Z[x], the >>>>>>> implementation is only "fast" up to degree about 250. It depends on >>>>>>> your definition of fast though. :-) In a later release of FLINT, GCD >>>>>>> will be asymptotically faster, and certainly by degree 1000 will be >>>>>>> many times faster than it currently is. >>> >>>>>> By "fast", I meant faster than my naive implementation of the sqrt. The >>>>>> ZZ >>>>>> gcd certainly seems fast enough for any problem I'm likely to have in the >>>>>> next 6 months :). The sage I'm using (3.0.6) is using the flint >>>>>> fmpz_poly_gcd. >>> >>>>> OK. That is currently quadratic, but with a low constant. :-) >>> >>>>>>> If I'm not mistaken, the linear algebra method is nearly cubic >>>>>>> complexity. The GCD method should be subquadratic (though only just). >>>>>>> The method using power series has the same complexity as Karatsuba >>>>>>> multiplication of polynomials, in terms of ring operations (see papers >>>>>>> of Paul Zimmermann), and the constant can be made very low by use of >>>>>>> the recursive middle product algorithm. >>> >>>>>> My algorithm is most definitely quadratic (as I had already thought from >>>>>> the >>>>>> simple nested for-loops). >>> >>>>> Yes I've had a look now, and I think you are right. >>> >>>>>> Doubling the degree produces a nearly perfect >>>>>> quadrupling of run-time. Note that it isn't a linear system of the >>>>>> coefficients and the solution is so obvious that it's easy to find one >>>>>> unknown at a time. I'm not sure what you are calling the "power series >>>>>> method" -- I wonder if it may be what I'm doing. >>> >>>>> No it definitely can't be, since the latter has much better asymptotic >>>>> complexity. I'm simply referring to thinking of the polynomial as a >>>>> power series, then essentially using Newton iteration to compute the >>>>> square root, though there are many variants of this basic idea which >>>>> can cut the implied constant by a very considerable factor. >>> >>>>> Another method is to take the logarithm of the power series, multiply >>>>> by 0.5 and then exponentiate. >>> >>>>>>> Note that if one requires the square root of an exact square of >>>>>>> polynomials, if one uses the power series method, one only needs to >>>>>>> work to a precision about half the length of the original polynomial >>>>>>> since that is how large the square root will be. I think so anyway, >>>>>>> unless I'm missing something silly here (altogether likely). >>> >>>>>> That makes sense. If you know that you have a perfect square you only >>>>>> need to >>>>>> think about half the coefficients of the square input. Of course, if you >>>>>> want to raise an error when something is not a perfect square, you'll >>>>>> need to >>>>>> test all the coefficients by comparing with the square of your tentative >>>>>> square root. >>> >>>>> It depends on how you implement it. As soon as you know one of the >>>>> coefficients is wrong for it to be a square, you can bail out and say >>>>> it isn't. It may not be necessary to jack the precision right up to n >>>>> for that. E.g. suppose n is just below a power of 2 and you use Newton >>>>> iteration to double the number of coefficients each time, you may >>>>> already know by the time you have just passed half way that you don't >>>>> have a square. For that matter, you may know at any time. >>> >>>>> It should also be possible to use a multimodular method to check it is >>>>> a square root in time something like n log (nm) where n is the degree >>>>> and m is the largest coefficient. This has a significant chance of >>>>> bailing out early if it isn't a square root and would certainly be >>>>> faster than multiplying out in full. I suspect on average O(1) primes >>>>> are needed on average to check if it isn't a square, and so this gives >>>>> an asymptotic improvement over multimodular multiplication, which >>>>> takes log m primes. >>> >>>>> Obviously over a general ring, you have to multiply using classical or >>>>> Karatsuba or multiplication, not all of which needs to be done to get >>>>> some of the coefficients. Of course that is a pointless observation in >>>>> special rings like Z, Z/pZ and Q where you don't have to use these >>>>> algorithms. >>> >>>>> And also note you only need a truncated product, since you know >>>>> already the top half of the terms are correct. Asymptotically that >>>>> doesn't save anything, but FLINT saves a small factor doing truncated >>>>> products anyway, since there is less real work to do, even when it >>>>> uses an FFT. You can combine this suggestion with the multimodular >>>>> suggestion. >>> >>>>> Another throw away observation is that squaring can be done faster >>>>> than a full product. FLINT does this automatically for sufficiently >>>>> large polynomials, i.e. length > 20 (FLINT 2.0 will do it for all >>>>> polynomials and I've already written some of the code for that). Again >>>>> you can use this with the multimodular and truncated product >>>>> suggestions. >>> >>>>> Anyhow, I've just outlined why I won't be doing this in FLINT >>>>> soon. :-) There's actually quite a lot of tricks to go through, and >>>>> there may even be an algorithm for detecting perfect squares which I >>>>> know nothing about. >>> >>>>> Bill. >>> >> > > > > --~--~---------~--~----~------------~-------~--~----~ To post to this group, send email to [email protected] To unsubscribe from this group, send email to [EMAIL PROTECTED] For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://www.sagemath.org -~----------~----~----~----~------~----~------~--~---
