Opened 7 years ago

Last modified 6 years ago

#18836 needs_work enhancement

Make refine_embedding into a method of number fields instead of stand-alone

Reported by: cremona Owned by:
Priority: minor Milestone: sage-6.8
Component: number fields Keywords:
Cc: Merged in:
Authors: John Cremona Reviewers: Vincent Delecroix
Report Upstream: N/A Work issues:
Branch: u/cremona/18836-refine_embedding (Commits, GitHub, GitLab) Commit: 94a1d567b17399c9e39517ab6943004210ee23ba
Dependencies: #20064 Stopgaps:

Status badges

Description

The stand-alone function {{{refine_embedding}} has been around for a while, right at the bottom of sage/rings/number_field/number_field.py, not imported into the global namespace and hence having to be imported whenever needed. More seriously, it is very easy for users and developers to miss its existence entirely! (Proof: today I learned that one of my students needed such a function and wrote his own since he did not know that I had written the existing function a few years ago!)

I am moving the function into the class NumberField?, so that one can now say

embx = K.refine_embedding(emb, prec)

instead of

from sage.rings.number_field.number_field import refine_embedding
embx = refine_embedding(emb,prec)

There are only about 3 places in the Sage library where this is used, which need small changes.

Change History (18)

comment:1 Changed 7 years ago by cremona

  • Branch set to u/cremona/18836-refine_embedding
  • Commit set to 87cb9bec69f77d1a3b4f33dacef1b7d6043c5062
  • Status changed from new to needs_review

New commits:

87cb9be#18836: mvoe refine_embedding function into NumberField class

comment:2 Changed 7 years ago by vdelecroix

  • Reviewers set to Vincent Delecroix
  • Status changed from needs_review to needs_work

Hello John,

  • Could you be more uniform in notations: you wrote RR and CC but ``AA`` and ``QQbar``.
  • if not self is e.domain() -> if self is not e.domain()
  • it is crazy that this function needs to compute all embeddings to refine one!! (might be for another ticket)
  • This is very obfuscated
            diffs = [(RC(ee(self.gen()))-old_root).abs() for ee in elist]
            return elist[min(izip(diffs,count()))[1]]
    
    Could you make it more readable like
       diff_min = Infinity
       ans = None
       for ee in elist:
           diff = (RC(ee(self.gen())) - old_root).abs()
           if diff < diff_min:
               ans = ee
               diff_min = diff
       return ans
    

Vincent

comment:3 Changed 7 years ago by cremona

OK, I am working on it. On your second point I am not sure that it is so wasteful to find all roots of a polynomial in a real / complex field than to refine one approximate root to higher precision, but I will experiment to see if that is possible. And of course, this code (mostly by me) has been in Sage for a few years already, all I did here was to move it!

comment:4 Changed 7 years ago by git

  • Commit changed from 87cb9bec69f77d1a3b4f33dacef1b7d6043c5062 to 87e8c14747eb0f5bf9b877b971b2ab2d5edbbffd

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

450209d#18836: mvoe refine_embedding function into NumberField class
accaa9f#18836: simplified code after review
87e8c14Merge branch 'u/cremona/18836-refine_embedding' of trac.sagemath.org:sage into refine_embedding

comment:5 Changed 7 years ago by cremona

I did what was asked for (almost) and rebased onto 6.8.beta7. There does not seem to be a way to ask for one root in a field close to a given approximation, so it still finds all roots in the higher precision field -- but only the roots, it does not contruct all teh associated embeddings, just the closest one. I think the code is less obfuscated now too.

There is still a problem,though: in testing I get two errors (from essentially the same test) in sage/schemes/elliptic_curves/ell_point.py and period_lattice.py. These occur when a complex embedding has been refined to infinite precision (i.e. to an embedding into QQbar), but an error is caused when an element of QQbar has its square root taken: its _value is an element of sage.rings.real_mpfi.RealIntervalFieldElement? which is surely wrong -- it should be the corresponding complex type -- and does not like its argument being asked for.

I have not been able to track this down, so I have not set the ticket to "needs review" as I know that these tests fail, but I am guessing that I have somehow stumbled across a bug which has nothing to do with this ticket. Help appreciated!

comment:6 Changed 6 years ago by git

  • Commit changed from 87e8c14747eb0f5bf9b877b971b2ab2d5edbbffd to ee94ab4b6da7a2bdf20b6fc50eefb2dba9190494

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

ee94ab4Merge branch 'master' (6.10) into refine_embedding

comment:7 Changed 6 years ago by cremona

I have merged this with 6.10 to make it possible to continue with it.

comment:8 Changed 6 years ago by cremona

This code was written while debugging the above. It runs fine unless you comment out the two lines near the end. I would say that this is a bug in QQbar.sqrt!

K.<i> = QuadraticField(-1)

# define a low-precision embedding from K to CC:                                                                

emb = K.embeddings(CC)[1]

# extend this to the closest embedding into QQbar:                                                              

old_gen = emb(K.gen())
rr = K.defining_polynomial().roots(QQbar, multiplicities=False)
diffs = [(CC(r)-old_gen).abs() for r in rr]
new_gen = rr[diffs.index(min(diffs))]
emb0 = K.hom([new_gen], check=False)

# Take a polynomial with 3 roots in K:                                                                          

e1 = -4+i
e2 = 1+i
e3 = 3-2*i
print("Original ei: %s with parent %s" % ([e1,e2,e3],parent(e1)))
x = polygen(K)
pol = (x-e1)*(x-e2)*(x-e3)

# Find the roots again in QQbar:                                                                                

pol0 = PolynomialRing(QQbar,'x')([emb0(c) for c in list(pol)])
e1, e2, e3 = pol0.roots(QQbar,multiplicities=False)
print("Roots ei: %s with parent %s" % ([e1,e2,e3],parent(e1)))

# Attempt to compute sqrt(e1-e2) from these:                                                                    

d = e1-e2
print("d=%s with parent %s" % (d,d.parent()))
# If the next 2 lines are commented out, an error is raised in the sqrt!                                        
s = d.imag().is_zero()
print("d.imag().is_zero()=%s" % s)
print("d=%s with parent %s" % (d,d.parent()))
d = d.sqrt()
print("d=%s" % d)

comment:9 Changed 6 years ago by nbruin

The fact that _value becomes a real interval element is probably due to this implementation (line 7431 of sage/rings/qqbar.py)

    def _interval_fast(self, prec):
        gen_val = self._generator._interval_fast(prec)
        v = self._value.polynomial()(gen_val)
        if self._exactly_real and is_ComplexIntervalFieldElement(v):
            return v.real()
        return v

It could well be that this is really the intention (and there could be other places where the interval is set to be a real thing!), in which case the bug is indeed in sqrt, which should avoid relying on "argument" being available.

Looking at it a little bit more, I think the error is in fact in __pow__ (sage/rings/qqbar.py:4233). I expect it's not invalid for _value to be a RIF element. That seems to happen quite a bit:

sage: type(QQbar(5)._value)
<type 'sage.rings.real_mpfi.RealIntervalFieldElement'>

It's just that elements where that happens are usually filtered out earlier. So in this case, d isn't recognized as a rational number yet when we enter, but then in the for loop on line 4304, we have that on line 4322 we execute:

               val = self._interval_fast(prec)

So actually, on the first pass val is a CIF element. We apparently get an error on the second pass, when val has been forced through _interval_fast.

So I see two solutions: either ensure that val is put back into CIF or ensure that "argument extraction" is done appropriately in cases where "val" is in RIF.

Last edited 6 years ago by nbruin (previous) (diff)

comment:10 Changed 6 years ago by cremona

That is good, I thought that you (Nils) would be able to help with this. I think we should fix this on #20068, which I made a dependency for this ticket, so I will copy your comments there.

comment:11 Changed 6 years ago by vdelecroix

Hello,

I followed the comment from #18333. It would actually make sense for elements in QQbar to always have their intervals being in some complex interval field. Since the zero is exact in interval arithmetic we can check if the number is real with

sage: a = CIF(5)
sage: a.imag().is_zero()
True

Though there is a lot of code shared between AA and QQbar and the change might be less trivial than it seems.

Vincent

comment:12 Changed 6 years ago by cremona

Thanks, I did not know of the meta-ticket #18333. If this issue (now at #20064) can be resolved as part of that, then good, but I hope that a quick resolution will be possible since #20064 is a real bug while it looks as if many of the issues at #18333 are less urgent (though worth doing, certainly).

comment:13 Changed 6 years ago by git

  • Commit changed from ee94ab4b6da7a2bdf20b6fc50eefb2dba9190494 to 94a1d567b17399c9e39517ab6943004210ee23ba

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

94a1d56Merge branch 'develop (7.1.beta3)' into 18836

comment:14 follow-up: Changed 6 years ago by cremona

I am about to see if the fix at #20064 deals with the remaining issues.

comment:15 in reply to: ↑ 14 Changed 6 years ago by cremona

  • Dependencies set to #20064
  • Status changed from needs_work to needs_review

Replying to cremona:

I am about to see if the fix at #20064 deals with the remaining issues.

It does, and now has a positive review, so I will set that as a dependency of this and ask for another review of this one.

comment:16 Changed 6 years ago by lftabera

  • Status changed from needs_review to needs_work

I got doctest failing for this one.

sage -t --long src/sage/schemes/elliptic_curves/ell_point.py # 1 doctest failed sage -t --long src/sage/schemes/elliptic_curves/period_lattice.py # 2 doctests failed

comment:17 Changed 6 years ago by cremona

Are those failures the ones which are fixed by #20064? Did you apply #20064 first -- it is merged but not yet in the latest beta, I think.

comment:18 Changed 6 years ago by lftabera

You are right, I thought #20064 was already merged, while it is not. Everything is fine now.

  • The code is ok, I miss extending embeddings from number to p-adic fields, but this issue is out of scope for this ticket.
  • Please document that precision cannot decrease.
sage: N=QQ[sqrt(2)]
sage: e=N.embeddings(RR)[0]
sage: N.refine_embedding(e,16) is e
True
sage: e
Ring morphism:
  From: Number Field in sqrt2 with defining polynomial x^2 - 2
  To:   Real Field with 53 bits of precision
  Defn: sqrt2 |--> -1.41421356237310

Also, it would be nice to note that, when the refinement is not unique, an "arbitrary" one is returned:

sage: N=NumberField(x^4 - 199999999/50000000*x^2 + 40000000400000001/10000000000000000,'a')
sage: e = N.embeddings(ComplexField(16))[1]
sage: e
Ring morphism:
  From: Number Field in a with defining polynomial x^4 - 199999999/50000000*x^2 + 40000000400000001/10000000000000000
  To:   Complex Field with 16 bits of precision
  Defn: a |--> 1.414
sage: N.refine_embedding(e,32)
Ring morphism:
  From: Number Field in a with defining polynomial x^4 - 199999999/50000000*x^2 + 40000000400000001/10000000000000000
  To:   Complex Field with 32 bits of precision
  Defn: a |--> 1.41421356 - 0.0000988882552*I
sage: len(N.embeddings(ComplexField(16)))
2
sage: len(N.embeddings(ComplexField(32)))
4
Note: See TracTickets for help on using tickets.