Allow finding minimal polynomial of RIF/CIF element?

user202729 commented 9 months ago

Problem Description

As in the title. I'm pretty sure this is not possible yet.

The motivation for this is that QQbar elements are printed as RIF/CIF printed form (1.2345678?), and it may be desirable to allow reconstructing the QQbar element from the printed form.

That having said, currently QQbar printed form does not necessarily print enough precision to reconstruct the number:

sage: (QQbar(sqrt(2))/10^19+1)
1.000000000000001?

Proposed Solution

Maybe a proof-of-concept implementation (I'm not sure if it works but it seem to give reasonable output most of the time):

from typing import Optional
def CIF_to_QQbar(x: CIF, d: Optional[int]=None)->Optional[QQbar]:
    if d is None:
        best_candidate=(Infinity, None)
        for d in (1..):
            y=CIF_to_QQbar(x, d)
            if y is None: continue
            p=y.minpoly()
            complexity=sum(len(str(abs(coefficient))) for coefficient in p)
            best_candidate=min(best_candidate, (complexity, y))
            if d>=best_candidate[0]: break
        return best_candidate[1]
    def compute_epsilon(l: list[RIF])->RR:
        l=[x for x in l if x!=0]
        if not l: return 1
        return max(x.center()/2^x.precision() if x.is_exact() else x.absolute_diameter() for x in l)
    powers=[x^i for i in range(1, d+1)]
    epsilon_real=compute_epsilon([a.real() for a in powers])
    epsilon_imag=compute_epsilon([a.imag() for a in powers])
    if epsilon_real==0: epsilon_real=min()
    if epsilon_imag==0: epsilon_imag=1  # possible if all numbers are exact
    discretized_powers_real=[round(1/epsilon_real)] + [round(a.real().center()/epsilon_real) for a in powers]
    discretized_powers_imag=[round(1/epsilon_imag)] + [round(a.imag().center()/epsilon_imag) for a in powers]
    global m
    m=matrix.identity(d+1).augment(vector(discretized_powers_real)).augment(vector(discretized_powers_imag)).LLL()
    error=sum(a*b for a, b in zip(m[0][:-2], [1]+powers))
    if error.contains_zero():
        ZZy.<y>=ZZ[]
        return min([a for a, b in ZZy([*m[0][:-2]]).roots(QQbar)],
          key=lambda u: abs(u-x)
          )
    return None

Alternatives Considered

-

Additional Information

No response

Is there an existing issue for this?

[X] I have searched the existing issues for a bug report that matches the one I want to file, without success.

mezzarobba commented 9 months ago

There is a function called algdep in sage.arith.misc that does essentially that. I agree that it would be nice to make it work directly on intervals (not just elements of RIF and CIF, but also of RBF and CBF) and to have a version that returns an element of QQbar or AA instead of the guessed minimal polynomial.

Note that one can print elements of QQbar in a format that allows reconstuction using sage_input.

user202729 commented 9 months ago

Ah, that's nice, evidently I'm reinventing the wheel.

It was too undiscoverable for me though. For discoverability maybe...

minpoly(<RIF element>) and minpoly(<RR element>) raises error message suggesting to use it
<RIF/RR element>.minpoly_approximate() works
QQbar(<RIF/RBF/CIF/CBF/RR/CC element>) works out of the box and/or raise an error message hinting to use the appropriate method.

For the last point: that's true, but the format is quite unwieldy I think (preferably the printed format should be something both readable "by human" and "by computer", to allow copy pasting from terminal output to assign to new variable). (or do you have any better workflow suggestion?)

sagemath / sage