Opened 6 years ago

Closed 6 years ago

# Bug in sqrt in QQbar

Reported by: Owned by: cremona major sage-7.1 numerical Nils Bruin John Cremona N/A 48c12ef 48c12efbf46cbd4143b1c605004bd6224461b3d8

### Description

See #18836. This bug is holding up that (and also #20028). The following code creates an elemnt d of QQbar and tries to do d.sqrt(). It fails unless you call d.imag().is_zero() first.

# 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:1 Changed 6 years ago by cremona

(Extracted from #18836 comment 9)

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.

### comment:2 follow-up: ↓ 3 Changed 6 years ago by cremona

Surely it is simplest to implement arg for RIF elements, returning 0 or Pi depending on sign provided that the interval has constant sign, and raising an error otherwise?

### comment:3 in reply to: ↑ 2 Changed 6 years ago by nbruin

Surely it is simplest to implement arg for RIF elements, returning 0 or Pi depending on sign provided that the interval has constant sign, and raising an error otherwise?

I would regard that as a "change in design/API", though. I think our real fields currently quite consistently do not have an "argument" method:

sage: CC(-4).argument()
3.14159265358979
sage: RR(-4).argument()
AttributeError: ...
sage: CDF(-4).argument()
3.141592653589793
sage: RDF(-4).argument()
AttributeError: ...

Your proposal is in line with the solution taken in #18337 (to put "real" and "imag" on RIF) so perhaps putting "argument" there is the simplest solution.

From an efficiency point of view: I'm not sure how hard we'd be hitting the coercion system by mixing RIF and CIF elements.

### comment:4 Changed 6 years ago by nbruin

• Branch set to u/nbruin/bug_in_sqrt_in_qqbar

### comment:5 Changed 6 years ago by nbruin

• Authors set to Nils Bruin
• Commit set to 48c12efbf46cbd4143b1c605004bd6224461b3d8
• Component changed from number fields to numerical
• Status changed from new to needs_review

New commits:

 ​48c12ef trac #20064: introduce "argument" on RealIntervalFieldElement

### comment:6 Changed 6 years ago by cremona

• Reviewers set to John Cremona
• Status changed from needs_review to positive_review

This looks good to me. I have checked that it deals with the original problems I had at #18836, so I am going to set this to positive_review, make that ticket depend on this and set that one to needs_review.

### comment:7 Changed 6 years ago by vbraun

• Branch changed from u/nbruin/bug_in_sqrt_in_qqbar to 48c12efbf46cbd4143b1c605004bd6224461b3d8
• Resolution set to fixed
• Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.