Opened 6 years ago
Closed 6 years ago
#15483 closed defect (fixed)
Simon 2descent gives RuntimeError for an elliptic curve over a quadratic field
Reported by:  pbruin  Owned by:  

Priority:  major  Milestone:  sage6.2 
Component:  elliptic curves  Keywords:  simon_two_descent 
Cc:  cremona  Merged in:  
Authors:  Peter Bruin  Reviewers:  Jeroen Demeyer 
Report Upstream:  Workaround found; Bug reported upstream.  Work issues:  
Branch:  b662aba (Commits)  Commit:  b662aba98ad4122a82fc135cf9b5dbf58de2ba0e 
Dependencies:  #11005  Stopgaps: 
Description (last modified by )
[See #15608 for a list of open simon_two_descent tickets]
The following tries to do 2descent on an elliptic curve over a quadratic field:
K.<s>=QuadraticField(229) c4=2173235*(1s)/2 c6=124369+15988*(1s)/2 E=EllipticCurve([c4/48,c6/864]) E.simon_two_descent()
It fails with a RuntimeError
:
Traceback (most recent call last): ... RuntimeError: *** at toplevel: ans=bnfellrank(K,[0,0,0, *** ^ *** in function bnfellrank: ...eqtheta,rnfeq,bbnf];rang= *** bnfell2descent_gen(b *** ^ *** in function bnfell2descent_gen: ...riv,r=nfsqrt(nf,norm(zc)) *** [1];if(DEBUGLEVEL_el *** ^ *** array index (1) out of allowed range [none]. An error occurred while running Simon's 2descent program
The problem is apparently that the square root of norm(zc)
(which is indeed a square) is not computed due to a precision bug; see comment:7 below.
Change History (23)
comment:1 followup: ↓ 2 Changed 6 years ago by
comment:2 in reply to: ↑ 1 Changed 6 years ago by
 Description modified (diff)
Replying to cremona:
There are many similar bug reports with separate trac tickets. I have no idea if this is a new bug or not. Anyway, the version Sage has of DS's scripts are not the most recent one so there is no point in reporting this back to him unless you check that it happens with the scripts on his web page and you email him without mentioning Sage, otherwise he will just assume that the problem lies with Sage (which is possible).
I searched for other tickets related to Simon's code and did not find reports of this particular problem. We should first upgrade to the latest version of his script and then see which of the various problems listed at #11041 persist. I added this ticket to that list for easier reference.
comment:3 Changed 6 years ago by
 Report Upstream changed from N/A to Not yet reported upstream; Will do shortly.
Lorenzini reported me a similar bug and I had a look. I now believe that this is a bug in Simon's script as I can reproduce it outside sage. I used the version of his file dated 06/04/2011, which is the newest I found on the web, while sage has one dated the 25/03/2009. So besides that bug, the file needs to be updated anyway. I tried both in sage's gp (2.5.5) and the current alpha (2.6.1).
I send an email to Denis and see if and how fast we can hope for a solution of this. Otherwise, I may at least update the script, later.
Lorenzini's example is a bit simpler:
R.<t> = QQ[] L.<g> = NumberField(t^3  9*t^2 + 13*t  4) E1 = EllipticCurve(L,[1g*(g1),g^2*(g1),g^2*(g1),0,0]) E1.rank()
comment:4 Changed 6 years ago by
 Description modified (diff)
comment:5 Changed 6 years ago by
 Milestone changed from sage6.1 to sage6.2
comment:6 Changed 6 years ago by
 Dependencies set to #11005
The problem still exists after #11005, which upgrades Simon's script to the latest version.
comment:7 Changed 6 years ago by
 Description modified (diff)
It seems that the bug is caused by a precision problem in subst
, called by nfrealsign
in ell.gp
. The following GP input demonstrates this (using the latest version of ell.gp
, see #11005):
\r src/ext/pari/simon/ell.gp b = 1554544300737274875964190134520312870631312460283689944298138572669148295776039072867720281361776956435252620954745928376624817557704277432961924925312*y + 23524523971732905757341977352314040726186200302188191824300117738073539522011689544444863977622786771332621915440577829842674416407299864303146477224320 a = Mod(b, y^2  229); \p 38 \\ default precision K = bnfinit(y^2  229); nfrealsign(K, a, 1) \\ result: 1 nfrealsign(K, a, 2) \\ result: 1 (wrong) subst(b, y, K.roots[2]) \\ result: 7.695704335233296721 E112, incorrect to this precision \p 500 \\ increase precision K = bnfinit(y^2  229); nfrealsign(K, a, 1) \\ result: 1 nfrealsign(K, a, 2) \\ result: 1 subst(b, y, K.roots[2]) \\ result: 55550556624985845118007242443189926820306719407.796355955086894963616120700422308457186069825115395206111444495243228840597412527828358943327267422898011087182852030228685210492253493260963313348652457717717703991690402543456279919008587884729307897541432624377818611055431792706380900043638904108307170236335925883131494247446074384269328376967381902948861789570537169455258630011202296209679102466188417
The value of K.roots
does not seem to be the problem; the roots are calculated correctly to the requested precision.
comment:8 Changed 6 years ago by
The basic problem seems to be that PARI is too optimistic about the resulting precision when adding real numbers:
? b = 1554544300737274875964190134520312870631312460283689944298138572669148295776039072867720281361776956435252620954745928376624817557704277432961924925312*y + 23524523971732905757341977352314040726186200302188191824300117738073539522011689544444863977622786771332621915440577829842674416407299864303146477224320 ? K = bnfinit(y^2  229); ? u=polcoeff(b,1)*K.roots[2] %167 = 2.3524523971732905757341977352314040726 E151 ? v=polcoeff(b,0)*1. %168 = 2.3524523971732905757341977352314040726 E151 ? u+v %169 = 7.695704335233296721 E112 ? precision(u) %170 = 38 ? precision(v) %171 = 38 ? precision(u+v) %172 = 19
comment:9 Changed 6 years ago by
 Report Upstream changed from Not yet reported upstream; Will do shortly. to Reported upstream. No feedback yet.
Emailed paridev
and Denis Simon about the precision issue, also suggesting a quick fix for ell.gp
.
comment:10 Changed 6 years ago by
The quick fix is to increase the bound in the precision check in nfrealsign
from 10 digits to 38 digits (128 bits). This at least makes the error in my example and that reported by Lorenzini disappear.
comment:11 Changed 6 years ago by
 Branch set to u/pbruin/15483two_descent_precision
 Commit set to d34a77ab7abe4072349e10e0154a598da94b7c31
 Status changed from new to needs_review
comment:12 Changed 6 years ago by
 Status changed from needs_review to needs_work
This needs a doctest...
comment:13 Changed 6 years ago by
 Commit changed from d34a77ab7abe4072349e10e0154a598da94b7c31 to b662aba98ad4122a82fc135cf9b5dbf58de2ba0e
Branch pushed to git repo; I updated commit sha1. New commits:
b662aba  add/fix doctests for simon_two_descent

comment:14 Changed 6 years ago by
 Status changed from needs_work to needs_review
comment:15 Changed 6 years ago by
 Report Upstream changed from Reported upstream. No feedback yet. to Workaround found; Bug reported upstream.
comment:16 followup: ↓ 17 Changed 6 years ago by
I quite like the suggestion by Karim Belabas to use nfsign()
together with nfnewprec
. I think that's a lot better than the current patch, which is really just a hack.
Maybe you should ask Denis Simon what he thinks (perhaps with a proofofconcept patch)? Ideally, he fixes his script "upstream".
comment:17 in reply to: ↑ 16 Changed 6 years ago by
Replying to jdemeyer:
I quite like the suggestion by Karim Belabas to use
nfsign()
together withnfnewprec
. I think that's a lot better than the current patch, which is really just a hack.Maybe you should ask Denis Simon what he thinks (perhaps with a proofofconcept patch)? Ideally, he fixes his script "upstream".
I completely agree with all of this and hope that we will soon hear Denis's opinion (I emailed him and paridev simultaneously). The current patch is indeed fundamentally a hack and will hopefully be superseded by a better solution.
comment:18 followup: ↓ 19 Changed 6 years ago by
Sorry that I have not had time to join this discussion  I am very interested in the outcome.
It would be good to get Denis to adapt his script. He has done this in the past several times. I have found that it does not work for someone else to revise his script and send him the result; he want to do it himself and remain sole author. But a proofofconcept patch might help him.
comment:19 in reply to: ↑ 18 Changed 6 years ago by
Replying to cremona:
It would be good to get Denis to adapt his script. He has done this in the past several times. I have found that it does not work for someone else to revise his script and send him the result; he want to do it himself and remain sole author. But a proofofconcept patch might help him.
Good to know. I also needed to make some changes to his scripts for the PARI2.7 port (#15767).
comment:20 followup: ↓ 21 Changed 6 years ago by
Did you get any answer from Denis Simon? If not, we should try to implement the nfsign()
solution ourselves.
comment:21 in reply to: ↑ 20 Changed 6 years ago by
Replying to jdemeyer:
Did you get any answer from Denis Simon?
Unfortunately not (yet).
If not, we should try to implement the
nfsign()
solution ourselves.
If you want to do it, go ahead, but to be honest, it seems a bit pointless to me at the moment.
Although I agreed above that the increase from 10 to 38 in the precision check was just a hack, there is actually a reason for the 38. (The following assumes a 64bit word length.) For any bound n < 38 in this test, a 1bit error in the result already has catastrophic consequences: a number that should be negative is taken to be a positive number with 64 bits of precision, due to the precision being only 1 bit, and that bit being affected by a roundoff error). For n >= 38, the roundoff errors have to accumulate to at least 65 bits of precision for the problem to occur (wrongly trusting the positivity of a 128bit number). Of course that can happen, but it is extremely unlikely and I would be surprised if you could find any example.
comment:22 Changed 6 years ago by
 Reviewers set to Jeroen Demeyer
 Status changed from needs_review to positive_review
comment:23 Changed 6 years ago by
 Branch changed from u/pbruin/15483two_descent_precision to b662aba98ad4122a82fc135cf9b5dbf58de2ba0e
 Resolution set to fixed
 Status changed from positive_review to closed
There are many similar bug reports with separate trac tickets. I have no idea if this is a new bug or not. Anyway, the version Sage has of DS's scripts are not the most recent one so there is no point in reporting this back to him unless you check that it happens with the scripts on his web page and you email him without mentioning Sage, otherwise he will just assume that the problem lies with Sage (which is possible).
Several of us spend hours or even days in 2011 sage Days (more than one) working on this sort of thing and gave up. We were trying to do more, namely rewrite the interface so that it can use a gp2ccompiled version of the scripts. One day I'll have another go at that.
Meanwhile, this entry is to make sure that noone wastes any time at all looking into the versions of the scripts which Sage currently includes; you should test against the current versions.