|
Re: cube root of a given number
Posted:
Aug 17, 2007 12:49 AM
|
|
On Aug 12, 2:33?pm, "sttscitr...@tesco.net" <sttscitr...@tesco.net> wrote: > On 12 Aug, 00:56, rich burge <r3...@aol.com> wrote: > > > > > > > On Aug 11, 4:11?pm, "sttscitr...@tesco.net" <sttscitr...@tesco.net> > > > \\ > > \\ > > \\ This script finds a solution to > > P_3(k;x,y,z)=x^3+ky^3+k^2z^3-3kxyz=1 for k=1000700. Here k,x,y,z are > > \\ positive integers. > > \\ > > > nm1(k,x)= > > { > > \\ > > \\ k is k. x is an integer triple, where x[3]/x[1] and x[2]/x[1] are > > the reduced approximations > > \\ to k^(2/3) and k^(1/3) found below. > > \\ > > return((x[3]^3+k*x[2]^3+k^2*x[1]^3-3*k*x[1]*x[2]*x[3])%modk); \\ > > return the norm mod k^6 > > \\ > > > } > > > lessthan(A,B)= > > \\ > > \\ This routine determines which of the two sets of column norms > > corresponding to the reduced approximation > > \\ matricies A and B is 'smallest'. Returns true iff A <* B. > > \\ > > \\ The idea here is this: given two matricies, A and B, the columns of > > which are reduced approximation > > \\ vectors, determine which one has the smallest maximum column norm. > > In the case where the largest column > > \\ norms are equal, compare the second largest, ect. > > \\ > > { > > local(i); > > > \\(re)compute the norms of the columns of A and B > > for(i=1,3,nms_A[i]=nm1(k,A[,i]~);nms_B[i]=nm1(k,B[,i]~)); \\ only 2 > > of these are new!!! > > > \\sortem, > > nms_A_srt=vecsort(nms_A); > > nms_B_srt=vecsort(nms_B); > > > \\ sort through the ties, > > i=3; > > while((i>1)&&nms_A_srt[i]==nms_B_srt[i],i=i-1); > > > \\and return the smallest one. > > return(nms_A_srt[i]<nms_B_srt[i]) > > > } > > > iter()= > > \\ > > \\ This routine acts on D and A. D[1,] is initially an integer > > approximation vector to (1,k^(1/3),k^(2/3)). A is > > \\ the identity matrix. Given D, two reduction matricies B_s and B > > are computed, along with the two possible > > \\ new values of A. lessthan is used to select the correct new A and > > B then D is reduced. > > \\ > > \\ The entries of A, if not controlled, grow rapidly and result in an > > increasing execution time. I.e if a large > > \\ number of iterations are needed to find a solution the process of > > determining what swap to make is painfully slow. > > \\ Note, however, that lessthan(A,B)=lessthan(A mod k^6,B mod k^6), > > given each column norm is bounded by k^6. > > \\ So the computational effort required to make the swap decision is > > bounded. > > \\ > > > { > > B_s=[1,-(D[1,1]\D[1,2]),0;0,-(D[1,3]\D[1,2]),1;0,1,0]; > > A_s=(A*(B_s)^(-1))%modk; > > \\ > > B=[-(D[1,2]\D[1,1]),1,0;-(D[1,3]\D[1,1]),0,1;1,0,0]; > > A=(A*B^(-1))%modk; > > \\ > > \\ this is the key: directed swaps > > \\ > > if(lessthan(A_s,A),A=A_s;B=B_s); > > \\ > > \\ now reduce D accordingly > > \\ > > D=D*B~; > > \\ > > \\ Every fifteenth iteration, add some more precision. 15 is a > > guess!! > > \\ > > if(itercnt%15==0,D=C10*D); > > \\ > > \\ Keep the unreduced values in A_act > > \\ > > A_act=A_act*B^(-1); > > \\ > > \\ (Re)compute the norm and stop if = 1, printing z where > > x=round(z*k^(2/3)) and y=round(z*k^(1/3)). > > \\ > > tmp=nm1(k,A[,3]~); > > if((tmp==1)&&(A_act[1,3]>10),print(vecmin(A_act[, > > 3]));return(0),return(1)); > > > } > > > find_p3()= > > { > > default(timer,1); > > k=1000700; > > modk=k^6; > > A=matid(3); > > A_act=A; > > nms_A=vector(3); > > nms_B=nms_A; > > \\ C is the matrix used to approximate [1,k^(1/3),k^(2/3)] > > \\ The entries are derived from the small rational solution > > [900630049/49,9004200/49,90021/49] to > > \\ P_3(k,x,y,z)=1. These small rational solutions are fairly easy to > > find. > > \\ > > C=[900630049,9004200,90021; > > k*90021,900630049,9004200; > > k*9004200,k*90021,900630049]~; > > C10=C; > > D=C; > > itercnt=0; > > while(iter(),itercnt=itercnt+1); > > print(); > > > } > > > find_p3; > > > \\ end script > > Very impressive. But I can't say I understand > how it all works at the moment. > What's the largest k you have found solutions for ?
The largest solution I have found so far is for k=1000007000 . It takes about 11 million iterations and can be done in about 90 minutes (w/o updating A_act). This requires a real good choice for C. I was hoping to find enough particular solutions to the family 10^(3*e) +7*10^e to find a parametric solution, but this does not seem feasible at this point.
Rich
|
|