F<t,a27> := FunctionField(Rationals(),2);

// Let Q(a27) be the function field of X_{27A^{0}-27a}.
// Let Q(t) be the function field of X_{0}(2). We use the model j = t^3/(t+16)
// from the work of Rouse, Zureick-Brown.

b9 := a27^3;
b3 := b9*(b9^2+9*b9+27);
j1 := (b3+3)^3*(b3+27)/b3;

j2 := t^3/(t+16);

pol := Numerator(j1-j2);
C := Curve(AffineSpace(Rationals(),2),pol);
C2 := ProjectiveClosure(C);
FF := FunctionField(C2);
jmap := (FF.1^3)/(FF.1+16);

// C2 is the fiber product of X_{27A^{0}-2a} and X_{0}(2) over X_{0}(1).
// It is singular.

printf "Genus of C2 is %o.\n",Genus(C2);
A, B, mp := IsHyperelliptic(C2);
// C2 is hyperelliptic
B2, mp2 := SimplifiedModel(B);
printf "Simplified model is %o.\n",B2;
printf "Rank bounds are %o.\n",RankBounds(Jacobian(B2));
// The Jacobian has rank 0.
pts := Chabauty0(Jacobian(B2));
printf "Rational points are %o.\n",pts;
// The set pts is a provably complete list of the rational points on B2.

mp2inv := Inverse(mp2);
Bpts := [ mp2inv(p) : p in pts ];
// Bpts is a list of the points on B.


// We find all the points on C2 by taking preimages of rational points on
// B.
Cpts := [];
for x in Bpts do
  X := x@@mp;
  ratpointsonX := RationalPoints(X);
  for p in ratpointsonX do
    if not C2!Eltseq(p) in Cpts then
      Append(~Cpts,C2!Eltseq(p));
    end if;
  end for;
end for;

printf "There are %o rational points on C2.\n",#Cpts;
for i in [1..#Cpts] do
  printf "Point %o is %o.\n",i,Cpts[i];
  if IsNonsingular(Cpts[i]) then
    printf "It's non-singular. It's image on the j-line is %o.\n",Evaluate(jmap,Cpts[i]);
  else
    plac := Places(Cpts[i]);
    printf "It's singular. There are %o places above it.\n",#plac;
    for j in [1..#plac] do
      if Degree(plac[j]) eq 1 then
        printf "Place %o has degree 1. It's image on the j-line is %o.\n",j,Evaluate(jmap,plac[j]);
      end if;
    end for;
  end if;
end for;