2009/6/21 William Stein <wst...@gmail.com>: > > On Sun, Jun 21, 2009 at 8:38 PM, gsw<georgswe...@googlemail.com> wrote: >> >> Hi John, >> >> On 21 Jun., 17:47, John Cremona <john.crem...@gmail.com> wrote: >>> This should be of interest to anyone who has ever had to manage >>> precision issues between Sage and pari real and complex numbers (e.g. >>> Alex Ghitza). Others can move on. >>> >>> In the conversion of a pari complex number back to Sage (in >>> sage/libs/pari/gen_py.py in the function python(z)), the precision of >>> the ComplexField in which the returned value should lie is taken to be >>> the *minimum* of the precisions of the real and imaginary parts of the >>> pari complex. Here I am assuming that both real and imaginary parts >>> are inexact, otherwise there is no issue; and note that a pari gen >>> object (type t_COMPLEX) can have different types, and also different >>> precisions, for the real and imaginary parts. And we *are* converting >>> properly between pari's precision measured in words and Sage's >>> measured in bits. >>> >>> I propose to change this to the *maximum* of the two precisions. >>> >>> I have been experimenting with the wrapping of the pari functions >>> ellwp0(), which essentially takes ((w1,w2),z) where (w1,w2) is a pair >>> of complex numbers defining a lattice L in C and z represents and >>> element of C/L, and returns the pair (P(z),P'(z)) where P() is the >>> Weierstrass P function. In some cases I find that the real and >>> imaginary parts returned by pari's function can have very different >>> precisions. >>> >>> Example (using code now developed, so this will not work in the >>> released version of Sage): >>> >>> sage: E = EllipticCurve([1,1,1,-8,6]) >>> sage: P = E(0,2) >>> sage: L = E.period_lattice() >>> sage: z = L.elliptic_logarithm(P, prec=129) >>> >>> Using the current min of the precisions (the exact value should be (0,2)): >>> sage: L.elliptic_exponential(z) >>> (-3.8805107275644938151041666666176877354e-11 - >>> 2.7369110631344083416479093423631174439e-48*I, >>> 2.0000000000194025536378224690755208333 + >>> 4.1053665947016125124718640135446761659e-48*I) >>> >>> Using the proposed max: >>> sage: L.elliptic_exponential(z) >>> (8.8162076311671563097655240291668425836e-39 - >>> 2.7369110631344083416479093423631174439e-48*I, >>> 2.0000000000000000000000000000000000000 + >>> 4.1053665947016125124718640135446761659e-48*I) >>> >>> I swear that the only change between those two is changing line 166 of >>> sage/libs/pari/gen_py.py from >>> prec = >>> min(gen.prec_words_to_bits(xprec),gen.prec_words_to_bits(yprec)) >>> to >>> prec = >>> max(gen.prec_words_to_bits(xprec),gen.prec_words_to_bits(yprec)) >>> >>> Note that in this example I started with 129 bits precision; at 128 bits >>> we get >>> >>> sage: z = L.elliptic_logarithm(P, prec=128) >>> sage: L.elliptic_exponential(z) >>> (8.6692708373143703712694319620140618739e-38, >>> 1.9999999999999999999999999999999999999) >>> (using min) >>> and the identical >>> sage: L.elliptic_exponential(z) >>> (8.6692708373143703712694319620140618739e-38, >>> 1.9999999999999999999999999999999999999) >>> (using max) >>> which are the same. Here 128 bits is an exact number of words, while >>> 129 forces pari to use an extra word most of which is garbage on >>> input, and somehow the effect is to catastrophically reduce the >>> precision of the output! >>> >>> I just checked (using testall) that this change has no effect at all >>> on any doctests, so I propose to include it in my upcoming patch >>> (which implements the complex elliptic exponential for elliptic >>> curves, by a simple call to the already wrapped pari.ellwp0(), but >>> which has taken ages to sort out the precision properly. >>> >>> John >> >> let me try to rephrase this, I'm not sure whether I did understand you >> correctly: >> >> a) On certain inputs a certain pari function returns complex values, >> whose real parts and imaginary parts it reports as being quite >> different (in terms of precision).
True. >> >> b) Closer inspection yields that the part with the lower precision >> reported actually holds more information than the reported precision >> would allow for, especially it is OK to regard this part to be of the >> (greater) precision of the other part. >> Hmm. In the case in question the imaginary part should have been 0 exactly if the inputs were known to infinite precision. But the value returned was approximately 0 with a precision of only 3 words while the real part had precision 7 words. >> c) If so, I would consider this behaviour a bug in pari, to be fixed >> in pari. I was planning to report this to pari, but it is not easy to give a reproducible example --precision is handled differently in pari's library mode than in gp, so it would mean writing a C program. I'm not sure how the pari people would like getting a bug report which relied on them installing Sage to verify! It does not help that the pari documentation is less than clear about what precision should be expected. (There is a prec parameter to ellwp0 but that parameter is frequently ignored in cases where the precision of the other parameters can be used to determine the best output precision). >> >> d) In the meantime, you propose as a workaround to alter the >> conversion behaviour from pari to Sage for complex numbers in general. >> >> If I got right --- I vote +1 for your proposal. >> Thanks, but I'll take William's suggestion on board too. > > If he got it right above, then I would suggest instead making a > function in gen.pyx that does what you want and using it directly. > You could also add a non-default flag to the conversion function so > that when that option is set the conversion function uses the max > instead of the min. I'll try that and see how it goes. thanks for the feedback, John > > William > > > > --~--~---------~--~----~------------~-------~--~----~ To post to this group, send email to sage-devel@googlegroups.com To unsubscribe from this group, send email to sage-devel-unsubscr...@googlegroups.com For more options, visit this group at http://groups.google.com/group/sage-devel URLs: http://www.sagemath.org -~----------~----~----~----~------~----~------~--~---
