# Schoof Algorithm : working on an example with SageMath

I am a noob amateur interested in Elliptic Cryptography and I am trying to work on Schoof Algorithm on a small example with the help of Sagemath

the algorithm description i found in a pdf called "Schoof's algorithm 18.783 Elliptic Curves, Lecture #9 Spring 2015"

and I am trying to work out example 2.2.10 from "Pairings for Beginners (Craig Costello) on SageMath

I want to verify that the trace of the Frobenius in the quotient ring with the division polynomials (given directly by Sagemath) is equal to whats written in the lecture (we should find t3 = 0, t5 = 1)

There must be an error in the implementation due to my misunderstanding of the formulas but I cannot see where.

#Example 2.2.10
q=13
A = 2
B = 1
E = EllipticCurve(GF(q),[A,B])
f = x^3+2*x+1

def addition_in_End_El(a1,b1,a2,b2):
if (a1 != a2):
r = (b1-b2) *inverse_mod(a1-a2,h)
else:
r = (3*a1*a1 +A) * inverse_mod(2*b1*f,h)

a3 = (r*r*f - a1 - a2) % h
b3 = (r*(a1-a3)-b1) % h
return a3,b3

def scalar_multiplication_in_End_El(k,a,b):
if k==0:
return (0,0)
elif k==1:
return (a,b)
else:
res = scalar_multiplication_in_End_El(k-1,a,b)
a2 = res[0]
b2 = res[1]
return addition_in_End_El(a,b,a2,b2)


I then look at l=3

l=3
ql = q % l
h=E.division_polynomial(l)
ql,h


(1, 3*x^4 + 12*x^2 + 12*x + 9)

R.<x> = PolynomialRing(GF(q),sparse=True)
S = R.quotient(x^3 + 2*x+1, 'a'); a = S.gen()
RR.<y> = PolynomialRing(R)

f = x^3+2*x+1
Q  = x^q
Q_sq = x^(q^2)

pil_x = Q % h
pil_y=(f^((q-1)/2))% h
pil_x,y*pil_y


(x^3 + 11*x^2 + 8*x + 12, (12*x^3 + 2*x^2 + 5*x + 12)*y)

pil_sq_x = Q_sq % h
pil_sq_y=(f^(((q^2)-1)/2))% h
pil_sq_x,y*pil_sq_y


(x, 12*y)

addition_in_End_El(pil_sq_x,pil_sq_y,x,1)
# i would expect at this stage (0,0) because t3 = 0


(x, 1)

for l=5,

l=5
ql = q % l
h=E.division_polynomial(l)
ql,h


(3, 5*x^12 + 7*x^10 + 3*x^9 + 9*x^8 + 12*x^7 + 12*x^6 + 11*x^5 + 10*x^4 + 9*x^3 + x^2 + 6*x + 7)

pil_x = Q % h
pil_y=(f^((q-1)/2)) % h
pil_x,y*pil_y


(9*x^11 + 2*x^10 + 6*x^9 + 8*x^8 + 8*x^7 + 3*x^6 + 11*x^5 + 6*x^4 + 5*x^3 + 4*x^2 + 9*x, (7*x^11 + 12*x^10 + 9*x^9 + 8*x^8 + x^7 + 7*x^6 + x^5 + x^4 + 11*x^2 + x + 12)*y)

pil_sq_x = Q_sq % h
pil_sq_y=(f^(((q^2)-1)/2)) % h
pil_sq_x,y*pil_sq_y


(5*x^11 + 12*x^10 + 6*x^9 + 11*x^8 + 2*x^7 + 5*x^6 + 7*x^5 + 2*x^4 + 2*x^3 + 3*x^2 + 11*x + 5, (6*x^11 + 2*x^10 + 6*x^9 + 5*x^8 + 3*x^7 + 5*x^6 + 10*x^5 + 2*x^4 + 9*x^3 + 9*x^2 + 3*x + 9)*y)

res1 = scalar_multiplication_in_End_El(3,x,1)
res1


(6*x^11 + 6*x^10 + 4*x^9 + 7*x^8 + 10*x^7 + 5*x^6 + 3*x^5 + 4*x^4 + 8*x^3 + 2*x^2 + 12*x + 8, 3*x^10 + 6*x^9 + 8*x^8 + 8*x^7 + 11*x^5 + 7*x^4 + 5*x^3 + 7*x + 10)

addition_in_End_El(pil_sq_x,pil_sq_y,res1[0],res1[1])
# i am expecting the same than pil_x,pil_y because t5 = 1


(9*x^11 + 2*x^10 + 6*x^9 + 8*x^8 + 8*x^7 + 3*x^6 + 11*x^5 + 6*x^4 + 5*x^3 + 4*x^2 + 9*x, 7*x^11 + 12*x^10 + 9*x^9 + 8*x^8 + x^7 + 7*x^6 + x^5 + x^4 + 11*x^2 + x + 12)

Example from the Lecture

From the algorithm description

[

Thank you for your help ! My goal is to understand the theorems and other statements by working out those small examples

By the way, are there some blocks of the code that already exist in Sage and can be replaced by existing functions ?

• i've been spending time to debug, and i noticed that my implementation doesn't even calculate $[1]_l=(a1=x,b1=1)$ and its multiples well. For l=3, I get $[2]_l = (x,12)$ but then $[3]_l$=(x,12)+(x,1) = (x,1) and $[3]_l$=(x,1)+(x,12)=(x,12) which is obviously wrong while I'm expecting both to give (0,0). However i feel i have respected the textbook Jun 12, 2020 at 6:54

## 1 Answer

I tried to understand the code in the OP, but it was not easy to do so. (It has too often pep8 conventions ignored, this makes the not highlighted code even less readable. A readable code should also have sound names for the variables, a structural approach - best respecting the structures already provided by sage - and some comments for what it is going on, else it is a once-only payment.)

Here is my try to reconstitute the Example 2.2.10 in [Costello] = Pairings for Beginners, Craig Costello .

We define the curve, and compute the order, just to know which is the number we try to exhibit using Schoof's algorithm.

sage: E = EllipticCurve(GF(13), [2, 1])
sage: E
Elliptic Curve defined by y^2 = x^3 + 2*x + 1
over Finite Field of size 13
sage: E.order()
8
sage: E.frobenius().minpoly()
x^2 - 6*x + 13


So indeed, we have so far the relation / the compatibility $$8=|E(\Bbb F_{13})|=13+1-t$$, where $$t=6$$ is the trace of the Frobenius morphism $$\Phi$$, extracted (Vieta) from the minimal polynomial above. So: $$\Phi^2 -6\Phi +13=0\ .$$ Explicitly, testing this on some point $$(x,y)\in E(\Bbb F_p)$$, $$p=q=13\ ,$$ it reads: $$(x^{p^2}, y^{p^2})-6(x^p,y^p)+p(x,y)=O\ .$$ We would like to detect this trace $$t=6$$ "from its pieces, first piece modulo $$l=3$$, and then modulo $$l=5$$".

For this, we test the above formula for points of $$l$$-torsion, for some "small" values of $$l$$, since our $$p$$ is... the huge prime $$13$$.

We consider $$l=3$$. For this prime, we compute the division polynomial $$\psi_3(x)$$...

sage: E.division_polynomial(3)
3*x^4 + 12*x^2 + 12*x + 9


and decide to switch to the field $$F$$, an extension of $$\Bbb F_p=\Bbb F_{13}$$, so that $$l=3$$-torsion points are available. So the above polynomial $$g=3x^4 + 12x^2 + 12x + 9$$ should have "roots" in $$F$$, well, at least one, and the equation $$y^2=x^3+2x+1$$ should also have solution in $$F$$ (for at least one root "$$x$$"). We factor $$g$$, and use one of the factors, $$g_1$$ say. For this reasons, we work in the field F constructed below... (Instead of using in the ring $$R_3$$, as specified in [Costello].)

sage: R.<x> = PolynomialRing(GF(13))
sage: g = E.division_polynomial(3)(x)
sage: g
3*x^4 + 12*x^2 + 12*x + 9
sage: g.factor()
(3) * (x^2 + 5*x + 5) * (x^2 + 8*x + 11)
sage: g1, mul1 = g.factor()[0]
sage: g1    # the first factor
x^2 + 5*x + 5
sage: K.<a> = GF(13^g1.degree(), modulus=g1)
sage: S.<y> = PolynomialRing(K)
sage: F.<b> = K.extension(y^2 - (a^3 + 2*a + 1))

sage: EF = E.base_extend(F)

sage: P = EF.point([a, b])
sage: FrobP = EF.point([a^13, b^13])
sage: FrobFrobP = EF.point([a^(13^2), b^(13^2)])

sage: # now we search for a good k (being t modulo 3)...
sage: for k in [0..2]:
....:    if FrobFrobP - k*FrobP + P == EF(0):
....:        print(f"t = {k} modulo 3")
t = 0 modulo 3


(The search is not optimal, i should have computed FrobFrobP + P separately, then break when a good k was found...)

The same in less lines for $$l=5$$. We compute $$\psi_5$$, factor it, use one of the factors $$g_1$$, take a root $$a$$ of it in a proper extension, and pass to an other extension $$F$$, so that there is some $$b\in F$$ with $$b^2 = a^3+2a+1$$. This makes $$(a,b)$$ be a point in $$E(F)$$, which is a $$5$$-torsion point. We test the relation $$\Phi^2 -k\Phi +p=0$$ on $$P$$, i.e. $$\Phi^2(P) -k\Phi(P) +3P=O\ ,$$ with $$3=p$$ modulo $$l$$ instead of $$p$$ by using the fact that $$P$$ has order $$l$$, so for instance $$13P=(5+5+3)P=3P$$, and let $$k$$ take in the search all values from $$0$$ (inclusively) to $$l=5$$ (exclusively).

The code would be:

p = 13
E = EllipticCurve(GF(p), [2, 1])

el = 5    # i will never write l for a variable
R.<x> = PolynomialRing(GF(p))
g = E.division_polynomial(el)
g1, mul1 = g.factor()[0]
K.<a> = GF(p^g1.degree(), modulus=g1)
S.<y> = PolynomialRing(K)
F.<b> = K.extension( E.defining_polynomial()(a, y, 1) )

EF = E.base_extend(F)
P = EF.point([a, b])
FrobP = EF.point([a^p, b^p])
FrobFrobP = EF.point([a^(p*p), b^(p*p)])

Q = FrobFrobP + (p % el)*P
t_mod_el = None
for k in [0..(el-1)]:
if k*FrobP == Q:
t_mod_el = k
break

print(f"t modulo {el} = {t_mod_el}")


This gives:

t modulo 5 = 1


Well, it turns out that for $$l=5$$ we have above $$a=10\in \Bbb F_{13}=\Bbb F_p$$, and that $$b$$ is the square root of $$a^3+2a+1$$, which lives in an extension of degree two of $$\Bbb F_p$$. The code above works also in such "simplified" cases.

Observing this, shows that it is always a good idea to sort the polynomials in the factorization by their degree. Factorization may be a time consuming issue. But we can improve the choice of the $$l$$-primes so that "small factors" are involved / insured.

So we are searching for a number $$t$$ between $$0$$ and $$3\cdot 5$$, which is:

• $$t_3=0$$, when taken modulo $$l=3$$, and
• $$t_5=1$$, when taken modulo $$l=5$$.

The number $$t=6$$ is the match.

From the code snippets above, i hope it is easy to exhibit a piece of code of slightly increased length starting from the last block. Above, i tested the relation $$\Phi^2-t\Phi+p=0$$ on exactly one point $$P$$ or order $$l$$, this is ok so from the mathematical point of view, it is easily supported by sage, and didactically i would also prefer to work over fields.

So i did not work in the endomorphism ring $$\operatorname{End}(E)$$, and i did not work over the ring $$R_l$$ from [Costello]. If these are main issues, then longer pieces of code are needed. Or other computing strategies. The following does not work:

p = 13
E = EllipticCurve(GF(p), [2, 1])
el = 5
g = E.division_polynomial(el)    # a "polynomial in x"
f = E.defining_polynomial()      # a "polynomial in x, y, z"

R.<x,y> = PolynomialRing(GF(p))
J = R.ideal( [g(x), f(x, y, 1)] )
RmodJ = R.quotient(J)
EJ = E.base_extend(RmodJ)

P         = EJ.point( [x      , y       ] )
FrobP     = EJ.point( [x^p    , y^p     ] )
FrobFrobP = EJ.point( [x^(p*p), y^(p*p) ] )

# !!! The following line does not work !!!

Q = FrobFrobP + (p % el)*P

t_mod_el = None
for k in [0..(el-1)]:
if k*FrobP == Q:
t_mod_el = k
break

print(f"t modulo {el} = {t_mod_el}")


Instead, addition and multiplication should be done manually, the components of $$\Phi^2(P)+pP$$ and $$k\Phi(P)=\Phi(kP)$$ have then to be compared explicitly modulo the ideal J above.

• thank you very much for your time and detailed response. I will definitely have a look at it. Given my beginner level, it will take some time to absorb. I just read quickly: are you saying there is a mistake in Costello pdf where he confuses a GF(p^n) and Z/p^n Z somewhere ? Jun 17, 2020 at 1:06
• If there are questions, please do not hesitate to ask. There is no mistake in the paper of Costello, but for the reason of a quick implementation i had to change slightly the way to proceed (for the two given primes). For instance, for $l=3$ Costello would work with points of the curve $E$ over the ring (not field)$$R_3=\Bbb F_{13}[x,y]\text{ modulo }(3x^4 + 12x^2 + 12x + 9,\ y^2-(x^3+2x+1))\ .$$This ring is a product of fields, and (**if we know** the factorization of the division polynomia) these fields can be constructed explicitly. In fact, we need only one of them to get $t_l$. Jun 17, 2020 at 9:14
• thank you for your explanations, but i was precisely trying to do the calculations in the ring $R_3$ using the formulas in End(E) given by the other pdf. And i'm going step by step (again) and i just can't make it work, and i'm losing my hair with it LOL Jun 17, 2020 at 14:01
• ok i finally got it, there is a little omission in the 2nd PDF, when giving the addition formula. using their notations, we need to distinguish a few more cases, namely (x1=x2, y1=y2) and (x1=x2, y1=-y2), and i'm getting the right result. Jun 17, 2020 at 15:52
• @dan_fulea Hello, can you kindly tell what changes is needed in your code to make it work over fields of order $p^n$ where $n>1$? It's working perfectly for $n=1$, but I am having difficulties for $n>1$ cases. Apr 27, 2023 at 13:54