4

I am trying to find an algorithm, similar to the Pell equation, that would solve for matricies of this form, that have a determinate of 1.

x  nz  ny
y   x  nz
z   y   x

Such a matrix would rapidly converge on the cube root of n. In the shown matrix, z is allowed to be fractional, divisible by some m, where m³ divides n². For example, the case of

n=12, x=55, y=24, z=10.5 works.

Some are quite large: n=1260, x=172226251, y=15945660, z=146337.5

wendy.krieger
  • 7,045
  • 1
    You are looking for norm-one units in $\Bbb Z[\sqrt[3]n]$. They are all powers of a fundamental unit. – Angina Seng Dec 18 '19 at 21:23
  • 1
    See the book "Pell's Equation" by Edward Barbeau Chapter 7, Page 106 for "A Procedure That Often, but Not Always, Works". Available online at the Springer website. – O. S. Dawg Dec 19 '19 at 01:08

1 Answers1

3

The standard method for solving this problem involves using built-in PARI/GP functions to calculate powers of fundamental units, as noted in the comments. However, a more experimental approach is to run the script below, which solves the cubic Pell equation $x^3 + Ny^3 + N^2z^3 - 3Nxyz = 1$ for $N=1260$. There are several noteworthy observations regarding this algorithm:

Rational Solutions: It is often straightforward to find rational matrices in the desired form. To do this, simply identify two $(x,y,z)$ triples with equal norms during the search for a determinant-one solution. Construct the special matrix form for each, then multiply the matrix of the larger triple by the inverse of the smaller. The result is a rational matrix with a determinant of one. Example: For $N=25$, iterations 3 and 13 yield the triples $(8,3,1)$ and $(38483,13161,4501)$, both with norm 12. Combining these produces the rational solution $(6079/4, 2079/4, 711/4)$ of norm one.

Efficiency & Precision: Although this is a two-dimensional, non-Markovian continued fraction algorithm, the memory requirements are lower than the code might suggest. Experiments indicate that initializing matrix $B$ to $B=1.0 \times A$ with a real precision of roughly ten digits yields identical results. This significantly reduces both storage and update costs.

Precomputation: Precomputing $N^{1/3}$ is unnecessary. Instead, you can replace vector $D$ in the script with powers of a matrix $D'$ that exhibits the convergence properties mentioned previously. Additional powers can be added via left multiplication as needed.

Bounds: Empirical evidence suggests that the norms of intermediate approximations remain below roughly $N^6$, a bound that can be used to speed up computation.

A_tmp=[[0, 0, 0; 1, 0, 0; 0, 1, 0], [1, 0, 0; 0, 0, 0; 0, 1, 0]];   B_tmp=A_tmp;    


A=matid(3); B_inv=A;    B=A;


alpha=vector(2);


default(realprecision,1000);


N=1260;


alpha[1]=N^(1/3);   alpha[2]=N^(2/3);


D=vector(3);


D[1]=10^990;


for(i=2,3,D[i]=round(D[1]*alpha[i-1]));


done=0;


while(!done,


for(i=1,2,
    for(j=1,3,A_tmp[i][j,3]=D[j]\D[i]);
    B_tmp[i]=A_tmp[i]^(-1);
);

B_szs=vector(2,i,norml2(B_tmp[i]*B));

B_szs_prm=vecsort(B_szs,,1);

idx=B_szs_prm[1];

A=A*(A_tmp[idx]);
    
B=(B_tmp[idx]*B);

D=D*(B_tmp[idx])~;

done=(1==A[3,3]^3+N*A[2,3]^3+N^2*A[1,3]^3-3*N*A[1,3]*A[2,3]*A[3,3]);


);


print(A[,3]~);

Added 1/16:

I did find the following reference which may be helpful to you:

C. L. E. Wolfe, On the indeterminate cubic equation x 3 + Dy 3 + D 2 z 3 − 3Dxyz = 1, Univ. California pub. Math. 1 (1923) no. 16, 359–369. Available at: University of Michigan Historical Math Collection http://quod.lib.umich.edu/u/umhistmath/ACH6090.0001.001/381

As a check of the algorithm given above I computed the first 100 solutions (in 672 ms) and compared the results to the solutions given in the above paper. It seems that the solutions given for 72 and 82 are wrong. For D=85, I find a smaller solution.

72: 1263601, 303738, 73011

82: 19481839755856966061390701, 4484272707726892450472280, 1032176733268676625192495

85: 658895013725266441, 149856842965183254, 34082931143344968

O. S. Dawg
  • 257