Bill Allombert on Sun, 30 Nov 2008 13:42:44 +0100

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

 Re: Difficulty with binary quadratic equation solver

```On Fri, Nov 28, 2008 at 09:36:47AM -0700, Kurt Foster wrote:
> Bill provides the example
>
> [a,b,c,d,e,f] = [1,0,-5,1,1,2]

(I'd like to note that http://www.alpertron.com.ar/QUAD.HTM
does not handle that equation either, so we cannot compare the results).

> for which the field discriminant N = 5, the index m is 2, U =
> (B.fu)^2, the smallest value of k is 3, and plowing through the
> calculations by hand shows that there are two "inequivalent" solutions
>
> [x,y] = [1,1] and [-5, -2]
>
> For a given positive fundamental discriminant N and index m > 1,
> determining the exact smallest value of k for which U^k is in O has
> proven to be rather exasperating, although certain multiples of it are
> easy to find.  I've dug into it a fair bit, but I'm sure this is a
> known calculation.  Can someone provide a reference?

Well, we can use the ray class number formula that give the ray class number
in term of the class number:

h_f=h*Phi(f)/phi(f)/k

So the following should compute k:

K=bnfinit(x^2-b*x+a*c);
B=bnrinit(K,K.index);
k=K.no*B.bid.no/(eulerphi(K.index)*B.no);

I join a new version of queq.gp that use this formula.

Cheers,
Bill.
```
```qeeval(Q,v)=local(a=Q[1],b=Q[2],c=Q[3],d=Q[4],e=Q[5],f=Q[6],x=v[1],y=v[2]);a*x^2+b*x*y+c*y^2+d*x+e*y+f

qered(Q)=
{
local(a=Q[1],b=Q[2],c=Q[3],d=Q[4],e=Q[5],f=Q[6],F,u,v);
F=[2*a,b;b,2*c]^-1*[-d,-e]~;u=F[1];v=F[2];
D=denominator(content(F));
[D^2*(u^2*a+v*u*b+v^2*c+u*d+v*e+f),D,D*u,D*v]
}

qesolve(Q)=
{
local(a,b,c,f,V,K,B,d,L,w,U,Ud,S=List(),idx);
local(R,D,M,u,v);
a=Q[1];b=Q[2];c=Q[3];
if (issquare(b^2-4*a*c),error("reducible equation"));
R=qered(Q); f=R[1]; D=R[2];
K=bnfinit(x^2-b*x+a*c);
L=bnfisintnorm(K,-f*a);
if (#K.fu,
idx=K.index; B=bnrinit(K,idx);
d=K.no*B.bid.no/(eulerphi(idx)*B.no);
U=K.fu[1]; if (norm(U)==-1,U=U^2);
,d=1; U=1);
for(i=1,#L,
for(j=0,K.tu[1]-1,
for(k=0,d-1,
w=U^k*L[i]*K.tu[2]^j;
w=lift(w); u=polcoeff(w,0); v=polcoeff(w,1);
if (denominator(u)!=1 || denominator(v)!=1,next);
if (u%a!=0,next,u\=a);
if ((R[3]-u)%D!=0,next,u=(R[3]-u)\D);
if ((R[4]-v)%D!=0,next,v=(R[4]-v)\D);
listput(S,[u,v]);
)));Vec(S)
}
```