//////////////////////////////////////////////////////////////////
// This code determines the rational points on X686
//////////////////////////////////////////////////////////////////


  A3<x,y,w> := AffineSpace(Rationals(),3);
  eqnsX686 := 
    [y^2 - (x^3 + x^2 - 3*x + 1),
     w^2 - 2*(y^2  -y*(-x + 1))*(x^2 - 2*x - 1)];
  X686 := Curve(A3, eqnsX686);



// The double cover of  X686 maps to the following hyperelliptic curve for d in {1,-1,2, -2} do
for d in {1,-1,2, -2} do
  A3<x,y,w> := AffineSpace(Rationals(),3);
  C1 := Curve(A3, 
  [y^2 - (x^3 + x^2 - 3*x + 1),
  d*w^2 - (x^2 - 2*x - 1)]
  );
  
  aut := AutomorphismGroup(C1);
  for i in [1..3] do 
    H, prC1H  :=   CurveQuotient(AutomorphismGroup(C1,[[aut.1, aut.2, aut.1*aut.2][i]]));
    if Genus(H) eq 2 then 
        r := RankBounds(Jacobian(H)); // returns 0 or 1
       // Do Chabauty on H
       if r eq 1 then
          J := Jacobian(H);
            //RankBound(J); // returns 1
          heightconst := HeightConstant(J : Effort:=2, Factor);
            pts := Points(J : Bound := 1000);
            min := Minimum({Height(P) : P in pts | Order(P) eq 0}); 
            bool := exists(PJ){pt : pt in pts | Height(pt) eq min};
            LogarithmicBound := Height(PJ) + heightconst;  // Bound on h(Q)
            AbsoluteBound := Ceiling(Exp(LogarithmicBound));
            PtsUpToAbsBound := Points(J : Bound:=AbsoluteBound );
            bas, ht := ReducedBasis([ pt : pt in PtsUpToAbsBound ]);
            ptsH :=    Chabauty(bas[1]);                 
       else
            ptsH :=    Chabauty0(Jacobian(H));
       end if;

            ptsC1 :=  {@qt  : qt in RationalPoints(pt @@ prC1H), pt in ptsH@};
            ptsC1; // all satisfy xy = 0 or x = -1.
    end if;      
  end for;
  
end for;
    
    
// Coordinates of x,y,w on C1 are the same as on the point coming from X686    
// all satisfy xy = 0 or x = -1.
  A2<x,y> := AffineSpace(Rationals(),2);
    [X686![0,a[1],a[2]] :  a in RationalPoints(Scheme(A2,[Evaluate(eqnsX686[i],[0,x,y])  : i in [1,2]]))];
    [X686![a[1],0,a[2]] :  a in RationalPoints(Scheme(A2,[Evaluate(eqnsX686[i],[x,0,y])  : i in [1,2]]))];
    [X686![-1,a[1],a[2]] :  a in RationalPoints(Scheme(A2,[Evaluate(eqnsX686[i],[-1,x,y])  : i in [1,2]]))];
    
// gives exactly the 3 known points.    
    