//easy_primes
Q<x> := PolynomialRing(Rationals());
label := <"121a1", "121b1", "121c1", "361a1",
		"1225h1", "1225h2", "1849a1", "4489a1">;
prime := <11, 11, 11, 19, 37, 37, 43, 67>;
S := [];
for i := 1 to 8 do
	E := EllipticCurve(label[i]);
	p := prime[i];
	T := Factorization(DivisionPolynomial(E,p) : DegreeLimit := p-1);
	S := Append(S, T);
	printf "Case E = %o:\n", label[i];
	for p in S[i] do
		Degree(p[1]);
	end for;
end for;
function f(K, S, label, l, r)
	for i := l to r do
		printf "Case E = %o:\n", label[i];
		for p in S[i] do
			HasRoot(p[1], K);
		end for;
	end for;
	return "done!";
end function;
//11
K := Subfields(CyclotomicField(25),5)[1][1];
f(K, S, label, 1, 3);
//19
K := Subfields(CyclotomicField(81),9)[1][1];
f(K, S, label, 4,4);
//37, deg = 6
F := Subfields(CyclotomicField(9),3)[1][1];
K := Compositum(F, NumberField(x^2 - 2));
f(K, S, label, 5, 5);
//37, deg = 18
F := Subfields(CyclotomicField(81),9)[1][1];
K := Compositum(F, NumberField(x^2 - 2));
f(K, S, label, 6, 6);
//43
F := Subfields(CyclotomicField(49),7)[1][1];
K := Subfields(CyclotomicField(9),3)[1][1];
K := Compositum(K, F);
f(K, S, label, 7, 7);
//67
F := Subfields(CyclotomicField(121),11)[1][1];
K := Subfields(CyclotomicField(9),3)[1][1];
K := Compositum(K, F);
f(K, S, label, 8, 8);

//p = 163
E := EllipticCurve("26569a1");
S := Factorization(DivisionPolynomial(E,163) : DegreeLimit := 162);
for p in S do
	Degree(p[1]);
	end for;
K := Subfields(CyclotomicField(9),3)[1][1];
p := ChangeRing(S[1][1],K);
IsIrreducible(p);

//10-torsion
A<u, v> := AffineSpace(Rationals(),2);
C := ProjectiveClosure(Curve(A,25*(v^2+10*v+5)^3-u^2*v^5-1728*v^5));
IsSingular(C);
DescentInformation(EllipticCurve(C,SingularPoints(C)[1]));
Points(C : Bound := 10);

//13-torsion
Q<x> := PolynomialRing(Rationals());
K := Subfields(CyclotomicField(9),3)[1][1];
K := Compositum(K, NumberField(x^4 - 4*x^2 + 2));
f := function(K, E)
	p := DivisionPolynomial(E, 13);
	fact := Factorization(p : DegreeLimit := 12);
	E := ChangeRing(E, K);
	for p in fact do
		r := Roots(p[1], K);
		for i in r do
			if #Points(E, i[1]) gt 0 then
				return true;
				end if;
			end for;
		end for;
	return false;
end function;
S := {};
for d in data do
	if f(K, EllipticCurve(d)) then
		S := S join {CremonaReference(EllipticCurve(d))};
	end if;
end for;
S;

//4+8-torsion
Q<x> := PolynomialRing(Rationals());
K := NumberField(x^8 - 8*x^6 + 20*x^4 - 16*x^2 + 2);
F := CyclotomicField(16);
L := CyclotomicField(32);
E := EllipticCurve("32a2");
TorsionSubgroup(ChangeRing(E,F));
TorsionSubgroup(ChangeRing(E,L));
RankBound(ChangeRing(E,K));
RankBound(ChangeRing(QuadraticTwist(E,-1),K));

//CM-case
Q<x> := PolynomialRing(Rationals());
J := [0, 54000, -12288000, 1728, 287496, -3375, 16581375, 8000, -32768,
	-884736, -884736000, -147197952000, -262537412640768000];
K := CyclotomicField(64);
for j in J do
	E := EllipticCurveWithjInvariant(j);
	p := DivisionPolynomial(E, 16) div DivisionPolynomial(E, 8);
	fact := Factorization(p : DegreeLimit := 32);
	printf "j = %o, #fact = %o\n", j, #fact;
	for i := 1 to #fact do
		printf "%o : %o\n", i, HasRoot(fact[i][1], K);
	end for;
end for;

//2+16-torsion
function aux_Id(m,s)
	if [Integers(2^s)!t : t in Eltseq(m)] eq [1,0,0,1]
		then return true;
		else return false;
	end if;
end function;
load "2primary_Ss.txt";
function Cuenta(G,N,s,d)
	M := Characteristic(BaseRing(G));
	GG := sub<GL(2,Integers(2^N)) | Generators(G),
		[[1,M,0,1],[1,0,M,1],[1+M,0,0,1],[1,0,0,1+M]]>;
	Odds := [[i,j] : i in [0..2^N-1], j in [0..2^N-1] | IsOdd(i) or IsOdd(j)]; 
	L := {};
	V := RSpace(GG);
	S := [];
	for v in Odds do
		Hv := Stabiliser(GG, V!v);
		Hvmod2s := [m : m in Hv | aux_Id(m,s)];
		HHv := sub<GL(2,Integers(2^N)) | Hvmod2s>;
		ord_s := #Hvmod2s;
		hv := Integers()!(Order(GG)/ord_s);
		if IsDivisibleBy(d,hv) then
			Kv := quo<GG | Core(GG,HHv)>;
			Gal := GroupName(Kv);
			if IsAbelian(Kv) eq true and Order(Kv) ge 16 then
				L := L join {Gal};
			end if;
		end if;
		end for;
	return L;
end function;
S := {};
for rzb in RZB do
	if rzb[3][2][4] le 16 then
		G := sub<GL(2,Integers(rzb[2])) | {m : m in rzb[4]}>;
	S := S join Cuenta(G,4,1,16);
	end if;
end for;
S;

//2-torsion
Q<x> := PolynomialRing(Rationals());
K := CyclotomicField(32);
f := function(E, K)
	r := Factorization(DivisionPolynomial(E,17) : DegreeLimit := 16);
	for p in r do if HasRoot(p[1], K) then return true; end if; end for;
	return false;
end function;
f(EllipticCurve("14450p1"), K);
f(EllipticCurve("14450p2"), K);