Bill Allombert on Mon, 06 Nov 2023 21:11:34 +0100


[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]

Re: Testing numbers to see if they look algebraic


On Mon, Nov 06, 2023 at 02:13:23PM -0500, Charles Greathouse wrote:
> It often happens that I'm looking to see if a number is algebraic, but I'm
> not sure what degree it might be or how complex the polynomial might be. I
> tend to run a loop like
> 
> \p1000
> for(n=2,100, if(#Str(algdep(x, n)) < 750, print("Possibly algebraic of
> degree ", n)))
> 
> adjusting the degrees, precision, and character tolerance
> (generrally 75-80% of the precision) as needed. If it finds something, I
> can eyeball the polynomial and replicate with higher precision (and
> possibly, depending on the definition of x, prove the relation).
> 
> Is there a better way to do this sort of exploratory analysis? (Of course a
> more complicated version with lindep can search for versions with known
> transcendental constants as well.)

To improve it, one first need to understand why your method works.

Obviously if x is the root of a polynomial of degree n with small coefficients,
it is also the root of a polynomial of degree m with small coefficients for 
all m>n. So why trying small n ?
>From a simple information-theoretic point of view, algdep can only guess
as much information as you give it, so if x is given with  N bits, and you
set the degree to n, you can expect to find polynomials with coefficients
with around N/(n+1) bits. So by decreasing n, you increase the size of the potential
coefficients. However N/(n+1) decreases slowly, so you do not need to test each
possible value of n, but only a geometric sequence (with a small ratio like 1.1).
You can do
my(n=2);for(i=1,37,n+=max(1,n\10);algdep(x, n))

Instead of algdep, I often use lindep(concat(powers(x,n),exp(1)))
and check wether the coefficient of exp(1) is indeed 0!

? \p38
   realprecision = 38 significant digits
? a=ellj((1+sqrt(-23))/2);
? lindep(concat(powers(a,3),exp(1)))
%2 = [-8358020684,858961064,-688165217,-197,-67040526]~ //obviously wrong
? \p100
   realprecision = 115 significant digits (100 digits displayed)
? a=ellj((1+sqrt(-23))/2);
? lindep(concat(powers(a,3),exp(1)))
%4 = [12771880859375,-5151296875,3491750,1,0]~ // very likely to be correct
?  polclass(-23)
%1 = x^3+3491750*x^2-5151296875*x+12771880859375

Cheers,
Bill.