upgrade eclib to v20190226

[cremona]

Release v20170815 introduced a modular symbol scaling bug (while sorting out modular symbol scaling essentially correctly) in functions not currently used by Sage. Other releases since the current one in Sage (v20180815) make some interface changes including precision handling (now using bits not digits), so some Sage code changes are likely to be necessary.

Updated source: ​http://homepages.warwick.ac.uk/staff/J.E.Cremona/ftp/eclib-20190226.tar.bz2

[cremona]

Thanks for adding the new keyword. I have built 8.7.beta5 and am starting to test.

[cremona]

The current tarball's configure script sets compiler flags which result in code requiring an optional FLINT module which Sage's FLINT does not have, which I have on my own machine and always forget to remove for Sage. I will replace it. (The module is hmod_mat which allows for linear algebra mod p using 32-bit integers.)

[cremona]

New commits:

​25d6db5

#27360 eclib upgrade to v20190226

[cremona]

I did the testing on a python3 build of Sage.

[jdemeyer]

Typo: procision

[cremona]

Replying to cremona:

I did the testing on a python3 build of Sage.

That is what I thought but it seems that the steps make distclean; make build took me back to python2 (where all is well with the patch uploaded).

I will build the branch from scratch again, this time making sure it is python3.

[jdemeyer]

Replying to cremona:

make distclean .. took me back to python2

Indeed. make distclean is meant to remove everything, including output generated by ./configure indicating that you wanted a Python 3 build.

[cremona]

Thanks -- I misunderstood Vincent's recent instructions on sage-devel. So now the situation is that this branch works fine under python2 but not under python3 where the output buffering issue arises. To be clear: after switching to this commit and putting the new source tarball into upstream/ I do

make distclean
./configure --with-python=3
MAKE='make -j20' make
./sage -t src/sage/libs/eclib

and see doctest failures in 2 files all of which (by eye) are caused by output buffering.

The same happens on a different machine and a python3 build of vanilla 8.7.beta5 so this issue is *not* related to the eclib upgrade at all. For that reason it might be better to open a new ticket for it, and not hold up this one.

[cremona]

As a *proof* that this is the issue, after

 export PYTHONUNBUFFERED=1

the test

./sage -t src/sage/libs/eclib

runs fine. I suppose that we do not want to turn off all python output buffering always?

[jdemeyer]

The last comment is surprising to me... how would disabling Python buffering fix buffering issues in eclib?

[jdemeyer]

Minor nitpick, but why this change?

-    r"""
-    Set the global NTL real number precision.  This has a massive
+    r""" Set the global NTL real number bit precision.  This has a massive

I don't know if it's right or wrong to format docstrings like that, but in any case this seems to go against the typical Sage docstring style. And if you really insist on starting the text on the line with the """ I would at least remove the leading space before Set the.

[jdemeyer]

Is -DNTL_ALL now superfluous? If so, you could remove the line

# distutils: extra_compile_args = -DNTL_ALL

from src/sage/libs/eclib/__init__.pxd

[cremona]

Replying to jdemeyer:

The last comment is surprising to me... how would disabling Python buffering fix buffering issues in eclib?

Well if you ignore my "proof" remark it seems to be a fact that setting that environment variable fixes the issue. Nothing has changed in eclib regarding buffering.

One fix would be to make changes in the Sage interface so that whenever an eclib function is called which outputs to stdout, the python stdout is flushed first. It would be possible to add flushing functions in eclib itself but that does not seem the right approach to me.

[cremona]

Replying to jdemeyer:

Is -DNTL_ALL now superfluous? If so, you could remove the line

# distutils: extra_compile_args = -DNTL_ALL

from src/sage/libs/eclib/__init__.pxd

You are right, and I will remove that line.

[cremona]

Replying to jdemeyer:

Minor nitpick, but why this change?

-    r"""
-    Set the global NTL real number precision.  This has a massive
+    r""" Set the global NTL real number bit precision.  This has a massive

I don't know if it's right or wrong to format docstrings like that, but in any case this seems to go against the typical Sage docstring style. And if you really insist on starting the text on the line with the """ I would at least remove the leading space before Set the.

Sorry that is an artefact introduced by emacs python mode which I usually change back manually. I'll do that.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​8c0f976

#27360 changes addressing reviewer suggestions

[jdemeyer]

I created #27383 for the buffering issues.

[jdemeyer]

There is still procision

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​36bf090

fix typo

[cremona]

Thanks -- I fixed the typo and assume that the buffering questions will be fixed in #27383.

[jdemeyer]

I have no more comments on the branch, now testing...

[jdemeyer]

On a 32-bit machine I get

sage -t --warn-long 43.4 src/sage/libs/eclib/mwrank.pyx
**********************************************************************
File "src/sage/libs/eclib/mwrank.pyx", line 365, in sage.libs.eclib.mwrank._Curvedata.cps_bound
Failed example:
    E.cps_bound()
Expected:
    0.11912451909250964
Got:
    0.11912451909250982
**********************************************************************
File "src/sage/libs/eclib/mwrank.pyx", line 404, in sage.libs.eclib.mwrank._Curvedata.height_constant
Failed example:
    E.height_constant()
Expected:
    0.11912451909250964
Got:
    0.11912451909250982
**********************************************************************
File "src/sage/libs/eclib/mwrank.pyx", line 1344, in sage.libs.eclib.mwrank._two_descent.regulator
Failed example:
    D2.regulator()
Expected:
    0.41714355875838
Got:
    0.417143558758384
**********************************************************************

and

sage -t --warn-long 43.4 src/sage/libs/eclib/interface.py
**********************************************************************
File "src/sage/libs/eclib/interface.py", line 601, in sage.libs.eclib.interface.mwrank_EllipticCurve.regulator
Failed example:
    E.regulator()
Expected:
    0.051111408239969
Got:
    0.05111140823996884
**********************************************************************
File "src/sage/libs/eclib/interface.py", line 1083, in sage.libs.eclib.interface.mwrank_MordellWeil.regulator
Failed example:
    E.regulator()
Expected:
    0.41714355875838
Got:
    0.417143558758384
**********************************************************************

[cremona]

Oh misery. I don't have a 32-bit machine to test on. But I also have no idea why there should be any difference. I don't think I have ever before had to deal with 32/64-bit issues with eclib in Sage. Sage uses NTL floating point (which used to require siwtching on with the NTL_ALL flag but now is on by default) and that is 32/64-independent.

[jdemeyer]

Could there be a problem with the precision handling? I made this change:

src/sage/libs/eclib/wrap.cpp

@@ -250 +250 @@
   bigfloat reg = t->regulator();
   ostringstream instore;
   instore << reg;
-  return stringstream_to_char(instore);
+  char *ptr = stringstream_to_char(instore);
+  printf("two_descent_regulator(): %s\n", ptr);
+  return ptr;
 }

to debug the bigfloat object t->regulator().

On 64-bit, I get

sage: E = mwrank_EllipticCurve([0, 0, 1, -1, 0]); E.regulator()
two_descent_regulator(): 0.051111408239968840235886099756942021609538202280852964249242761521097051665601926126276280166145013158726830156913137573334269246681032033729207105828
0.05111140823996884

On 32-bit, I get

sage: E = mwrank_EllipticCurve([0, 0, 1, -1, 0]); E.regulator()
two_descent_regulator(): 0.0511114082399688402358860997569420216095382
0.05111140823996884

[jdemeyer]

I applied this patch to debug:

src/sage/libs/eclib/wrap.cpp

@@ -250 +250 @@
   bigfloat reg = t->regulator();
   ostringstream instore;
   instore << reg;
-  return stringstream_to_char(instore);
+  char *ptr = stringstream_to_char(instore);
+  printf("two_descent_regulator():\n%s\nprecision = %ld\nas double = %a\n", ptr, reg.precision(), to_double(reg));
+  return ptr;
 }

What I see when I run the doctests on 64-bit:

sage -t src/sage/libs/eclib/interface.py                                                                                                                                
**********************************************************************
File "src/sage/libs/eclib/interface.py", line 601, in sage.libs.eclib.interface.mwrank_EllipticCurve.regulator
Failed example:
    E.regulator()
Expected:
    0.051111408239969
Got:
    two_descent_regulator():
    0.051111408239969
    precision = 53
    as double = 0x1.a2b4645afb41dp-5
    0.051111408239969
**********************************************************************

So the precision is only 53 bits here. So one of the tests preceding it apparently changes the precision from 150 to 53.

[jdemeyer]

Another obvious thing to ask is why does the eclib interface uses strings for conversion? For the regulator, an NTL RR is converted to a string and then to a Python float. It would be more efficient and more reliable to use the NTL conversion function to_double.

(Edit: I don't think that it would fix the doctest failures, it's just something that I saw while reading the code)

[jdemeyer]

The problem is obvious: you have a doctest sage: mwrank_set_precision(50) which means that all following doctests are run with a precision much lower than default.

[cremona]

Replying to jdemeyer:

Another obvious thing to ask is why does the eclib interface uses strings for conversion? For the regulator, an NTL RR is converted to a string and then to a Python float. It would be more efficient and more reliable to use the NTL conversion function to_double.

I agree. I didn't write the wrapping code -- I have only tweaked it so take changes to eclib into account.

Do you know how to convert an NTL RR into an element of RealField?() or some other Sage real type?

(Edit: I don't think that it would fix the doctest failures, it's just something that I saw while reading the code)

[cremona]

Replying to jdemeyer:

The problem is obvious: you have a doctest sage: mwrank_set_precision(50) which means that all following doctests are run with a precision much lower than default.

Of course, how stupid. This hardly needs a doctest; it would be silly to use the default here; perhaps the doctest should be expanded to save the current precision, set it to a new value and then rest it. I will do that.

[cremona]

@jdemeyer: after editing the file wrap.cpp what exactly needs to be done to rebuild sage? I have the impression that neither "sage -b" nor "make" caused that file to be recompiled when it changed.

[jdemeyer]

Yes, I noticed that too and created a new ticket #27389. To solve this problem concretely, remove all copies of wrap.cpp in the src/build directory and try again.

[cremona]

Replying to jdemeyer:

Yes, I noticed that too and created a new ticket #27389. To solve this problem concretely, remove all copies of wrap.cpp in the src/build directory and try again.

Thanks for that tip.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​d63916a

#27360: fix flushing, avoid some stringy interfaces

[cremona]

The latest commit d63916a address some of the issues discovered.

1. By adding a couple of "cou<<flush;" statememts in the c++ wrapper file wrap.cpp I think I have dealt with that problem under python3. At the same time I deleted some flushing which was there already but redundant. The result works fine in python3.
2. Fixed an interface bug where I called eclib's set_precision() instead of set_bit_precision. The former still takes a decimal precision input.
3. Fixed a doctesting bug where one test left the preciaion at a different value than the usual default.
4. Avoid converting an NTL:RR to a float via strings by using iNTL's to_double() as suggested, returning a double.

All that left a lot of computed number ending in some different digits. I fixed those doctests in the stupid way, i.e. changing the expected digits to what I now get. But I only tested on a 64-bit machine. If there are problems on 32-bit I will probably get around that by replacing some trailing digits by "...".

Putting this back to "needs review", and the testing and debugging help already received are much appreciated.

[jdemeyer]

There is still a mwrank_set_precision(50) doctest in src/sage/libs/eclib/interface.py, shouldn't that be fixed too?

[jdemeyer]

Replying to cremona:

At the same time I deleted some flushing which was there already but redundant.

As explained in #27383, it's flushing the wrong buffers anyway (the Python buffers, not the C buffers used by eclib).

[jdemeyer]

Speaking of getting/setting precision, I wonder why there are get_precision and set_precision functions in both src/sage/libs/eclib/mwrank.pyx and src/sage/libs/eclib/interface.py. There is probably a historical reason for this, but it's confusing.

I suggest to remove the get_precision and set_precision functions from src/sage/libs/eclib/interface.py and then adjust the imports in src/sage/libs/eclib/all.py.

Note that the docstring for these functions in src/sage/libs/eclib/mwrank.pyx and src/sage/libs/eclib/interface.py are different, so they should be merged.

[jdemeyer]

I'll fix this.

[cremona]

Replying to jdemeyer:

I'll fix this.

Thanks -- but I also could do it.

The eclib wrapping happened in 2007 or 2006, by William --I made many changes to eclib at the time to make an interface even possible, but did not do the wrapping as at that time I had never written a single line of python. After that William often suggested that it needed to be redone from scratch but the intersection of people who had the expertise to do that was always empty.

[git]

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

​72f217b	#27360 eclib upgrade to v20190226
​2defe8a	#27360: fix flushing, avoid some stringy interfaces
​69d2f64	Allow any kind of iterable for ainvs
​ff2e9dd	Remove get/set_precision from eclib/interface.py
​17635d8	Do not import eclib at startup
​29b2288	Remove superfluous "cdef double" declarations

[jdemeyer]

Rebased to 8.7.beta6 and partially squashed. The last 4 commits are by me and fix various minor things. Each of these commits is more or less independent.

[cremona]

Thanks. The diffs look good and I will do my own testing once the rebuild is done. As soon as we are both happy I think we can set this to Positive Review though I do not know what the strange colour of the patchbot icon means.

[jdemeyer]

Replying to cremona:

I do not know what the strange colour of the patchbot icon means.

That it's a package upgrade and the patchbot doesn't (or can't) know how to deal with those.

[jdemeyer]

This passes all tests in src/sage/libs/eclib/ and src/sage/schemes/elliptic_curves/ both on 32-bit and 64-bit Linux.

[cremona]

When I test this on a python3 build, then there are some doctest failures in elliptic_curves most of which are unrelated unresolved python3 things, except perhaps one in height.py:

File "src/sage/schemes/elliptic_curves/height.py", line 1934, in sage.schemes.elliptic_curves.height.EllipticCurveCanonicalHeight.min_gr
Failed example:
    H.min_gr(.1,5,True) # long time (~22s)
Expected:
    B_1(1) = 1540.199246369678
    ...
    halving mu to 0.25 and increasing n_max to 6
    ...
    halving mu to 0.001953125 and increasing n_max to 13
    doubling mu to 0.0078125
    doubling mu to 0.015625
    height bound in [0.0078125, 0.015625] using n_max = 13
    ...
    height bound in [0.0120485220735, 0.0131390064883] using n_max = 13
    0.012048522073499539
Got:
    B_1(1) = 1540.199246369678
    B_2(1) = 621360.723701339
    B_3(1) = 13686380726.84389
    B_4(1) = 1.6459300096812694e+16
    B_5(1) = 4.322868545515137e+22
    halving mu to 0.25 and increasing n_max to 6
    halving mu to 0.125 and increasing n_max to 7
    halving mu to 0.0625 and increasing n_max to 8
    halving mu to 0.03125 and increasing n_max to 9
    halving mu to 0.015625 and increasing n_max to 10
    halving mu to 0.0078125 and increasing n_max to 11
    halving mu to 0.00390625 and increasing n_max to 12
    halving mu to 0.001953125 and increasing n_max to 13
    doubling mu to 0.0078125
    doubling mu to 0.015625
    height bound in [0.0078125, 0.015625] using n_max = 13
    height bound in [0.011048543456039806, 0.015625000000000003] using n_max = 13
    height bound in [0.011048543456039806, 0.01313900648833929] using n_max = 13
    height bound in [0.012048522073499539, 0.01313900648833929] using n_max = 13
    0.012048522073499539
**********************************************************************

where the penultimate line (only) is slightly different. Did that test really pass for you?

[jdemeyer]

I tested Python 2, not Python 3. Did you test both Python 2 and Python 3 or only Python 3?

I would find it quite surprising that you get different numerical noise on Python 3 compared to Python 2. More likely, it depends on the machine and not on the Python version.

[cremona]

Replying to jdemeyer:

I tested Python 2, not Python 3. Did you test both Python 2 and Python 3 or only Python 3?

Only Python 3

I would find it quite surprising that you get different numerical noise on Python 3 compared to Python 2. More likely, it depends on the machine and not on the Python version.

Me too. The only numerical computation with output other than halving or doubling is line 1975 of height.py where math.sqrt() is used on a float. Perhaps mu and eps should be in RDF from the start, as the output is coerced in there at the end. But this does not explain the discrepancy, it must be a machine difference.

[git]

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

​0e3efa2	#27360 eclib upgrade to v20190226
​424678c	#27360: fix flushing, avoid some stringy interfaces
​dc98893	Allow any kind of iterable for ainvs
​00754ba	Remove get/set_precision from eclib/interface.py
​c2fcfc3	Do not import eclib at startup
​3ff19c1	Remove superfluous "cdef double" declarations

[jdemeyer]

I tested this again (on x86_64 laptop with Sage 8.7.beta7 under Python 2) and all long tests in src/sage/schemes/elliptic_curves/ do pass for me.

[git]

Branch pushed to git repo; I updated commit sha1. New commits:

​8936611

Numerical noise

[jdemeyer]

This should work now.

[cremona]

What is the difference between the commits in comment 56 and those in comment 45?

[jdemeyer]

Replying to cremona:

What is the difference between the commits in comment 56 and those in comment 45?

I just rebased to 8.7.beta7 (which includes #27389)

[cremona]

So should I rebuild my python3 using your rebase and test again what I saw at comment 53? Or should we leave this to be dealt with later when other python3 issues are finally fixed?

[jdemeyer]

I built and tested this on Python 3 myself. I can indeed reproduce the failure from 52. I'll try to figure out why it's happening, if only because I'm honestly curious.

[jdemeyer]

OK, the difference between Python 2 and 3 is that the same numbers are printed differently. So this suggests a different fix: use %r which always print floats exactly.

[git]

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

​901eaf8

Print numbers consistently on Python 2 and 3

[cremona]

Well spotted, and thanks for fixing. I'm still building the previous version but will continue.

[jdemeyer]

I made a few more minor fixes in that function, in particular replacing /= 2 by *= 0.5 (dividing an integer by an integer has different semantics in Python 3, so this might fix a not yet discovered Python 3 bug). I could have written /= 2.0 but decided that the multiplication would be faster too (not that it matters here).

[cremona]

When I pull your commit 901eaf8 on top of 8936611 I get a merge conflict in height.py. That seems to be because the parent of your latest commit is not 8936611 but 3ff16c1. Did you mean to get rid of commit 8936611?

[jdemeyer]

Don't pull. Just checkout my branch.

[cremona]

Replying to jdemeyer:

Don't pull. Just checkout my branch.

Of course. For some reason I had created a new branch tracking trac/u/jdemeyer/eclib. Still, I think it may be bad practice to upload a branch with same name as an old one but which does not have the old one as parent...

[cremona]

OK, so testing all the elliptic_curves directory with a python3 build the only failures I get are ones which have nothing to do with this fix.

For the record there is one failure in each of ell_rational_field, ell_number_field, and padic_lseries. The one in ell_rational_field has been fixed in #27432.

[jdemeyer]

Replying to cremona:

Still, I think it may be bad practice to upload a branch with same name as an old one but which does not have the old one as parent...

I think you're right in principle. On the other hand, sometimes it is more practical to just rebase a branch. For example, having one commit and then another one which reverts that is just adding noise to the git history. In that case, I'd rather just remove the commit.

[cremona]

Replying to jdemeyer:

Replying to cremona:

Still, I think it may be bad practice to upload a branch with same name as an old one but which does not have the old one as parent...

I think you're right in principle. On the other hand, sometimes it is more practical to just rebase a branch. For example, having one commit and then another one which reverts that is just adding noise to the git history. In that case, I'd rather just remove the commit.

*and* as a courtesy to others, add a comment on trac to say that that is what you have done!

[jdemeyer]

Good news! This ticket took a bit longer than expected...

[cremona]

I would like to publicly thank Jeroen for his work on this. It's not the first time that eclib's inclusion in Sage has led to improvements in the library itself.

The less good news is that some of my code it not converging well with the new precision handling (decimal to bits) but not in a part which affects Sage, so there will not be another ticket to upgrade again. Well, not this month.