version 1.1, 2001/10/02 11:17:12 |
version 1.2, 2002/09/11 07:27:06 |
Line 13 Check the License for details. You should have receive |
|
Line 13 Check the License for details. You should have receive |
|
with the package; see the file 'COPYING'. If not, write to the Free Software |
with the package; see the file 'COPYING'. If not, write to the Free Software |
Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ |
Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. */ |
|
|
|
#include "pari.h" |
|
|
/* Thue equation solver. In all the forthcoming remarks, "paper" |
/* Thue equation solver. In all the forthcoming remarks, "paper" |
* designs the paper "Thue Equations of High Degree", by Yu. Bilu and |
* designs the paper "Thue Equations of High Degree", by Yu. Bilu and |
* G. Hanrot, J. Number Theory (1996). Note that numbering of the constants |
* G. Hanrot, J. Number Theory (1996). Note that numbering of the constants |
* is that of Hanrot's thesis rather than that of the paper. |
* is that of Hanrot's thesis rather than that of the paper. |
* The last part of the program (bnfisintnorm) was written by K. Belabas. |
* The last part of the program (bnfisintnorm) was written by K. Belabas. */ |
*/ |
|
#include "pari.h" |
|
|
|
static int curne,r,s,t,deg,Prec,ConstPrec,numroot; |
|
static GEN c1,c2,c3,c4,c5,c6,c7,c8,c10,c11,c12,c13,c14,c15,B0,x1,x2,x3,x0; |
|
static GEN halpha,eps3,gdeg,A,MatNE,roo,MatFU,Delta,Lambda,delta,lambda; |
|
static GEN Vect2,SOL,uftot; |
|
|
|
GEN bnfisintnorm(GEN, GEN); |
|
|
|
/* Check whether tnf is a valid structure */ |
/* Check whether tnf is a valid structure */ |
static int |
static int |
checktnf(GEN tnf) |
checktnf(GEN tnf) |
{ |
{ |
|
int n, R, S, T; |
if (typ(tnf)!=t_VEC || (lg(tnf)!=8 && lg(tnf)!=3)) return 0; |
if (typ(tnf)!=t_VEC || (lg(tnf)!=8 && lg(tnf)!=3)) return 0; |
if (typ(tnf[1])!=t_POL) return 0; |
if (typ(tnf[1]) != t_POL) return 0; |
if (lg(tnf) != 8) return 1; /* s=0 */ |
if (lg(tnf) != 8) return 1; /* S=0 */ |
|
|
deg=degpol(tnf[1]); |
n = degpol(tnf[1]); |
if (deg<=2) err(talker,"invalid polynomial in thue (need deg>2)"); |
if (n <= 2) err(talker,"invalid polynomial in thue (need n>2)"); |
s=sturm((GEN)tnf[1]); t=(deg-s)>>1; r=s+t-1; |
S = sturm((GEN)tnf[1]); T = (n-S)>>1; R = S+T-1; |
(void)checkbnf((GEN)tnf[2]); |
(void)checkbnf((GEN)tnf[2]); |
if (typ(tnf[3]) != t_COL || lg(tnf[3]) != deg+1) return 0; |
if (typ(tnf[3]) != t_COL || lg(tnf[3]) != n+1) return 0; |
if (typ(tnf[4]) != t_COL || lg(tnf[4]) != r+1) return 0; |
if (typ(tnf[4]) != t_COL || lg(tnf[4]) != R+1) return 0; |
if (typ(tnf[5]) != t_MAT || lg(tnf[5]) != r+1 |
if (typ(tnf[5]) != t_MAT || lg(tnf[5])!=R+1 || lg(gmael(tnf,5,1)) != n+1) return 0; |
|| lg(gmael(tnf,5,1)) != deg+1) return 0; |
if (typ(tnf[6]) != t_MAT || lg(tnf[6])!=R+1 || lg(gmael(tnf,6,1)) != R+1) return 0; |
if (typ(tnf[6]) != t_MAT || lg(tnf[6])!=r+1 |
|
|| lg(gmael(tnf,6,1)) != r+1) return 0; |
|
if (typ(tnf[7]) != t_VEC || lg(tnf[7]) != 7) return 0; |
if (typ(tnf[7]) != t_VEC || lg(tnf[7]) != 7) return 0; |
return 1; |
return 1; |
} |
} |
|
|
static GEN distoZ(GEN z) |
static GEN |
|
distoZ(GEN z) |
{ |
{ |
GEN p1=gfrac(z); |
GEN p1 = gfrac(z); |
return gmin(p1, gsub(gun,p1)); |
return gmin(p1, gsub(gun,p1)); |
} |
} |
|
|
/* Check whether a solution has already been found */ |
/* Compensates rounding errors for computation/display of the constants. |
static int |
* Round up if dir > 0, down otherwise */ |
_thue_new(GEN zz) |
|
{ |
|
int i; |
|
for (i=1; i<lg(SOL); i++) |
|
if (gegal(zz,(GEN)SOL[i])) return 0; |
|
return 1; |
|
} |
|
|
|
/* Compensates rounding errors for computation/display of the constants */ |
|
static GEN |
static GEN |
myround(GEN Cst, GEN upd) |
myround(GEN x, long dir) |
{ |
{ |
return gmul(Cst, gadd(gun, gmul(upd,gpuigs(stoi(10),-10)))); |
GEN eps = gpowgs(stoi(dir > 0? 10: -10), -10); |
|
return gmul(x, gadd(gun, eps)); |
} |
} |
|
|
/* Returns the index of the largest element in a vector */ |
/* Returns the index of the largest element in a vector */ |
static int |
static GEN |
Vecmax(GEN Vec, int size) |
_Vecmax(GEN Vec, int *ind) |
{ |
{ |
int k, tno = 1; |
int k, tno = 1, l = lg(Vec); |
GEN tmax = (GEN)Vec[1]; |
GEN tmax = (GEN)Vec[1]; |
for (k=2; k<=size; k++) |
for (k = 2; k < l; k++) |
if (gcmp((GEN)Vec[k],tmax)==1) { tmax=(GEN)Vec[k]; tno=k; } |
if (gcmp((GEN)Vec[k],tmax) > 0) { tmax = (GEN)Vec[k]; tno = k; } |
return tno; |
if (ind) *ind = tno; |
|
return tmax; |
} |
} |
|
|
/* Performs basic computations concerning the equation: buchinitfu, c1, c2 */ |
static GEN |
static void |
Vecmax(GEN v) { return _Vecmax(v, NULL); } |
inithue(GEN poly, long flag) |
|
{ |
|
GEN roo2, tmp, gpmin, de; |
|
int k,j; |
|
|
|
x0=gzero; x1=gzero; deg=degpol(poly); gdeg=stoi(deg); |
static int |
|
Vecmaxind(GEN v) { int i; (void)_Vecmax(v, &i); return i; } |
|
|
if (!uftot) |
static GEN |
{ |
tnf_get_roots(GEN poly, long prec, int S, int T) |
uftot=bnfinit0(poly,1,NULL,Prec); |
{ |
if (uftot != checkbnf_discard(uftot)) |
GEN R0 = roots(poly, prec), R = cgetg(lg(R0), t_COL); |
err(talker,"non-monic polynomial in thue"); |
int k; |
if (flag) certifybuchall(uftot); |
|
s=itos(gmael3(uftot,7,2,1)); |
|
t=itos(gmael3(uftot,7,2,2)); |
|
} |
|
/* Switch the roots to get the usual order */ |
|
roo=roots(poly, Prec); roo2=cgetg(deg+1,t_COL); |
|
for (k=1; k<=s; k++) roo2[k]=roo[k]; |
|
for (k=1; k<=t; k++) |
|
{ |
|
roo2[k+s]=roo[2*k-1+s]; |
|
roo2[k+s+t]=lconj((GEN)roo2[k+s]); |
|
} |
|
roo=roo2; |
|
|
|
r=s+t-1; halpha=gun; |
for (k=1; k<=S; k++) R[k] = lreal((GEN)R0[k]); |
for (k=1; k<=deg; k++) |
/* swap roots to get the usual order */ |
halpha = gmul(halpha,gmax(gun,gabs((GEN)roo[k],Prec))); |
for (k=1; k<=T; k++) |
halpha=gdiv(glog(halpha,Prec),gdeg); |
|
|
|
de=derivpol(poly); c1=gabs(poleval(de,(GEN)roo[1]),Prec); |
|
for (k=2; k<=s; k++) |
|
{ |
{ |
tmp=gabs(poleval(de,(GEN)roo[k]),Prec); |
R[k+S] = R0[2*k+S-1]; |
if (gcmp(tmp,c1)==-1) c1=tmp; |
R[k+S+T]= R0[2*k+S]; |
} |
} |
|
return R; |
|
} |
|
|
c1=gdiv(gpui(gdeux,gsub(gdeg,gun),Prec),c1); c1=myround(c1,gun); |
/* Computation of the logarithmic height of x (given by embeddings) */ |
c2=gabs(gsub((GEN)roo[1],(GEN)roo[2]),Prec); |
static GEN |
|
LogHeight(GEN x, long prec) |
for (k=1; k<=deg; k++) |
{ |
for (j=k+1; j<=deg; j++) |
int i, n = lg(x)-1; |
{ |
GEN LH = gun; |
tmp=gabs(gsub((GEN)roo[j],(GEN)roo[k]),Prec); |
for (i=1; i<=n; i++) LH = gmul(LH, gmax(gun, gabs((GEN)x[i], prec))); |
if (gcmp(c2,tmp)==1) c2=tmp; |
return gdivgs(glog(LH,prec), n); |
} |
|
|
|
c2=myround(c2,stoi(-1)); |
|
if (t==0) x0=gun; |
|
else |
|
{ |
|
gpmin=gabs(poleval(de,(GEN)roo[s+1]),Prec); |
|
for (k=2; k<=t; k++) |
|
{ |
|
tmp=gabs(poleval(de,(GEN)roo[s+k]),Prec); |
|
if (gcmp(tmp,gpmin)==-1) gpmin=tmp; |
|
} |
|
|
|
/* Compute x0. See paper, Prop. 2.2.1 */ |
|
x0=gpui(gdiv(gpui(gdeux,gsub(gdeg,gun),Prec), |
|
gmul(gpmin, |
|
gabs((GEN)gimag(roo)[Vecmax(gabs(gimag(roo),Prec),deg)],Prec))), |
|
ginv(gdeg),Prec); |
|
} |
|
|
|
if (DEBUGLEVEL >=2) fprintferr("c1 = %Z\nc2 = %Z\n", c1, c2); |
|
|
|
} |
} |
|
|
/* Computation of M, its inverse A and check of the precision (see paper) */ |
/* x^(1/n) */ |
static void |
static GEN |
T_A_Matrices(void) |
mpsqrtn(GEN x, long n) |
{ |
{ |
GEN mask, eps1, eps2, nia, m1, IntM; |
if (typ(x) != t_REAL) err(typeer,"mpsqrtn"); |
int i,j; |
return mpexp(divrs(mplog(x), n)); |
|
} |
|
|
m1=glog(gabs(MatFU,Prec),Prec); mask=gsub(gpui(gdeux,stoi(r),Prec),gun); |
/* |x|^(1/n), x t_INT */ |
m1=matextract(m1,mask,mask); |
static GEN |
|
absisqrtn(GEN x, long n, long prec) { |
A=invmat(m1); IntM=gsub(gmul(A,m1),idmat(r)); |
GEN r = itor(x,prec); |
|
if (signe(r) < 0) r = negr(r); |
eps2=gzero; |
return mpsqrtn(r, n); |
for (i=1; i<=r; i++) |
|
for (j=1; j<=r; j++) |
|
if (gcmp(eps2,gabs(gcoeff(IntM,i,j),20)) == -1) |
|
eps2=gabs(gcoeff(IntM,i,j),20); |
|
|
|
eps1=gpui(gdeux,stoi((Prec-2)*32),Prec); eps2=gadd(eps2,ginv(eps1)); |
|
|
|
nia=gzero; |
|
for (i=1; i<=r; i++) |
|
for (j=1; j<=r; j++) |
|
if (gcmp(nia,gabs(gcoeff(A,i,j),20)) == -1) |
|
nia = gabs(gcoeff(A,i,j),20); |
|
|
|
/* Check for the precision in matrix inversion. See paper, Lemma 2.4.2. */ |
|
if (gcmp(gmul(gadd(gmul(stoi(r),gmul(nia,eps1)),eps2), |
|
gmul(gdeux,stoi(r))),gun) == -1) |
|
err(precer,"thue"); |
|
|
|
eps3=gmul(gdeux,gmul(gmul(gsqr(stoi(r)),nia), |
|
gadd(gmul(stoi(r),gdiv(nia,eps1)),eps2))); |
|
myround(eps3,gun); |
|
|
|
if (DEBUGLEVEL>=2) fprintferr("epsilon_3 -> %Z\n",eps3); |
|
} |
} |
|
|
/* Computation of the constants c5, c7, c10, c12, c15 */ |
static GEN |
static void |
get_emb(GEN x, GEN r) |
ComputeConstants(void) |
|
{ |
{ |
int k; |
int l = lg(r), i, tx; |
|
GEN e, y = cgetg(l, t_COL); |
|
|
GEN Vect; |
tx = typ(x); |
|
if (tx != t_INT && tx != t_POL) err(typeer,"get_emb"); |
Vect=cgetg(r+1,t_COL); for (k=1; k<=r; k++) Vect[k]=un; |
for (i=1; i<l; i++) |
if (numroot <= r) Vect[numroot]=lsub(gun,gdeg); |
|
|
|
Delta=gmul(A,Vect); |
|
|
|
c5=(GEN)(gabs(Delta,Prec)[Vecmax(gabs(Delta,Prec),r)]); |
|
c5=myround(c5,gun); c7=mulsr(r,c5); |
|
c10=divsr(deg,c7); c13=divsr(deg,c5); |
|
|
|
|
|
if (DEBUGLEVEL>=2) |
|
{ |
{ |
fprintferr("c5 = %Z\n",c5); |
e = poleval(x, (GEN)r[i]); |
fprintferr("c7 = %Z\n",c7); |
if (gcmp0(e)) return NULL; |
fprintferr("c10 = %Z\n",c10); |
y[i] = (long)e; |
fprintferr("c13 = %Z\n",c13); |
|
} |
} |
|
return y; |
} |
} |
|
|
/* Computation of the constants c6, c8, c9 */ |
/* Computation of the conjugates (sigma_i(v_j)), and log. heights, of elts of v */ |
static void |
static GEN |
ComputeConstants2(GEN poly, GEN rhs) |
Conj_LH(GEN v, GEN *H, GEN r, long prec) |
{ |
{ |
GEN Vect1, tmp; |
int j, l = lg(v); |
int k; |
GEN e, M = (GEN)cgetg(l,t_MAT); |
|
|
Vect1=cgetg(r+1,t_COL); |
(*H) = cgetg(l,t_COL); |
for (k=1; k<=r; k++) Vect1[k]=un; |
for (j = 1; j < l; j++) |
Vect1=gmul(gabs(A,ConstPrec),Vect1); |
{ |
|
if (! (e = get_emb((GEN)v[j], r)) ) return NULL; /* FAIL */ |
|
M[j] = (long)e; |
|
(*H)[j]= (long)LogHeight(e, prec); |
|
} |
|
return M; |
|
} |
|
|
|
static GEN abslog(GEN x, long prec) { return gabs(glog(x,prec), prec); } |
|
static GEN logabs(GEN x, long prec) { return glog(gabs(x,prec), prec); } |
|
|
Vect2=cgetg(r+1,t_COL); |
/* Computation of M, its inverse A and precision check (see paper) */ |
for (k=1; k<=r; k++) |
static GEN |
{ if (k!=numroot) |
T_A_Matrices(GEN MatFU, int r, GEN *eps5, long prec) |
{ Vect2[k]=llog(gabs(gdiv(gsub((GEN)roo[numroot],(GEN)roo[k]), |
{ |
gcoeff(MatNE,k,curne)),Prec),Prec); } |
GEN A, p1, m1, IntM, nia, eps3, eps2, eps1 = shifti(gun, bit_accuracy(prec)); |
else { Vect2[k]=llog(gabs(gdiv(rhs,gmul(poleval(derivpol(poly) |
|
,(GEN)roo[numroot]), |
|
gcoeff(MatNE,k,curne))),Prec),Prec); |
|
} |
|
} |
|
Lambda=gmul(A,Vect2); |
|
|
|
tmp=(GEN)Vect1[Vecmax(Vect1,r)]; |
m1 = rowextract_i(vecextract_i(MatFU, 1,r), 1,r); /* minor order r */ |
x2=gmax(x1,gpui(mulsr(10,mulrr(c4,tmp)),ginv(gdeg),ConstPrec)); |
m1 = logabs(m1,prec); |
c14=mulrr(c4,tmp); |
|
|
|
c6=gabs((GEN)Lambda[Vecmax(gabs(Lambda,ConstPrec),r)],ConstPrec); |
A = invmat(m1); |
c6=addrr(c6,dbltor(0.1)); c6=myround(c6,gun); |
IntM = gsub(gmul(A,m1), idmat(r)); |
|
|
c8=addrr(dbltor(1.23),mulsr(r,c6)); |
eps2 = gadd(vecmax(gabs(IntM,prec)), ginv(eps1)); |
c11=mulrr(mulsr(2,c3),gexp(divrr(mulsr(deg,c8),c7),ConstPrec)); |
nia = vecmax(gabs(A,prec)); |
|
|
tmp=gexp(divrr(mulsr(deg,c6),c5),ConstPrec); |
/* Check for the precision in matrix inversion. See paper, Lemma 2.4.2. */ |
c12=mulrr(mulsr(2,c3),tmp); |
p1 = gadd(gmulsg(r, gmul(nia,eps1)), eps2); |
c15=mulsr(2,mulrr(c14,tmp)); |
if (gcmp(gmulgs(p1, 2*r), gun) < 0) err(precer,"thue"); |
|
|
if (DEBUGLEVEL>=2) |
p1 = gadd(gmulsg(r, gdiv(nia,eps1)), eps2); |
{ |
eps3 = gmul(gmulsg(2*r*r,nia), p1); |
fprintferr("c6 = %Z\n",c6); |
eps3 = myround(eps3, 1); |
fprintferr("c8 = %Z\n",c8); |
|
fprintferr("c11 = %Z\n",c11); |
if (DEBUGLEVEL>1) fprintferr("epsilon_3 -> %Z\n",eps3); |
fprintferr("c12 = %Z\n",c12); |
*eps5 = mulsr(r, eps3); return A; |
fprintferr("c14 = %Z\n",c14); |
|
fprintferr("c15 = %Z\n",c15); |
|
} |
|
} |
} |
|
|
/* Computation of the logarithmic heights of the units */ |
/* Performs basic computations concerning the equation. |
|
* Returns a "tnf" structure containing |
|
* 1) the polynomial |
|
* 2) the bnf (used to solve the norm equation) |
|
* 3) roots, with presumably enough precision |
|
* 4) The logarithmic heights of units |
|
* 5) The matrix of conjugates of units |
|
* 6) its inverse |
|
* 7) a few technical constants */ |
static GEN |
static GEN |
Logarithmic_Height(int s) |
inithue(GEN P, GEN bnf, long flag, long prec) |
{ |
{ |
int i; |
GEN MatFU, x0, tnf, tmp, gpmin, dP, csts, ALH, eps5, ro, c1, c2; |
GEN LH=gun,mat; |
int k,j, s,t, n = degpol(P); |
|
|
mat=gabs(MatFU,Prec); |
if (!bnf) |
for (i=1; i<=deg; i++) |
{ |
LH = gmul(LH,gmax(gun,gabs(gcoeff(mat,i,s),Prec))); |
bnf = bnfinit0(P,1,NULL,prec); |
return gmul(gdeux,gdiv(glog(LH,Prec),gdeg)); |
if (bnf != checkbnf_discard(bnf)) err(talker,"non-monic polynomial in thue"); |
} |
if (flag) (void)certifybuchall(bnf); |
|
} |
|
nf_get_sign(checknf(bnf), (long*)&s, (long*)&t); |
|
ro = tnf_get_roots(P, prec, s, t); |
|
MatFU = Conj_LH(gmael(bnf,8,5), &ALH, ro, prec); |
|
if (!MatFU) return NULL; /* FAIL */ |
|
|
/* Computation of the matrix ((\sigma_i(\eta_j)) used for M_1 and |
dP = derivpol(P); |
the computation of logarithmic heights */ |
c1 = NULL; /* min |P'(r_i)|, i <= s */ |
static void |
for (k=1; k<=s; k++) |
Compute_Fund_Units(GEN uf) |
|
{ |
|
int i,j; |
|
MatFU=cgetg(r+1,t_MAT); |
|
|
|
for (i=1; i<=r; i++) |
|
MatFU[i]=lgetg(deg+1,t_COL); |
|
for (i=1; i<=r; i++) |
|
{ |
{ |
if (typ(uf[i])!=t_POL) err(talker,"incorrect system of units"); |
tmp = gabs(poleval(dP,(GEN)ro[k]),prec); |
for (j=1; j<=deg; j++) |
if (!c1 || gcmp(tmp,c1) < 0) c1 = tmp; |
coeff(MatFU,j,i)=(long)poleval((GEN)uf[i],(GEN)roo[j]); |
|
} |
} |
} |
c1 = gdiv(shifti(gun,n-1), c1); |
|
c1 = gprec_w(myround(c1, 1), DEFAULTPREC); |
|
|
/* Computation of the conjugates of the solutions to the norm equation */ |
c2 = NULL; /* max |r_i - r_j|, i!=j */ |
static void |
for (k=1; k<=n; k++) |
Conj_Norm_Eq(GEN ne, GEN *Hmu) |
for (j=k+1; j<=n; j++) |
{ |
{ |
int p,q,nesol,x; |
tmp = gabs(gsub((GEN)ro[j],(GEN)ro[k]), prec); |
|
if (!c2 || gcmp(c2,tmp) > 0) c2 = tmp; |
|
} |
|
c2 = gprec_w(myround(c2, -1), DEFAULTPREC); |
|
|
nesol=lg(ne); MatNE=(GEN)cgetg(nesol,t_MAT); (*Hmu)=cgetg(nesol,t_COL); |
if (t==0) x0 = gun; |
|
else |
for (p=1; p<nesol; p++) { MatNE[p]=lgetg(deg+1,t_COL); (*Hmu)[p]=un; } |
|
for (q=1; q<nesol; q++) |
|
{ |
{ |
x=typ(ne[q]); |
gpmin = NULL; /* min |P'(r_i)|, i > s */ |
if (x!=t_INT && x!=t_POL) |
for (k=1; k<=t; k++) |
err(talker,"incorrect solutions of norm equation"); |
|
for (p=1; p<=deg; p++) |
|
{ |
{ |
coeff(MatNE,p,q)=(long)poleval((GEN)ne[q],(GEN)roo[p]); |
tmp = gabs(poleval(dP,(GEN)ro[s+k]), prec); |
/* Logarithmic height of the solutions of the norm equation */ |
if (!gpmin || gcmp(tmp,gpmin) < 0) gpmin = tmp; |
(*Hmu)[q]=lmul((GEN)(*Hmu)[q],gmax(gun,gabs(gcoeff(MatNE,p,q),Prec))); |
|
} |
} |
|
gpmin = gprec_w(gpmin, DEFAULTPREC); |
|
|
|
/* Compute x0. See paper, Prop. 2.2.1 */ |
|
x0 = gmul(gpmin, Vecmax(gabs(gimag(ro), prec))); |
|
x0 = mpsqrtn(gdiv(shifti(gun,n-1), x0), n); |
} |
} |
for (q=1; q<nesol; q++) |
if (DEBUGLEVEL>1) fprintferr("c1 = %Z\nc2 = %Z\n", c1, c2); |
(*Hmu)[q]=ldiv(glog((GEN)(*Hmu)[q],Prec),gdeg); |
|
|
ALH = gmul2n(ALH, 1); |
|
tnf = cgetg(8,t_VEC); csts = cgetg(7,t_VEC); |
|
tnf[1] = (long)P; |
|
tnf[2] = (long)bnf; |
|
tnf[3] = (long)ro; |
|
tnf[4] = (long)ALH; |
|
tnf[5] = (long)MatFU; |
|
tnf[6] = (long)T_A_Matrices(MatFU, s+t-1, &eps5, prec); |
|
tnf[7] = (long)csts; |
|
csts[1] = (long)c1; csts[2] = (long)c2; csts[3] = (long)LogHeight(ro, prec); |
|
csts[4] = (long)x0; csts[5] = (long)eps5; csts[6] = (long)stoi(prec); |
|
return tnf; |
} |
} |
|
|
/* Compute Baker's bound c11, and B_0, the bound for the b_i's .*/ |
typedef struct { |
static void |
GEN c10, c11, c13, c15, bak, NE, ALH, hal, MatFU, ro, Hmu; |
Baker(GEN ALH, GEN Hmu) |
GEN delta, lambda, errdelta; |
|
int r, iroot; |
|
} baker_s; |
|
|
|
/* Compute Baker's bound c9 and B_0, the bound for the b_i's. See Thm 2.3.1 */ |
|
static GEN |
|
Baker(baker_s *BS) |
{ |
{ |
GEN c9=gun, gbak, hb0=gzero; |
const long prec = DEFAULTPREC; |
int k,i1,i2; |
GEN tmp, B0, hb0, c9 = gun, ro = BS->ro, ro0 = (GEN)ro[BS->iroot]; |
|
int k, i1, i2, r = BS->r; |
|
|
gbak = gmul(gmul(gdeg,gsub(gdeg,gun)),gsub(gdeg,gdeux)); |
switch (BS->iroot) { |
|
case 1: i1=2; i2=3; break; |
switch (numroot) { |
case 2: i1=1; i2=3; break; |
case 1: i1=2; i2=3; break; |
default: i1=1; i2=2; break; |
case 2: i1=1; i2=3; break; |
|
case 3: i1=1; i2=2; break; |
|
default: i1=1; i2=2; break; |
|
} |
} |
|
|
/* For the symbols used for the computation of c11, see paper, Thm 2.3.1 */ |
|
/* Compute h_1....h_r */ |
/* Compute h_1....h_r */ |
for (k=1; k<=r; k++) |
for (k=1; k<=r; k++) |
{ |
{ |
c9=gmul(c9,gmax((GEN)ALH[k], |
tmp = gdiv(gcoeff(BS->MatFU,i1,k), gcoeff(BS->MatFU,i2,k)); |
gmax(ginv(gbak), |
tmp = gmax(gun, abslog(tmp,prec)); |
gdiv(gabs(glog(gdiv(gcoeff(MatFU,i1,k), |
c9 = gmul(c9, gmax((GEN)BS->ALH[k], gdiv(tmp, BS->bak))); |
gcoeff(MatFU,i2,k)), |
|
Prec),Prec),gbak)))); |
|
} |
} |
|
|
/* Compute a bound for the h_0 */ |
/* Compute a bound for the h_0 */ |
hb0=gadd(gmul(stoi(4),halpha),gadd(gmul(gdeux,(GEN)Hmu[curne]), |
hb0 = gadd(gmul2n(BS->hal,2), gmul(gdeux, gadd(BS->Hmu,mplog2(prec)))); |
gmul(gdeux,glog(gdeux,Prec)))); |
tmp = gdiv(gmul(gsub(ro0, (GEN)ro[i2]), (GEN)BS->NE[i1]), |
hb0=gmax(hb0,gmax(ginv(gbak), |
gmul(gsub(ro0, (GEN)ro[i1]), (GEN)BS->NE[i2])); |
gdiv(gabs(glog(gdiv(gmul(gsub((GEN)roo[numroot], |
tmp = gmax(gun, abslog(tmp, prec)); |
(GEN)roo[i2]), |
hb0 = gmax(hb0, gdiv(tmp, BS->bak)); |
gcoeff(MatNE,i1,curne)), |
c9 = gmul(c9,hb0); |
gmul(gsub((GEN)roo[numroot], |
|
(GEN)roo[i1]), |
|
gcoeff(MatNE,i2,curne))), |
|
Prec),Prec),gbak))); |
|
c9=gmul(c9,hb0); |
|
/* Multiply c9 by the "constant" factor */ |
/* Multiply c9 by the "constant" factor */ |
c9=gmul(c9,gmul(gmul(gmul(stoi(18),mppi(Prec)),gpui(stoi(32),stoi(4+r),Prec)), |
c9 = gmul(c9, gmul(mulri(mulsr(18,mppi(prec)), gpowgs(stoi(32),4+r)), |
gmul(gmul(mpfact(r+3), gpowgs(gmul(gbak,stoi(r+2)), r+3)), |
gmul(gmul(mpfact(r+3), gpowgs(mulis(BS->bak,r+2), r+3)), |
glog(gmul(gdeux,gmul(gbak,stoi(r+2))),Prec)))); |
glog(mulis(BS->bak,2*(r+2)),prec)))); |
c9=myround(c9,gun); |
c9 = gprec_w(myround(c9, 1), DEFAULTPREC); |
/* Compute B0 according to Lemma 2.3.3 */ |
/* Compute B0 according to Lemma 2.3.3 */ |
B0=gmax(gexp(gun,Prec), |
B0 = mulsr(2, divrr(addrr(mulrr(c9,mplog(divrr(c9,BS->c10))), mplog(BS->c11)), |
mulsr(2,divrr(addrr(mulrr(c9,glog(divrr(c9,c10),ConstPrec)), |
BS->c10)); |
glog(c11,ConstPrec)),c10))); |
B0 = gmax(B0, dbltor(2.71828183)); |
|
|
if (DEBUGLEVEL>=2) fprintferr("Baker -> %Z\nB0 -> %Z\n",c9,B0); |
if (DEBUGLEVEL>1) { |
|
fprintferr(" B0 = %Z\n",B0); |
|
fprintferr(" Baker = %Z\n",c9); |
|
} |
|
return B0; |
} |
} |
|
|
/* Compute delta and lambda */ |
/* || x d ||, x t_REAL, d t_INT */ |
static void |
static GEN |
Create_CF_Values(int i1, int i2, GEN *errdelta) |
errnum(GEN x, GEN d) |
{ |
{ |
GEN eps5; |
GEN dx = mulir(d, x); |
|
return mpabs(subri(dx, ground(dx))); |
if (r>1) |
|
{ |
|
delta=divrr((GEN)Delta[i2],(GEN)Delta[i1]); |
|
eps5=mulrs(eps3,r); |
|
*errdelta=mulrr(addsr(1,delta), |
|
divrr(eps5,subrr(gabs((GEN)Delta[i1],ConstPrec),eps5))); |
|
|
|
lambda=gdiv(gsub(gmul((GEN)Delta[i2],(GEN)Lambda[i1]), |
|
gmul((GEN)Delta[i1],(GEN)Lambda[i2])), |
|
(GEN)Delta[i1]); |
|
} |
|
else |
|
{ |
|
GEN Pi2 = gmul2n(mppi(Prec),1); |
|
delta=gdiv(gcoeff(MatFU,2,1),gcoeff(MatFU,3,1)); |
|
delta=gdiv(garg(delta,Prec),Pi2); |
|
*errdelta=gdiv(gpui(gdeux,stoi(-(Prec-2)*32+1),Prec), |
|
gabs(gcoeff(MatFU,2,1),Prec)); |
|
lambda=gmul(gdiv(gsub((GEN)roo[numroot],(GEN)roo[2]), |
|
gsub((GEN)roo[numroot],(GEN)roo[3])), |
|
gdiv(gcoeff(MatNE,3,curne),gcoeff(MatNE,2,curne))); |
|
lambda=gdiv(garg(lambda,Prec),Pi2); |
|
} |
|
if (DEBUGLEVEL>=2) outerr(*errdelta); |
|
} |
} |
|
|
/* Try to reduce the bound through continued fractions; see paper. */ |
/* Try to reduce the bound through continued fractions; see paper. */ |
static int |
static int |
CF_First_Pass(GEN kappa, GEN errdelta) |
CF_1stPass(GEN *B0, GEN kappa, baker_s *BS) |
{ |
{ |
GEN q,ql,qd,l0; |
GEN q, ql, qd, l0, denbound = mulri(*B0, kappa); |
|
|
if (gcmp(gmul(dbltor(0.1),gsqr(mulri(B0,kappa))),ginv(errdelta))==1) |
if (gcmp(gmul(dbltor(0.1),gsqr(denbound)), ginv(BS->errdelta)) > 0) return -1; |
{ |
|
return -1; |
|
} |
|
|
|
q=denom(bestappr(delta, mulir(kappa, B0))); ql=mulir(q,lambda); |
q = denom( bestappr(BS->delta, denbound) ); |
qd=gmul(q,delta); ql=gabs(subri(ql, ground(ql)),Prec); |
qd = errnum(BS->delta, q); |
qd=gabs(mulrr(subri(qd,ground(qd)),B0),Prec); |
ql = errnum(BS->lambda,q); |
l0=subrr(ql,addrr(qd,divri(dbltor(0.1),kappa))); |
|
if (signe(l0) <= 0) |
l0 = subrr(ql, addrr(mulrr(qd, *B0), divri(dbltor(0.1),kappa))); |
|
if (signe(l0) <= 0) return 0; |
|
|
|
if (BS->r > 1) |
|
*B0 = divrr(mplog(divrr(mulir(q,BS->c15), l0)), BS->c13); |
|
else |
{ |
{ |
if (DEBUGLEVEL>=2) |
l0 = mulrr(l0,Pi2n(1, DEFAULTPREC)); |
fprintferr("CF_First_Pass failed. Trying again with larger kappa\n"); |
*B0 = divrr(mplog(divrr(mulir(q,BS->c11), l0)), BS->c10); |
return 0; |
|
} |
} |
|
if (DEBUGLEVEL>1) fprintferr(" B0 -> %Z\n",*B0); |
|
return 1; |
|
} |
|
|
if (r>1) |
/* Check whether a solution has already been found */ |
B0=divrr(glog(divrr(mulir(q,c15),l0),ConstPrec),c13); |
static int |
else |
new_sol(GEN z, GEN S) |
B0=divrr(glog(divrr(mulir(q,c11),mulrr(l0,gmul2n(mppi(ConstPrec),1))),ConstPrec),c10); |
{ |
|
int i, l = lg(S); |
|
for (i=1; i<l; i++) |
|
if (gegal(z,(GEN)S[i])) return 0; |
|
return 1; |
|
} |
|
|
if (DEBUGLEVEL>=2) |
static void |
fprintferr("CF_First_Pass successful !!\nB0 -> %Z\n",B0); |
add_sol(GEN *pS, GEN x, GEN y) |
return 1; |
{ |
} |
GEN u = cgetg(3,t_VEC); u[1] = (long)x; u[2] = (long)y; |
|
if (new_sol(u, *pS)) *pS = concatsp(*pS, _vec(u)); |
|
} |
|
|
/* Check whether a "potential" solution is truly a solution. */ |
|
static void |
static void |
Check_Solutions(GEN xmay1, GEN xmay2, GEN poly, GEN rhs) |
check_sol(GEN x, GEN y, GEN P, GEN rhs, GEN *pS) |
{ |
{ |
GEN xx1, xx2, yy1, yy2, zz, u; |
if (gcmp0(y)) |
|
{ if (gegal(gpowgs(x,degpol(P)), rhs)) add_sol(pS, x, gzero); } |
|
else |
|
{ if (gegal(poleval(rescale_pol(P,y),x), rhs)) add_sol(pS, x, y); } |
|
} |
|
|
yy1=ground(greal(gdiv(gsub(xmay2,xmay1), gsub((GEN)roo[1],(GEN)roo[2])))); |
/* Check whether a potential solution is a true solution. Return 0 if |
yy2=gneg(yy1); |
* truncation error (increase precision) */ |
|
static int |
|
CheckSol(GEN *pS, GEN z1, GEN z2, GEN P, GEN rhs, GEN ro) |
|
{ |
|
GEN x, y, ro1 = (GEN)ro[1], ro2 = (GEN)ro[2]; |
|
long e; |
|
|
xx1=greal(gadd(xmay1,gmul((GEN)roo[1],yy1))); |
y = grndtoi(greal(gdiv(gsub(z2,z1), gsub(ro1,ro2))), &e); |
xx2=gneg(xx1); |
if (e > 0) return 0; |
|
x = gadd(z1, gmul(ro1, y)); |
if (gcmp(distoZ(xx1),dbltor(0.0001))==-1) |
x = grndtoi(greal(x), &e); |
|
if (e > 0) return 0; |
|
if (e <= -13) |
{ |
{ |
xx1=ground(xx1); |
check_sol( x , y, P,rhs,pS); |
if (!gcmp0(yy1)) |
check_sol(gneg(x), gneg(y), P,rhs,pS); |
{ |
|
if (gegal(gmul(poleval(poly,gdiv(xx1,yy1)), |
|
gpui(yy1,gdeg,Prec)),(GEN)rhs)) |
|
{ |
|
zz=cgetg(2,t_VEC); u=cgetg(3,t_VEC); |
|
u[1]=(long)xx1; u[2]=(long)yy1; zz[1]=(long)u; |
|
if (_thue_new(u)) SOL=concatsp(SOL,zz); |
|
} |
|
} |
|
else if (gegal(gpui(xx1,gdeg,Prec),(GEN)rhs)) |
|
{ |
|
zz=cgetg(2,t_VEC); u=cgetg(3,t_VEC); |
|
u[1]=(long)xx1; u[2]=(long)gzero; zz[1]=(long)u; |
|
if (_thue_new(u)) SOL=concatsp(SOL,zz); |
|
} |
|
|
|
xx2=ground(xx2); |
|
if (!gcmp0(yy2)) |
|
{ |
|
if (gegal(gmul(poleval(poly,gdiv(xx2,yy2)), |
|
gpui(yy2,gdeg,Prec)),(GEN)rhs)) |
|
{ |
|
zz=cgetg(2,t_VEC); u=cgetg(3,t_VEC); |
|
u[1]=(long)xx2; u[2]=(long)yy2; zz[1]=(long)u; |
|
if (_thue_new(u)) SOL=concatsp(SOL,zz); |
|
} |
|
} |
|
else if (gegal(gpui(xx2,gdeg,Prec),(GEN)rhs)) |
|
{ |
|
zz=cgetg(2,t_VEC); u=cgetg(3,t_VEC); |
|
u[1]=(long)xx2; u[2]=(long)gzero; zz[1]=(long)u; |
|
if (_thue_new(u)) SOL=concatsp(SOL,zz); |
|
} |
|
} |
} |
|
return 1; |
} |
} |
|
|
|
/* find q1,q2,q3 st q1 + b q2 + c q3 ~ 0 */ |
static GEN |
static GEN |
GuessQi(GEN *q1, GEN *q2, GEN *q3) |
GuessQi(GEN b, GEN c, GEN *eps) |
{ |
{ |
GEN C, Lat, eps; |
GEN Q, Lat, C = gpowgs(stoi(10), 10); |
|
|
C=gpui(stoi(10),stoi(10),Prec); |
Lat = idmat(3); |
Lat=idmat(3); mael(Lat,1,3)=(long)C; mael(Lat,2,3)=lround(gmul(C,delta)); |
mael(Lat,1,3) = (long)C; |
mael(Lat,3,3)=lround(gmul(C,lambda)); |
mael(Lat,2,3) = lround(gmul(C,b)); |
|
mael(Lat,3,3) = lround(gmul(C,c)); |
|
|
Lat=lllint(Lat); |
Q = (GEN)lllint(Lat)[1]; |
*q1=gmael(Lat,1,1); *q2=gmael(Lat,1,2); *q3=gmael(Lat,1,3); |
if (gcmp0((GEN)Q[3])) return NULL; /* FAIL */ |
|
|
eps=gabs(gadd(gadd(*q1,gmul(*q2,delta)),gmul(*q3,lambda)),Prec); |
*eps = gadd(gadd((GEN)Q[1], gmul((GEN)Q[2],b)), gmul((GEN)Q[3],c)); |
return(eps); |
*eps = mpabs(*eps); return Q; |
} |
} |
|
|
static void |
|
TotRat(void) |
|
{ |
|
err(bugparier,"thue (totally rational case)"); |
|
} |
|
|
|
/* Check for solutions under a small bound (see paper) */ |
/* Check for solutions under a small bound (see paper) */ |
static void |
static GEN |
Check_Small(int bound, GEN poly, GEN rhs) |
SmallSols(GEN S, int Bx, GEN poly, GEN rhs, GEN ro) |
{ |
{ |
GEN interm, xx, zz, u, maxr, tmp, ypot, xxn, xxnm1, yy; |
gpmem_t av = avma, lim = stack_lim(av, 1); |
long av = avma, lim = stack_lim(av,1); |
const long prec = DEFAULTPREC; |
int x, j, bsupy, y; |
GEN Hpoly, interm, X, Xn, Xnm1, Y, sqrtnRHS; |
|
int x, y, j, By, n = degpol(poly); |
double bndyx; |
double bndyx; |
|
|
maxr=gabs((GEN)roo[1],Prec); for(j=1; j<=deg; j++) |
if (DEBUGLEVEL>1) fprintferr("* Checking for small solutions\n"); |
{ if(gcmp(tmp=gabs((GEN)roo[j],Prec),maxr)==1) maxr=tmp; } |
|
|
sqrtnRHS = absisqrtn(rhs, n, prec); |
|
bndyx = gtodouble(gadd(sqrtnRHS, Vecmax(gabs(ro,prec)))); |
|
|
bndyx=gtodouble(gadd(gpui(gabs(rhs,Prec),ginv(gdeg),Prec),maxr)); |
/* x = 0 first */ |
|
Y = ground(sqrtnRHS); |
for (x=-bound; x<=bound; x++) |
if (gegal(gpowgs(Y,n), rhs)) add_sol(&S, Y, gzero); |
|
Y = negi(Y); |
|
if (gegal(gpowgs(Y,n), rhs)) add_sol(&S, Y, gzero); |
|
/* x != 0 */ |
|
for (x = -Bx; x <= Bx; x++) |
{ |
{ |
if (low_stack(lim,stack_lim(av,1))) |
if (!x) continue; |
|
|
|
X = stoi(x); Xn = gpowgs(X,n); |
|
interm = subii(rhs, mulii(Xn, (GEN)poly[2])); |
|
Xnm1 = Xn; j = 2; |
|
while (gcmp0(interm)) |
{ |
{ |
if(DEBUGMEM>1) err(warnmem,"Check_small"); |
Xnm1 = divis(Xnm1, x); |
SOL = gerepilecopy(av, SOL); |
interm = mulii((GEN)poly[++j], Xnm1); |
} |
} |
if (x==0) |
/* y | interm */ |
{ |
|
ypot=gmul(stoi(gsigne(rhs)),ground(gpui(gabs(rhs,0),ginv(gdeg),Prec))); |
|
if (gegal(gpui(ypot,gdeg,0),rhs)) |
|
{ |
|
zz=cgetg(2,t_VEC); u=cgetg(3,t_VEC); |
|
u[1]=(long)ypot; u[2]=(long)gzero; zz[1]=(long)u; |
|
if (_thue_new(u)) SOL=concatsp(SOL,zz); |
|
} |
|
if (gegal(gpui(gneg(ypot),gdeg,0),rhs)) |
|
{ |
|
zz=cgetg(2,t_VEC); u=cgetg(3,t_VEC); |
|
u[1]=(long)gneg(ypot); u[2]=(long)gzero; zz[1]=(long)u; |
|
if (_thue_new(u)) SOL=concatsp(SOL,zz); |
|
} |
|
} |
|
else |
|
{ |
|
bsupy=(int)(x>0?bndyx*x:-bndyx*x); |
|
|
|
xx=stoi(x); xxn=gpui(xx,gdeg,Prec); |
|
interm=gsub((GEN)rhs,gmul(xxn,(GEN)poly[2])); |
|
|
|
/* Verifier ... */ |
Hpoly = rescale_pol(poly,X); /* homogenize */ |
xxnm1=xxn; j=2; |
By = (int)(x>0? bndyx*x: -bndyx*x); |
while(gcmp0(interm)) |
if (gegal(gmul((GEN)poly[2],Xn),rhs)) add_sol(&S, gzero, X); /* y = 0 */ |
{ |
for(y = -By; y <= By; y++) |
xxnm1=gdiv(xxnm1,xx); |
{ |
interm=gmul((GEN)poly[++j],xxnm1); |
if (!y || smodis(interm, y)) continue; |
} |
Y = stoi(y); |
|
if (gegal(poleval(Hpoly, Y), rhs)) add_sol(&S, Y, X); |
|
} |
|
|
for(y=-bsupy; y<=bsupy; y++) |
if (low_stack(lim,stack_lim(av,1))) |
{ |
{ |
yy=stoi(y); |
if(DEBUGMEM>1) err(warnmem,"Check_small"); |
if(y==0) { |
S = gerepilecopy(av, S); |
if (gegal(gmul((GEN)poly[2],xxn),rhs)) |
} |
{ |
|
zz=cgetg(2,t_VEC); u=cgetg(3,t_VEC); |
|
u[1]=(long)gzero; u[2]=(long)xx; zz[1]=(long)u; |
|
if (_thue_new(u)) SOL=concatsp(SOL,zz); |
|
} |
|
} |
|
else if (gcmp0(gmod(interm,yy))) |
|
if(gegal(poleval(poly,gdiv(yy,xx)),gdiv(rhs,xxn))) |
|
/* Remplacer par un eval *homogene* */ |
|
{ |
|
zz=cgetg(2,t_VEC); u=cgetg(3,t_VEC); |
|
u[1]=(long)yy; u[2]=(long)xx; zz[1]=(long)u; |
|
if (_thue_new(u)) SOL=concatsp(SOL,zz); |
|
} |
|
} |
|
} |
|
} |
} |
|
return S; |
} |
} |
|
|
/* Computes [x]! */ |
/* Computes [x]! */ |
static double |
static double |
fact(double x) |
fact(double x) |
{ |
{ |
double ft=1.0; |
double ft = 1.0; |
x=(int)x; while (x>1) { ft*=x; x--; } |
x = (int)x; while (x>1) { ft *= x; x--; } |
return ft ; |
return ft ; |
} |
} |
|
|
/* From a polynomial and optionally a system of fundamental units for the |
/* From a polynomial and optionally a system of fundamental units for the |
* field defined by poly, computes all the relevant constants needed to solve |
* field defined by poly, computes all relevant constants needed to solve |
* the equation P(x,y)=a whenever the solutions of N_{ K/Q }(x)=a. Returns a |
* the equation P(x,y)=a given the solutions of N_{K/Q}(x)=a (see inithue). */ |
* "tnf" structure containing |
|
* 1) the polynomial |
|
* 2) the bnf (used to solve the norm equation) |
|
* 3) roots, with presumably enough precision |
|
* 4) The logarithmic heights of units |
|
* 5) The matrix of conjugates of units |
|
* 6) its inverse |
|
* 7) a few technical constants |
|
*/ |
|
GEN |
GEN |
thueinit(GEN poly, long flag, long prec) |
thueinit(GEN poly, long flag, long prec) |
{ |
{ |
GEN thueres,ALH,csts,c0; |
GEN tnf, bnf = NULL; |
ulong av = avma; |
gpmem_t av = avma; |
long k,st; |
long k, s; |
double d,dr; |
|
|
|
uftot = NULL; |
if (checktnf(poly)) { bnf = checkbnf((GEN)poly[2]); poly = (GEN)poly[1]; } |
if (checktnf(poly)) { uftot=(GEN)poly[2]; poly=(GEN)poly[1]; } |
if (typ(poly)!=t_POL) err(notpoler,"thueinit"); |
else |
if (degpol(poly) <= 2) err(talker,"invalid polynomial in thue (need deg>2)"); |
if (typ(poly)!=t_POL) err(notpoler,"thueinit"); |
|
if (degpol(poly)<=2) err(talker,"invalid polynomial in thue (need deg>2)"); |
|
|
|
if (!gisirreducible(poly)) err(redpoler,"thueinit"); |
if (!gisirreducible(poly)) err(redpoler,"thueinit"); |
st=sturm(poly); |
s = sturm(poly); |
if (st) |
if (s) |
{ |
{ |
dr=(double)((st+lgef(poly)-5)>>1); |
long PREC, n = degpol(poly); |
d=(double)degpol(poly); d=d*(d-1)*(d-2); |
double d, dr, dn = (double)n; |
/* Try to guess the precision by approximating Baker's bound. |
|
* Note that the guess is most of the time pretty generous, |
dr = (double)((s+n-2)>>1); /* s+t-1 */ |
* ie 10 to 30 decimal digits above what is *really* necessary. |
d = dn*(dn-1)*(dn-2); |
* Note that the limiting step is the reduction. See paper. |
/* Guess precision by approximating Baker's bound. The guess is most of |
*/ |
* the time not sharp, ie 10 to 30 decimal digits above what is _really_ |
Prec=3 + (long)((5.83 + (dr+4)*5 + log(fact(dr+3)) + (dr+3)*log(dr+2) + |
* necessary. Note that the limiting step is the reduction. See paper. */ |
|
PREC = 3 + (long)((5.83 + (dr+4)*5 + log(fact(dr+3)) + (dr+3)*log(dr+2) + |
(dr+3)*log(d) + log(log(2*d*(dr+2))) + (dr+1)) / 10.); |
(dr+3)*log(d) + log(log(2*d*(dr+2))) + (dr+1)) / 10.); |
ConstPrec=4; |
if (PREC < prec) PREC = prec; |
if (Prec<prec) Prec = prec; |
for (;;) |
if (!checktnf(poly)) inithue(poly,flag); |
{ |
|
if (( tnf = inithue(poly, bnf, flag, PREC) )) break; |
|
PREC = (PREC<<1)-2; |
|
if (DEBUGLEVEL>1) err(warnprec,"thueinit",PREC); |
|
bnf = NULL; avma = av; |
|
} |
|
} |
|
else |
|
{ |
|
GEN c0 = gun, ro = roots(poly, DEFAULTPREC); |
|
for (k=1; k<lg(ro); k++) c0 = gmul(c0, gimag((GEN)ro[k])); |
|
c0 = ginv( mpabs(c0) ); |
|
tnf = cgetg(3,t_VEC); |
|
tnf[1] = (long)poly; |
|
tnf[2] = (long)c0; |
|
} |
|
return gerepilecopy(av,tnf); |
|
} |
|
|
thueres=cgetg(8,t_VEC); |
static GEN |
thueres[1]=(long)poly; |
get_B0(int i1, GEN Delta, GEN Lambda, GEN eps5, long prec, baker_s *BS) |
thueres[2]=(long)uftot; |
{ |
thueres[3]=(long)roo; |
GEN delta, lambda, errdelta, B0 = Baker(BS); |
Compute_Fund_Units(gmael(uftot,8,5)); |
int i2, r = BS->r; |
ALH=cgetg(r+1,t_COL); |
|
for (k=1; k<=r; k++) ALH[k]=(long)Logarithmic_Height(k); |
|
thueres[4]=(long)ALH; |
|
thueres[5]=(long)MatFU; |
|
T_A_Matrices(); |
|
thueres[6]=(long)A; |
|
|
|
csts=cgetg(7,t_VEC); |
i2 = (i1 == 1)? 2: 1; |
csts[1]=(long)c1; csts[2]=(long)c2; csts[3]=(long)halpha; |
for(;;) /* i2 from 1 to r unless r = 1 [then i2 = 2] */ |
csts[4]=(long)x0; csts[5]=(long)eps3; |
{ |
csts[6]=(long)stoi(Prec); |
if (r > 1) |
|
{ |
|
delta = divrr((GEN)Delta[i2],(GEN)Delta[i1]); |
|
lambda = gdiv(gsub(gmul((GEN)Delta[i2],(GEN)Lambda[i1]), |
|
gmul((GEN)Delta[i1],(GEN)Lambda[i2])), |
|
(GEN)Delta[i1]); |
|
errdelta = mulrr(addsr(1,delta), |
|
divrr(eps5, subrr(mpabs((GEN)Delta[i1]),eps5))); |
|
} |
|
else |
|
{ /* r == 1, single fundamental unit (i1 = s = t = 1) */ |
|
GEN p1, Pi2 = Pi2n(1, prec); |
|
GEN fu = (GEN)BS->MatFU[1], ro = BS->ro; |
|
|
thueres[7]=(long)csts; |
p1 = gdiv((GEN)fu[2], (GEN)fu[3]); |
return gerepilecopy(av,thueres); |
delta = divrr(garg(p1,prec), Pi2); |
} |
|
|
|
thueres=cgetg(3,t_VEC); c0=gun; Prec=4; |
p1 = gmul(gdiv(gsub((GEN)ro[1], (GEN)ro[2]), |
roo=roots(poly,Prec); |
gsub((GEN)ro[1], (GEN)ro[3])), |
for (k=1; k<lg(roo); k++) |
gdiv((GEN)BS->NE[3], (GEN)BS->NE[2])); |
c0=gmul(c0, gimag((GEN)roo[k])); |
lambda = divrr(garg(p1,prec), Pi2); |
c0=ginv(gabs(c0,Prec)); |
|
thueres[1]=(long)poly; thueres[2]=(long)c0; |
errdelta = gdiv(gmul2n(gun, 1 - bit_accuracy(prec)), |
return gerepilecopy(av,thueres); |
gabs((GEN)fu[2],prec)); |
|
} |
|
BS->delta = delta; |
|
BS->lambda= lambda; BS->errdelta = errdelta; |
|
if (DEBUGLEVEL>1) fprintferr(" errdelta = %Z\n",errdelta); |
|
|
|
if (DEBUGLEVEL>1) fprintferr(" Entering CF...\n"); |
|
/* Reduce B0 as long as we make progress: newB0 < oldB0 - 0.1 */ |
|
for (;;) |
|
{ |
|
GEN oldB0 = B0, kappa = stoi(10); |
|
int cf; |
|
|
|
for (cf = 0; cf < 10; cf++, kappa = mulis(kappa,10)) |
|
{ |
|
int res = CF_1stPass(&B0, kappa, BS); |
|
if (res < 0) return NULL; /* prec problem */ |
|
if (res) break; |
|
if (DEBUGLEVEL>1) fprintferr("CF failed. Increasing kappa\n"); |
|
} |
|
if (cf == 10) |
|
{ /* Semirational or totally rational case */ |
|
GEN Q, ep, q, l0, denbound; |
|
|
|
if (! (Q = GuessQi(delta, lambda, &ep)) ) break; |
|
|
|
denbound = gadd(B0, absi((GEN)Q[2])); |
|
q = denom( bestappr(delta, denbound) ); |
|
l0 = subrr(errnum(delta, q), ep); |
|
if (signe(l0) <= 0) break; |
|
|
|
B0 = divrr(mplog(divrr(mulir((GEN)Q[3], BS->c15), l0)), BS->c13); |
|
if (DEBUGLEVEL>1) fprintferr("Semirat. reduction: B0 -> %Z\n",B0); |
|
} |
|
/* if no progress, stop */ |
|
if (gcmp(oldB0, gadd(B0,dbltor(0.1))) <= 0) return gmin(oldB0, B0); |
|
} |
|
i2++; if (i2 == i1) i2++; |
|
if (i2 > r) break; |
|
} |
|
err(bugparier,"thue (totally rational case)"); |
|
return NULL; /* not reached */ |
} |
} |
|
|
/* Given a tnf structure as returned by thueinit, and a right-hand-side and |
static GEN |
* optionally the solutions to the norm equation, returns the solutions to |
LargeSols(GEN tnf, GEN rhs, GEN ne, GEN *pro, GEN *pS) |
* the Thue equation F(x,y)=a |
|
*/ |
|
GEN |
|
thue(GEN thueres, GEN rhs, GEN ne) |
|
{ |
{ |
GEN Kstart,Old_B0,ALH,errdelta,Hmu,c0,poly,csts,bd; |
GEN Vect, P, ro, bnf, MatFU, A, csts, dP, vecdP; |
GEN xmay1,xmay2,b,zp1,tmp,q1,q2,q3,ep; |
GEN c1,c2,c3,c4,c10,c11,c13,c14,c15, x0, x1, x2, x3, b, zp1, tmp, eps5; |
long Nb_It=0,upb,bi1,i1,i2,i, flag,cf,fs; |
int iroot, ine, n, i, r,s,t; |
ulong av = avma; |
long upb, bi1, Prec, prec; |
|
baker_s BS; |
|
gpmem_t av = avma; |
|
|
if (!checktnf(thueres)) err(talker,"not a tnf in thue"); |
bnf = (GEN)tnf[2]; |
|
if (!ne) ne = bnfisintnorm(bnf, rhs); |
|
if (lg(ne)==1) return NULL; |
|
nf_get_sign(checknf(bnf), (long*)&s, (long*)&t); |
|
BS.r = r = s+t-1; |
|
P = (GEN)tnf[1]; n = degpol(P); |
|
ro = (GEN)tnf[3]; |
|
BS.ALH = (GEN)tnf[4]; |
|
MatFU = (GEN)tnf[5]; |
|
A = (GEN)tnf[6]; |
|
csts = (GEN)tnf[7]; |
|
c1 = (GEN)csts[1]; c1 = gmul(absi(rhs), c1); |
|
c2 = (GEN)csts[2]; |
|
BS.hal = (GEN)csts[3]; |
|
x0 = (GEN)csts[4]; |
|
eps5 = (GEN)csts[5]; |
|
Prec = gtolong((GEN)csts[6]); |
|
BS.MatFU = MatFU; |
|
BS.bak = mulss(n, (n-1)*(n-2)); /* safe */ |
|
*pS = cgetg(1, t_VEC); |
|
|
SOL = cgetg(1,t_VEC); |
if (t) x0 = gmul(x0, absisqrtn(rhs, n, Prec)); |
upb = 0; ep = NULL; /* gcc -Wall */ |
tmp = divrr(c1,c2); |
|
c3 = mulrr(dbltor(1.39), tmp); |
|
c4 = mulsr(n-1, c3); |
|
x1 = gmax(x0, mpsqrtn(mulsr(2,tmp),n)); |
|
|
if (lg(thueres)==8) |
Vect = cgetg(r+1,t_COL); for (i=1; i<=r; i++) Vect[i]=un; |
|
Vect = gmul(gabs(A,DEFAULTPREC), Vect); |
|
c14 = mulrr(c4, Vecmax(Vect)); |
|
x2 = gmax(x1, mpsqrtn(mulsr(10,c14), n)); |
|
if (DEBUGLEVEL>1) { |
|
fprintferr("x1 -> %Z\n",x1); |
|
fprintferr("x2 -> %Z\n",x2); |
|
fprintferr("c14 = %Z\n",c14); |
|
} |
|
|
|
dP = derivpol(P); |
|
vecdP = cgetg(s+1, t_VEC); |
|
for (i=1; i<=s; i++) vecdP[i] = (long)poleval(dP, (GEN)ro[i]); |
|
|
|
zp1 = dbltor(0.01); |
|
x3 = gmax(x2, mpsqrtn(mulsr(2,divrr(c14,zp1)),n)); |
|
|
|
b = cgetg(r+1,t_COL); |
|
for (iroot=1; iroot<=s; iroot++) |
{ |
{ |
poly=(GEN)thueres[1]; deg=degpol(poly); gdeg=stoi(deg); |
GEN Delta, MatNE, Hmu, c5, c7; |
uftot=(GEN)thueres[2]; |
|
roo=gcopy((GEN)thueres[3]); |
|
ALH=(GEN)thueres[4]; |
|
MatFU=gcopy((GEN)thueres[5]); |
|
A=gcopy((GEN)thueres[6]); |
|
csts=(GEN)thueres[7]; |
|
|
|
if (!ne) ne = bnfisintnorm(uftot, rhs); |
Vect = cgetg(r+1,t_COL); for (i=1; i<=r; i++) Vect[i] = un; |
if (lg(ne)==1) { avma=av; return cgetg(1,t_VEC); } |
if (iroot <= r) Vect[iroot] = lstoi(1-n); |
|
Delta = gmul(A,Vect); |
|
|
c1=gmul(gabs(rhs,Prec), (GEN)csts[1]); c2=(GEN)csts[2]; |
c5 = Vecmax(gabs(Delta,Prec)); |
halpha=(GEN)csts[3]; |
c5 = myround(gprec_w(c5,DEFAULTPREC), 1); |
if (t) |
c7 = mulsr(r,c5); |
x0 = gmul((GEN)csts[4],gpui(gabs(rhs,Prec), ginv(gdeg), Prec)); |
c10 = divsr(n,c7); BS.c10 = c10; |
eps3=(GEN)csts[5]; |
c13 = divsr(n,c5); BS.c13 = c13; |
Prec=gtolong((GEN)csts[6]); |
if (DEBUGLEVEL>1) { |
b=cgetg(r+1,t_COL); |
fprintferr("* real root no %ld/%ld\n", iroot,s); |
tmp=divrr(c1,c2); |
fprintferr(" c10 = %Z\n",c10); |
x1=gmax(x0,gpui(mulsr(2,tmp),ginv(gdeg),ConstPrec)); |
fprintferr(" c13 = %Z\n",c13); |
if(DEBUGLEVEL >=2) fprintferr("x1 -> %Z\n",x1); |
} |
|
|
c3=mulrr(dbltor(1.39),tmp); |
prec = Prec; |
c4=mulsr(deg-1,c3); |
for (;;) |
|
{ |
|
if (( MatNE = Conj_LH(ne, &Hmu, ro, prec) )) break; |
|
prec = (prec<<1)-2; |
|
if (DEBUGLEVEL>1) err(warnprec,"thue",prec); |
|
ro = tnf_get_roots(P, prec, s, t); |
|
} |
|
BS.ro = ro; |
|
BS.iroot = iroot; |
|
|
for (numroot=1; numroot<=s; numroot++) |
for (ine=1; ine<lg(ne); ine++) |
{ |
{ |
ComputeConstants(); |
GEN Lambda, B0, c6, c8; |
Conj_Norm_Eq(ne,&Hmu); |
GEN NE = (GEN)MatNE[ine], Vect2 = cgetg(r+1,t_COL); |
for (curne=1; curne<lg(ne); curne++) |
int k, i1; |
|
|
|
if (DEBUGLEVEL>1) fprintferr(" - norm sol. no %ld/%ld\n",ine,lg(ne)-1); |
|
for (k=1; k<=r; k++) |
{ |
{ |
ComputeConstants2(poly,rhs); |
if (k == iroot) |
Baker(ALH,Hmu); |
tmp = gdiv(rhs, gmul((GEN)vecdP[k], (GEN)NE[k])); |
|
else |
i1=Vecmax(gabs(Delta,Prec),r); |
tmp = gdiv(gsub((GEN)ro[iroot],(GEN)ro[k]), (GEN)NE[k]); |
if (i1!=1) i2=1; else i2=2; |
Vect2[k] = llog(gabs(tmp,prec), prec); |
do |
} |
{ |
Lambda = gmul(A,Vect2); |
fs=0; |
|
Create_CF_Values(i1,i2,&errdelta); |
|
if (DEBUGLEVEL>=2) fprintferr("Entering CF\n"); |
|
Old_B0=gmul(B0,gdeux); |
|
|
|
/* Try to reduce B0 while |
c6 = addrr(dbltor(0.1), Vecmax(gabs(Lambda,DEFAULTPREC))); |
* 1) there was less than 10 reductions |
c6 = myround(c6, 1); |
* 2) the previous reduction improved the bound of more than 0.1. |
c8 = addrr(dbltor(1.23), mulsr(r,c6)); |
*/ |
c11= mulrr(mulsr(2,c3) , mpexp(divrr(mulsr(n,c8),c7))); |
while (Nb_It<10 && gcmp(Old_B0,gadd(B0,dbltor(0.1))) == 1 && fs<2) |
c15= mulrr(mulsr(2,c14), mpexp(divrr(mulsr(n,c6),c5))); |
{ |
|
cf=0; Old_B0=B0; Nb_It++; Kstart=stoi(10); |
|
while (!fs && CF_First_Pass(Kstart,errdelta) == 0 && cf < 8 ) |
|
{ |
|
cf++; |
|
Kstart=mulis(Kstart,10); |
|
} |
|
if ( CF_First_Pass(Kstart,errdelta) == -1 ) |
|
{ ne = gerepilecopy(av, ne); Prec+=10; |
|
if(DEBUGLEVEL>=2) fprintferr("Increasing precision\n"); |
|
thueres=thueinit(thueres,0,Prec); |
|
return(thue(thueres,rhs,ne)); } |
|
|
|
if (cf==8 || fs) /* Semirational or totally rational case */ |
if (DEBUGLEVEL>1) { |
{ |
fprintferr(" c6 = %Z\n",c6); |
if (!fs) |
fprintferr(" c8 = %Z\n",c8); |
{ ep=GuessQi(&q1,&q2,&q3); } |
fprintferr(" c11 = %Z\n",c11); |
bd=gmul(denom(bestappr(delta,gadd(B0,gabs(q2,Prec)))) |
fprintferr(" c15 = %Z\n",c15); |
,delta); |
} |
bd=gabs(gsub(bd,ground(bd)),Prec); |
BS.c11 = c11; |
if(gcmp(bd,ep)==1 && !gcmp0(q3)) |
BS.c15 = c15; |
{ |
BS.NE = NE; |
fs=1; |
BS.Hmu = (GEN)Hmu[ine]; |
B0=gdiv(glog(gdiv(gmul(q3,c15),gsub(bd,ep)), Prec),c13); |
|
if (DEBUGLEVEL>=2) |
|
fprintferr("Semirat. reduction ok. B0 -> %Z\n",B0); |
|
} |
|
else fs=2; |
|
} |
|
else fs=0; |
|
Nb_It=0; B0=gmin(Old_B0,B0); upb=gtolong(gceil(B0)); |
|
} |
|
i2++; if (i2==i1) i2++; |
|
} |
|
while (fs == 2 && i2 <= r); |
|
|
|
if (fs == 2) TotRat(); |
|
|
|
if (DEBUGLEVEL >=2) fprintferr("x2 -> %Z\n",x2); |
i1 = Vecmaxind(gabs(Delta,prec)); |
|
if (! (B0 = get_B0(i1, Delta, Lambda, eps5, prec, &BS)) ) goto PRECPB; |
|
|
/* For each possible value of b_i1, compute the b_i's |
/* For each possible value of b_i1, compute the b_i's |
* and 2 conjugates of x-\alpha y. Then check. |
* and 2 conjugates of z = x - alpha y. Then check. */ |
*/ |
upb = gtolong(gceil(B0)); |
zp1=dbltor(0.01); |
for (bi1=-upb; bi1<=upb; bi1++) |
x3=gmax(x2,gpui(mulsr(2,divrr(c14,zp1)),ginv(gdeg),ConstPrec)); |
{ |
|
GEN z1, z2; |
|
for (i=1; i<=r; i++) |
|
{ |
|
b[i] = ldiv(gsub(gmul((GEN)Delta[i], stoi(bi1)), |
|
gsub(gmul((GEN)Delta[i],(GEN)Lambda[i1]), |
|
gmul((GEN)Delta[i1],(GEN)Lambda[i]))), |
|
(GEN)Delta[i1]); |
|
if (gcmp(distoZ((GEN)b[i]), zp1) > 0) break; |
|
} |
|
if (i <= r) continue; |
|
|
for (bi1=-upb; bi1<=upb; bi1++) |
z1 = z2 = gun; |
{ |
for(i=1; i<=r; i++) |
flag=1; |
{ |
xmay1=gun; xmay2=gun; |
GEN c = ground((GEN)b[i]); |
for (i=1; i<=r; i++) |
z1 = gmul(z1, powgi(gcoeff(MatFU,1,i), c)); |
{ |
z2 = gmul(z2, powgi(gcoeff(MatFU,2,i), c)); |
b[i]=(long)gdiv(gsub(gmul((GEN)Delta[i],stoi(bi1)), |
} |
gsub(gmul((GEN)Delta[i],(GEN)Lambda[i1]), |
z1 = gmul(z1, (GEN)NE[1]); |
gmul((GEN)Delta[i1],(GEN)Lambda[i]))), |
z2 = gmul(z2, (GEN)NE[2]); |
(GEN)Delta[i1]); |
if (!CheckSol(pS, z1,z2,P,rhs,ro)) goto PRECPB; |
if (gcmp(distoZ((GEN)b[i]),zp1)==1) { flag=0; break; } |
|
} |
|
|
|
if(flag) |
|
{ |
|
for(i=1; i<=r; i++) |
|
{ |
|
b[i]=lround((GEN)b[i]); |
|
xmay1=gmul(xmay1,gpui(gcoeff(MatFU,1,i),(GEN)b[i],Prec)); |
|
xmay2=gmul(xmay2,gpui(gcoeff(MatFU,2,i),(GEN)b[i],Prec)); |
|
} |
|
xmay1=gmul(xmay1,gcoeff(MatNE,1,curne)); |
|
xmay2=gmul(xmay2,gcoeff(MatNE,2,curne)); |
|
Check_Solutions(xmay1,xmay2,poly,rhs); |
|
} |
|
} |
|
} |
} |
} |
} |
if (DEBUGLEVEL>=2) fprintferr("Checking for small solutions\n"); |
|
Check_Small((int)(gtodouble(x3)),poly,rhs); |
|
return gerepilecopy(av,SOL); |
|
} |
} |
|
*pro = ro; return x3; |
|
|
/* Case s=0. All the solutions are "small". */ |
PRECPB: |
Prec=DEFAULTPREC; poly=(GEN)thueres[1]; deg=degpol(poly); |
ne = gerepilecopy(av, ne); |
gdeg=stoi(deg); c0=(GEN)thueres[2]; |
prec += 5 * (DEFAULTPREC-2); |
roo=roots(poly,Prec); |
if (DEBUGLEVEL>1) err(warnprec,"thue",prec); |
Check_Small(gtolong(ground(gpui(gmul((GEN)poly[2],gdiv(gabs(rhs,Prec),c0)), |
tnf = inithue(P, bnf, 0, prec); |
ginv(gdeg),Prec) )), poly, rhs); |
return LargeSols(tnf, rhs, ne, pro, pS); |
return gerepilecopy(av,SOL); |
|
} |
} |
|
|
|
/* Given a tnf structure as returned by thueinit, a RHS and |
|
* optionally the solutions to the norm equation, returns the solutions to |
|
* the Thue equation F(x,y)=a |
|
*/ |
|
GEN |
|
thue(GEN tnf, GEN rhs, GEN ne) |
|
{ |
|
gpmem_t av = avma; |
|
GEN P, ro, x3, S; |
|
|
|
if (!checktnf(tnf)) err(talker,"not a tnf in thue"); |
|
if (typ(rhs) != t_INT) err(typeer,"thue"); |
|
|
|
P = (GEN)tnf[1]; |
|
if (lg(tnf) == 8) |
|
{ |
|
x3 = LargeSols(tnf, rhs, ne, &ro, &S); |
|
if (!x3) { avma = av; return cgetg(1,t_VEC); } |
|
} |
|
else |
|
{ /* Case s=0. All solutions are "small". */ |
|
GEN c0 = (GEN)tnf[2]; /* t_REAL */ |
|
S = cgetg(1,t_VEC); |
|
ro = roots(P, DEFAULTPREC); |
|
x3 = mpsqrtn(mulir(constant_term(P), divir(absi(rhs),c0)), degpol(P)); |
|
} |
|
return gerepilecopy(av, SmallSols(S, gtolong(x3), P, rhs, ro)); |
|
} |
|
|
static GEN *Relations; /* primes above a, expressed on generators of Cl(K) */ |
static GEN *Relations; /* primes above a, expressed on generators of Cl(K) */ |
static GEN *Partial; /* list of vvectors, Partial[i] = Relations[1..i] * u[1..i] */ |
static GEN *Partial; /* list of vvectors, Partial[i] = Relations[1..i] * u[1..i] */ |
static GEN *gen_ord; /* orders of generators of Cl(K) given in bnf */ |
static GEN *gen_ord; /* orders of generators of Cl(K) given in bnf */ |
Line 853 test_sol(long i) |
|
Line 811 test_sol(long i) |
|
|
|
if (Partial) |
if (Partial) |
{ |
{ |
long av=avma; |
gpmem_t av=avma; |
for (k=1; k<lg(Partial[1]); k++) |
for (k=1; k<lg(Partial[1]); k++) |
if ( signe(modii( (GEN)Partial[i][k], gen_ord[k] )) ) |
if ( signe(modii( (GEN)Partial[i][k], gen_ord[k] )) ) |
{ avma=av; return; } |
{ avma=av; return; } |
Line 872 test_sol(long i) |
|
Line 830 test_sol(long i) |
|
|
|
for (k=1; k<=i; k++) sol[k]=u[k]; |
for (k=1; k<=i; k++) sol[k]=u[k]; |
for ( ; k<=Nprimes; k++) sol[k]=0; |
for ( ; k<=Nprimes; k++) sol[k]=0; |
if (DEBUGLEVEL > 2) |
if (DEBUGLEVEL>2) |
{ |
{ |
fprintferr("sol = %Z\n",sol); |
fprintferr("sol = %Z\n",sol); |
if (Partial) fprintferr("Partial = %Z\n",Partial); |
if (Partial) fprintferr("Partial = %Z\n",Partial); |
Line 882 test_sol(long i) |
|
Line 840 test_sol(long i) |
|
static void |
static void |
fix_Partial(long i) |
fix_Partial(long i) |
{ |
{ |
long k, av = avma; |
long k; |
|
gpmem_t av = avma; |
for (k=1; k<lg(Partial[1]); k++) |
for (k=1; k<lg(Partial[1]); k++) |
addiiz( |
addiiz( |
(GEN) Partial[i-1][k], |
(GEN) Partial[i-1][k], |
Line 1062 bnfisintnorm(GEN bnf, GEN a) |
|
Line 1021 bnfisintnorm(GEN bnf, GEN a) |
|
{ |
{ |
GEN nf,pol,res,unit,x,id, *Primes; |
GEN nf,pol,res,unit,x,id, *Primes; |
long sa,i,j,norm_1; |
long sa,i,j,norm_1; |
ulong av = avma; |
gpmem_t av = avma; |
|
|
bnf = checkbnf(bnf); nf = (GEN)bnf[7]; pol = (GEN)nf[1]; |
bnf = checkbnf(bnf); nf = (GEN)bnf[7]; pol = (GEN)nf[1]; |
if (typ(a)!=t_INT) |
if (typ(a)!=t_INT) |
err(talker,"expected an integer in bnfisintnorm"); |
err(talker,"expected an integer in bnfisintnorm"); |
sa = signe(a); |
sa = signe(a); |
if (sa == 0) { res=cgetg(2,t_VEC); res[1]=zero; return res; } |
if (sa == 0) return _vec(gzero); |
if (gcmp1(a)) { res=cgetg(2,t_VEC); res[1]=un; return res; } |
if (gcmp1(a)) return _vec(gun); |
|
|
get_sol_abs(bnf, absi(a), &Primes); |
get_sol_abs(bnf, absi(a), &Primes); |
|
|
Line 1077 bnfisintnorm(GEN bnf, GEN a) |
|
Line 1036 bnfisintnorm(GEN bnf, GEN a) |
|
norm_1 = 0; /* gcc -Wall */ |
norm_1 = 0; /* gcc -Wall */ |
for (i=1; i<=sindex; i++) |
for (i=1; i<=sindex; i++) |
{ |
{ |
id = gun; x = normsol[i]; |
x = normsol[i]; |
for (j=1; j<=Nprimes; j++) /* compute prod Primes[i]^u[i] */ |
id = idealpow(nf,Primes[1], stoi(x[1])); |
if (x[j]) |
for (j=2; j<=Nprimes; j++) /* compute prod Primes[i]^u[i] */ |
{ |
id = idealmulpowprime(nf,id, Primes[j], stoi(x[j])); |
GEN id2 = Primes[j]; |
|
if (x[j] != 1) id2 = idealpow(nf,id2, stoi(x[j])); |
|
id = idealmul(nf,id,id2); |
|
} |
|
x = (GEN) isprincipalgenforce(bnf,id)[2]; |
x = (GEN) isprincipalgenforce(bnf,id)[2]; |
x = gmul((GEN)nf[7],x); /* x possible solution */ |
x = gmul((GEN)nf[7],x); /* x possible solution */ |
if (signe(gnorm(gmodulcp(x,(GEN)nf[1]))) != sa) |
if (signe(gnorm(gmodulcp(x,(GEN)nf[1]))) != sa) |
Line 1093 bnfisintnorm(GEN bnf, GEN a) |
|
Line 1048 bnfisintnorm(GEN bnf, GEN a) |
|
if (norm_1) x = gmul(unit,x); |
if (norm_1) x = gmul(unit,x); |
else |
else |
{ |
{ |
if (DEBUGLEVEL > 2) |
if (DEBUGLEVEL > 2) fprintferr("%Z eliminated because of sign\n",x); |
fprintferr("%Z eliminated because of sign\n",x); |
|
x = NULL; |
x = NULL; |
} |
} |
} |
} |