/* WARNING: I haven't given an algorithm to compute the geometric component
   group at non-semistable primes. */


no_curve := [0,0];

visdat:=
[<3, [681,2], [681,3]>,
 <5, [1058,4], [1058,3]>,
 <5, [1246,2], [1246,3]>,
 <5, [1664,11], [1664,14]>,
 <3, [1913,2], [1913,1]>,
 <3, [2006,5], [2006,4]>,
 <3, [2366,4], [2366,5]>, 
 <5, [2366,6], [2366,5]>,  
 <3, [2429,2], [2429,4]>,  
 <3, [2534,5], [2534,7]>,  
 <3, [2534,6], [2534,7]>,  
 <3, [2541,4], [2541,3]>,  
 <5, [2574,4], [2574,7]>,
 <3, [2601,8], [2601,12]>,
 <3, [2674,2], [2674,1]>,
 <3, [2710,3], [2710,2]>,
 <3, [2718,4], [2710,2]>, 
 <3, [2768,3], [2768,2]>,
 <5, [2834,4], [2834,3]>,
 <3, [2849,1], no_curve>, 
 <5, [2900,4], [2900,3]>,
 <3, [2932,1], no_curve>,
 <3, [2955,2], [2955, 3]>,
 <3, [3054,1], [3054,3]>,
 <5, [3185,3], [3185,2]>,
 <3, [3306,2], [1102,1]>,
 <7, [3364,3], no_curve>,
 <5, [3384,1], [3384,3]>,
 <3, [3536,8], [3536,7]>,
 <3, [3555,5], [3555,4]>,
 <3, [3712,10], [3712,9]>,
 <3, [3879,5], [3879,4]>,
 <3, [3933,1], [3933,2]>,
 <5, [3952,3], [3952,5]>,
 <3, [3954,3], [3954, 4]>,
 <5, [4092,1], [4092,2]>,
 <3, [4229,1], no_curve>,
 <3, [4343,3], no_curve>,
 <5, [4592,4], [4592,7]>,
 <3, [4592,6], [4592,3]>,
 <3, [4592,6], [4592,7]>,
 <3, [4606,2], [4606,3]>,
 <3, [4675,10], [4675,9]>,
 <3, [4914,14], no_curve>,
 <3, [4963,3], [4963,4]>,
 <3, [5046,8], [5046,10]>,
 <3, [5054,3], no_curve>,
 <3, [5073,4], no_curve>,
 <5, [5082,3], [5082, 4]>,
 <3, [5136,2], [1712,4]>,
 <3, [5389,1], no_curve>
/* <3, [5499,5], [5499,6]> */  // not in database.
];


function E(label)  // label = [N,iso].
   if label eq no_curve then
      error "No such curve.";
   end if;
   return EllipticCurve(CremonaDatabase(),label[1],label[2],1);
end function;

function verify_ranks()
   for dat in visdat do
      if dat[3] ne no_curve then
         e := E(dat[2]);
         f := E(dat[3]);
/*         if Rank(e) ne 0 then
            error "Lie in ", dat, ": rank is not zero.";
         end if;
*/
         if Rank(f) ne 2 then
            error "Lie in ", dat, ": rank is not two.";
         end if;
         dat, "passed.";
      end if;
   end for; 
end function;


function cpbar(e,p)
   assert Type(e) eq CrvEll;
   assert Type(p) eq RngIntElt;
   assert IsPrime(p);
 
   N := Conductor(e);
   if Valuation(N,p) eq 0 then
      return 1;
   end if;
   if Valuation(N,p) eq 1 then
      return -Valuation(Discriminant(e),p);
   end if;
   return TamagawaNumber(e,p);   // this is wrong, but probably not by much!
end function;


function Product_cpbar(e)
   assert Type(e) eq CrvEll;
   return &*[cpbar(e,p[1]) : p in Factorization(Conductor(e))];
end function;


function Product_cp(e)
   assert Type(e) eq CrvEll;
   return &*[TamagawaNumber(e,p[1]) : p in Factorization(Conductor(e))];
end function;


procedure verify_visible()
   for dat in visdat do
      if dat[3] ne no_curve then
         e := E(dat[2]);
         f := E(dat[3]);
         p := dat[1];
         N := dat[2][1];
         bad := false;
         if N mod p eq 0 then
            print N,":  p =", p, " divides the level.";
            bad := true;
         end if;
         if Product_cpbar(e) mod p eq 0 then
            print N,":  p =", p, " divides the prod cpbar(E).";
            bad := true;
         end if;
         if Product_cp(f) mod p eq 0 then
            print N,":  p =", p, " divides the prod cp(F).";
            bad := true;
         end if;
         if not (N eq SquareFree(N)) then
            print N,":  additive reduction.";
            bad := true;
         end if;
         if not bad then
            print N,":  OK";
         end if;
      end if;
   end for; 
end procedure;


