Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Allow finding minimal polynomial of RIF/CIF element? #36798

Open
1 task done
user202729 opened this issue Dec 2, 2023 · 3 comments · May be fixed by #39017
Open
1 task done

Allow finding minimal polynomial of RIF/CIF element? #36798

user202729 opened this issue Dec 2, 2023 · 3 comments · May be fixed by #39017

Comments

@user202729
Copy link
Contributor

user202729 commented Dec 2, 2023

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?

  • I have searched the existing issues for a bug report that matches the one I want to file, without success.
@mezzarobba
Copy link
Member

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
Copy link
Contributor Author

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?)

@user202729
Copy link
Contributor Author

The funny part is currently minpoly(SR([element of RR])) actually works:

sage: minpoly(1.2)
[error]
sage: minpoly(1.2*x/x)
x - 6/5
sage: minpoly(SR(1.2))
x - 6/5

but doesn't always work as expected:

sage: minpoly(SR(RR(sqrt(2))))
x - 6369051672525773/4503599627370496

@user202729 user202729 linked a pull request Nov 21, 2024 that will close this issue
5 tasks
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

Successfully merging a pull request may close this issue.

2 participants