Numerical glitch causes error in computation of fundamental group of curves

[mmarco]

In Sage 8.3 we get the following (wrong) answer:

sage: var('t')
sage: wp = QQ[t](t^2 + t + 1).roots(QQbar)[0][0]
sage: Kw.<wp> = NumberField(wp.minpoly(), embedding=wp)
sage: R.<x,y> = Kw[]
sage: z = 1
sage: f = y * (y + (-wp - 1)) * x * (x - 1) * (x - y) * (x + (-wp - 1)*y - 1) * (x + (-wp - 1)*y + (wp))
sage: from sage.schemes.curves import zariski_vankampen as zvk
sage: zvk.fundamental_group(f)
Finitely presented group < x0, x1, x2, x3, x5 | x2*x1*x2^-1*x1^-1, x1*x0*x1^-1*x0^-1, x5*x2*x5^-1*x2^-1, x0^-1*x3*x0*x3^-1, x0^-1*x2*x0*x2^-1, x1*x3*x1^-1*x3^-1, x5^-1*x0*x5*x0^-1, x3^-1*x5*x3*x5^-1, x1*x5*x1^-1*x5^-1 >

It must be wrong because the abelianized must have rank 7, so it must have at least 7 generators.

Note that Sage 8.1 used to give right answers.

Anyways, I think I did find the problem. When computing the braids corresponding to each segment, the endpoints are represented by floating point numbers, which may cause some errors. Using algebraic numbers for this seems to solve the issue.

[slelievre]

Trying to see if Miguel Marco's branch is still there.

---

New commits:

​d0ee394	Use algebraic numbers to pin the beginning and end of braids
​cb4c624	Added docytest
​7b55efd	Merge branch 'fixzvk' into t/26503/numerical_glitch_causes_error_in_computation_of_fundamental_group_of_curves

[tscrim]

So this change will likely slow things down, but it will make things more numerically stable (as evidenced by this example).

One thing that should be changed:

-    TESTS::
+    TESTS:
 
-        # Check that # 26503 is fixed
+    Check that :trac:`26503` is fixed::

[slelievre]

Pushing the suggestion by Travis Scrimshaw a bit further, I would change

    TESTS::

        # Check that # 26503 is fixed
        sage: var('t')
        sage: wp = QQ[t](t^2+t+1).roots(QQbar)[0][0]
        sage: Kw.<wp>=NumberField(wp.minpoly(), embedding=wp)
        sage: R.<x,y> = Kw[]
        sage: z = 1
        sage: f = y * (y + (-wp - 1)) * x * (x - 1) * (x - y) * (x + (-wp - 1)*y - 1) * (x + (-wp - 1)*y + (wp))
        sage: from sage.schemes.curves import zariski_vankampen as zvk
        sage: g = f.subs({x:x+y})
        sage: g = g.subs({x:x+y})
        sage: disc = zvk.discrim(g)
        sage: segs = zvk.segments(disc)
        sage: segs[16]
        (0.577350269189626*I, 0.500000000000000 + 0.288675134594813*I)
        sage: zvk.braid_in_segment(g, *segs[16])  # optional - sirocco
        s0*s1*s3*s5*s1^-1*s0^-1*s3^-1

to the following, if you agree it means the same (please double check):

    TESTS:

    Check that :trac:`26503` is fixed::
        
        sage: wp = QQ['t']([1, 1, 1]).roots(QQbar)[0][0]
        sage: Kw.<wp> = NumberField(wp.minpoly(), embedding=wp)
        sage: R.<x, y> = Kw[]
        sage: z = -wp - 1
        sage: f = y*(y + z)*x*(x - 1)*(x - y)*(x + z*y - 1)*(x + z*y + wp)
        sage: from sage.schemes.curves import zariski_vankampen as zvk
        sage: g = f.subs({x: x + 2*y})
        sage: disc = zvk.discrim(g)
        sage: segs = zvk.segments(disc)
        sage: segs[16]
        (0.577350269189626*I, 0.500000000000000 + 0.288675134594813*I)
        sage: zvk.braid_in_segment(g, *segs[16])  # optional - sirocco
        s0*s1*s3*s5*s1^-1*s0^-1*s3^-1

This combines some pep8 fixes, getting rid of a symbolic variable, using z to clarify f, and combining two .subs into a single one.

[mmarco]

Replying to tscrim:

So this change will likely slow things down, but it will make things more numerically stable (as evidenced by this example).

Yes, but on the other hand, the other minor change I made that computes radical of the polynomial before computing the discriminant would make things faster, so overall the change could actually be faster.

Samuel: I think your change is better. I will push a commit with it tomorrow.

[git]

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

​5ba2b90

Improve doctests with reviewers suggestions

[slelievre]

A blank line is needed between :: and the code block.

Also, please consider simplifying

        sage: g = f.subs({x:x+y})
        sage: g = g.subs({x:x+y})

into

        sage: g = f.subs({x: x + 2*y})

[git]

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

​6e96677

Minor fix in the doctest

[chapoton]

I have made some trivial changes, and I allow myself to turn this to positive.

---

New commits:

​fd8864e	Merge branch 'u/mmarco/numerical_glitch_causes_error_in_computation_of_fundamental_group_of_curves' in 8.5.b5
​67f1188	trac 26503 some minor details

[vbraun]

On 32-bit:

Doctesting 3781 files using 2 threads.
**********************************************************************
File "src/sage/schemes/curves/zariski_vankampen.py", line 333, in sage.schemes.curves.zariski_vankampen.braid_in_segment
Failed example:
    segs[16]
Expected:
    (0.577350269189626*I, 0.500000000000000 + 0.288675134594813*I)
Got:
    (3.82181265801140e-17 + 0.577350269189626*I,
     -0.500000000000000 + 0.288675134594813*I)
**********************************************************************
1 item had failures:
   1 of  15 in sage.schemes.curves.zariski_vankampen.braid_in_segment
    [35 tests, 1 failure, 0.91 s]
----------------------------------------------------------------------
sage -t --long src/sage/schemes/curves/zariski_vankampen.py  # 1 doctest failed
----------------------------------------------------------------------

[git]

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

​1269f8c

#26503 Add abs tol to one doctest

[slelievre]

Added # abs tol 1e-16.

[chapoton]

ok

[vbraun]

On 32-bit:

File "src/sage/schemes/curves/zariski_vankampen.py", line 333, in sage.schemes.curves.zariski_vankampen.braid_in_segment
Failed example:
    segs[16]  # abs tol 1e-16
Expected:
    (0.577350269189626*I, 0.500000000000000 + 0.288675134594813*I)
Got:
    (3.82181265801140e-17 + 0.577350269189626*I,
     -0.500000000000000 + 0.288675134594813*I)
**********************************************************************
1 item had failures:
   1 of  15 in sage.schemes.curves.zariski_vankampen.braid_in_segment
    [35 tests, 1 failure, 0.79 s]

[chapoton]

There is a sign problem, not only a matter of tolerance..

[slelievre]

Seems the order of the elements in disc can vary.

Maybe exactifying the elements of disc and sorting them would fix that:

        sage: wp = QQ['t']([1, 1, 1]).roots(QQbar)[0][0]
        sage: Kw.<wp> = NumberField(wp.minpoly(), embedding=wp)
        sage: R.<x, y> = Kw[]
        sage: z = -wp - 1
        sage: f = y*(y + z)*x*(x - 1)*(x - y)*(x + z*y - 1)*(x + z*y + wp)
        sage: from sage.schemes.curves import zariski_vankampen as zvk
        sage: g = f.subs({x: x + 2*y})
        sage: disc = zvk.discrim(g)
        sage: for a in disc: a.exactify()
        sage: disc = sorted(disc)
        sage: segs = zvk.segments(disc)
        sage: segs[16]  # abs tol 1e-16
        (0.577350269189626*I, 0.500000000000000 + 0.288675134594813*I)

[mmarco]

Replying to slelievre:

Seems the order of the elements in disc can vary.

Maybe exactifying the elements of disc and sorting them would fix that:

        sage: wp = QQ['t']([1, 1, 1]).roots(QQbar)[0][0]
        sage: Kw.<wp> = NumberField(wp.minpoly(), embedding=wp)
        sage: R.<x, y> = Kw[]
        sage: z = -wp - 1
        sage: f = y*(y + z)*x*(x - 1)*(x - y)*(x + z*y - 1)*(x + z*y + wp)
        sage: from sage.schemes.curves import zariski_vankampen as zvk
        sage: g = f.subs({x: x + 2*y})
        sage: disc = zvk.discrim(g)
        sage: for a in disc: a.exactify()
        sage: disc = sorted(disc)
        sage: segs = zvk.segments(disc)
        sage: segs[16]  # abs tol 1e-16
        (0.577350269189626*I, 0.500000000000000 + 0.288675134594813*I)

Yes, that is what I think too. The failing doctest relies on the order of the roots of the discriminant being fixed, which apparently might not be depending on the architecture (as the failing doctest shows). Samuel's proposal makes more sense.

[git]

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

​98c921d

Sort to make braid_in_segment doctest pass on 32-bit

[slelievre]

Pushed a commit to add this exactifying and sorting. Can someone test on 32-bit?

[tscrim]

Patchbot seems to have trouble with it, plus one of them is getting an error with sirocco:

File "src/sage/schemes/curves/zariski_vankampen.py", line 337, in sage.schemes.curves.zariski_vankampen.braid_in_segment
Failed example:
    zvk.braid_in_segment(g, *segs[16])  # optional - sirocco
Expected:
    s0*s1*s3*s5*s1^-1*s0^-1*s3^-1
Got:
    s5^-1*s3*s0*s1*s0*s3^-1*s4^-1*s3^-1*s2^-1

[tmonteil]

Let me run a 32-bit patchbot on it.

[tmonteil]

My 32-bit patchbot is unhappy :(

[slelievre]

Order seems to depend on the platform, so we can't rely on segs[16] and we need a better doctest.

[mmarco]

Can you check if the discrepancy on the order appears on disc or segs ?

Anyways, maybe the right thing to do about this doctest is to be explicit and give the specific segment to test?

[tscrim]

What precisely is that part of the doctest suppose to be testing?

[mmarco]

It is testing the computation of the braid along the segment. The segments are created from the Voronoi diagram of the points in disc. We compute them all, and pick a specific one (where the bug shows up) and test it. The problem seems to be that the orderof the computed segments seems to deppend on the platform, so we are picking a different one on 32 bits, and get different results.

[tscrim]

So if the bug is showing up in that specific segment, perhaps the thing to do is isolate that segment, maybe by checking if that segment is in segs and then passing it to braid_in_segments?

[tscrim]

Ah, wait, that won't work outright because of the numerical noise.

[tscrim]

Perhaps instead of that specific segment, we check

s0*s1*s3*s5*s1^-1*s0^-1*s3^-1 in [zvk.braid_in_segment(g,s) for s in segs]

I guess we could then get the index of that braid in the corresponding list and then check the segment itself (if you also want to show that output).

[git]

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

​a348797	Use algebraic numbers to pin the beginning and end of braids
​b3cd564	Added docytest
​446a5fb	Improve doctests with reviewers suggestions
​666941e	Minor fix in the doctest
​581cf35	trac 26503 some minor details
​2f081a6	#26503 Add abs tol to one doctest
​5c74d7a	Sort to make braid_in_segment doctest pass on 32-bit
​1a1af79	Explicit path in doctest for braid_in_segment
​5012a7a	Merge branch 'public/ticket/26503' of git://trac.sagemath.org/sage into t/26503/public/ticket/26503

[mmarco]

I changed the doctest by giving explicitely the endpoints of the path. It still detects the bug.

I didn't test it in a 32 bits system though.

[git]

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

​bfa164b

Use left normal form of braid in doctest

[tscrim]

The changes are good with me, but you need to mark 2 doctests with # optional - sirocco (see pathcbot). Thierry, can you run your 32-bit patchbot on this?

[git]

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

​0a95947	Merge branch 'public/ticket/26503' in 8.7.b7
​b91818d	adding the optional - sirocco tags

[chapoton]

I have added the 2 missing sirocco tags

[tscrim]

Sorry, this fell off my radar by the time I got to a place where I could have made those modifications. Thank you Frédéric for the fix and the reminder.

Can someone with a 32-bit machine retest this for completeness?

[chapoton]

@tmonteil: Thierry, please launch your 32 bit patchbot on this ticket

[tmonteil]

OK. I have to reboostrap one patchbot, so do not expect some result before this evening.

[tscrim]

Replying to tmonteil:

OK. I have to reboostrap one patchbot, so do not expect some result before this evening.

No problem. Thank you for doing that.

[tmonteil]

This looks good now on 32bit, see ​https://patchbot.sagemath.org/log/26503/debian/9.7/i686/4.9.0-7-686/tmonteil-debian-stretch-32/2019-03-21%2013:03:07.

There is an error for the file src/sage/tests/books/computational-mathematics-with-sagemath/graphique_doctest.py but it also exists for ticket 0, see ​https://patchbot.sagemath.org/ticket/0/ so it is unrelated to this ticket.

[tmonteil]

The error is actually tracked at #27250.

[tscrim]

Thank you. I am now setting this to positive review.

[slelievre]

There's been a lot of trial and error on this, do we want to keep all these commits? We could squash to a single commit instead -- I'm fine either way.

[embray]

Ticket retargeted after milestone closed (if you don't believe this ticket is appropriate for the Sage 8.8 release please retarget manually)