Isogenies 2: So... what are isogenies again? (in Sage this time)

March 23, 2024

Heya! Welcome to my second isogeny post. In this blog post, I will demonstrate claims I made in the previous post with code in Sage. Personally, I find experimenting in Sage incredibly useful, especially since there are a few isogenists actively working on improving its functionalities. The main target audience is for beginner isogenists, but on the other hand, if you are an experienced isogenist with less Sage experience, maybe this blog will help too! If you haven’t used Sage before, you can try it out online. Let’s get started!

Introduction

Let’s begin by reviewing how we define basic objects. We construct finite fields with GF:

sage: GF(13)
Finite Field of size 13
sage: GF(13^5)
Finite Field in z5 of size 13^5
sage: GF((13, 5))
Finite Field in z5 of size 13^5
sage: GF(13^5, name="t")
Finite Field in t of size 13^5
sage: K.<t> = GF(13^5); K
Finite Field in t of size 13^5
sage: t.minpoly()
x^5 + 4*x + 11
sage: GF(13) is GF(13)
True

We can also specify the modulus:

sage: K.<t> = GF(13^5, modulus=x^5+4*x-2); K
Finite Field in t of size 13^5
sage: t.minpoly()
x^5 + 4*x + 11

Next, we can define elliptic curves with EllipticCurve, passing in either \([a_4, a_6]\) or \([a_1, a_2, a_3, a_4, a_6]\) for the second argument:

sage: set_random_seed(1337)
sage: K = GF(103)
sage: E = EllipticCurve(K, [3, 5]); E
Elliptic Curve defined by y^2 = x^3 + 3*x + 5 over Finite Field of size 103
sage: E.order()
96
sage: E.random_point()
(40 : 33 : 1)
sage: E(0) == E(0, 1, 0) == E.zero()
True

We can also construct them by their \(j\)-invariant, though we won’t need this for now:

sage: E = EllipticCurve(K, j=17); E
Elliptic Curve defined by y^2 = x^3 + 20*x + 87 over Finite Field of size 103

Here is a fun one. When \(E\) is defined over a finite field \(K\), it’ll be a finite abelian group, in which case Sage provides a really nice method, called .abelian_group, that allows you work with the abelian group \(E(K)\) directly. For example, we can compute the group structure of the \(k\)-torsion points by taking the \(k\)-torsion subgroup of \(E(K)\). Let’s start with the following curve:

sage: E = EllipticCurve(GF(103), [1, 0])
sage: G = E.abelian_group(); G
Additive abelian group isomorphic to Z/104 embedded in Abelian group of points on Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 103
sage: G.torsion_subgroup(5)
Trivial group embedded in Abelian group of points on Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 103
sage: G.torsion_subgroup(2)
Additive abelian group isomorphic to Z/2 embedded in Abelian group of points on Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 103

As we can see, it has no points of order \(5\), since \(|E(\mathbb{F}_{103})| = 104\) (it is supersingular, as described below). Despite \(2^2\) dividing \(104\), the curve still doesn’t have full \(2\)-torsion, which we recall should look like \((\mathbb{Z} / 2\mathbb{Z})^2\). This is because we have to take a quadratic extension:

sage: E = E.change_ring(GF(103^2))
sage: G = E.abelian_group(); G
Additive abelian group isomorphic to Z/104 + Z/104 embedded in Abelian group of points on Elliptic Curve defined by y^2 = x^3 + x over Finite Field in z2 of size 103^2
sage: G.torsion_subgroup(5)
Trivial group embedded in Abelian group of points on Elliptic Curve defined by y^2 = x^3 + x over Finite Field in z2 of size 103^2
sage: G.torsion_subgroup(2)
Additive abelian group isomorphic to Z/2 + Z/2 embedded in Abelian group of points on Elliptic Curve defined by y^2 = x^3 + x over Finite Field in z2 of size 103^2
sage: P, Q = _.gen(0).element(), _.gen(1).element(); P, Q
((96*z2 + 55 : 0 : 1), (7*z2 + 48 : 0 : 1))
sage: P.weil_pairing(Q, 2)
102

What are Isogenies?

In the first post, we defined isogenies as a kind of group homomorphism, but didn’t go into much details about their properties or how to compute them. I hope that experimenting with the concept in Sage would help.

Let’s begin by investing our curves \(E : y^2 = x^3 + x\) and \(E’ : y^2 = x^3 - 4x\), say over \(K = \mathbb{Q}\). In the post, I told you that there is a degree-\(2\) isogeny between them. Currently, we cannot directly construct an isogeny between two isogenous curves in Sage - well, this is the isogeny path problem, which is meant to be hard I believe it’s easy in this case, since we are working over global fields, but nevermind. ! However, Sage implements the method .isogenies_prime_degree that enumerate all separable prime degree This requirement will be removed in #35949 isogenies from a curve. We can extract the desired isogeny and compute its properties.

sage: E0 = EllipticCurve(QQ, [1, 0])
sage: E0.isogenies_prime_degree(2)
[Isogeny of degree 2 from Elliptic Curve defined by y^2 = x^3 + x over Rational Field to Elliptic Curve defined by y^2 = x^3 - 4*x over Rational Field]
sage: phi = _[0]; phi.domain(), phi.codomain()
(Elliptic Curve defined by y^2 = x^3 + x over Rational Field,
 Elliptic Curve defined by y^2 = x^3 - 4*x over Rational Field)
sage: phi.degree()
2
sage: phi.rational_maps()
((x^2 + 1)/x, (x^2*y - y)/x^2)

We can verify the rational maps indeed map onto the codomain curve:

sage: fx, fy = _
sage: R.<x, y> = ZZ[]
sage: (fy^2 - (fx^3 - 4*fx)).numerator() / (y^2 - (x^3 + x))
x^4 - 2*x^2 + 1

As there is no remainder, we see that \(f_y^2 - (f_x^3 - 4f_x) \in \langle y^2 - (x^3 + x) \rangle\), verifying the claim.

Evaluating an isogeny is easy, we just call phi(P) where P is the point:

sage: phi(E(0, 0))
(0 : 1 : 0)

We can also compute the dual isogeny, which goes the opposite way:

sage: phi_hat = phi.dual(); phi_hat
Isogeny of degree 2 from Elliptic Curve defined by y^2 = x^3 - 4*x over Rational Field to Elliptic Curve defined by y^2 = x^3 + x over Rational Field
sage: phi_hat.domain(), phi_hat.codomain()
(Elliptic Curve defined by y^2 = x^3 - 4*x over Rational Field,
 Elliptic Curve defined by y^2 = x^3 + x over Rational Field)
sage: phi_hat.degree()
2
sage: phi_hat.rational_maps()
((1/4*x^2 - 1)/x, (1/8*x^2*y + 1/2*y)/x^2)

As noted in the first post, these rational maps are not isomorphisms, so \(\varphi\) and \(\widehat{\varphi}\) cannot be inverses. What do they compose to then? As it turns out, we have \(\varphi \circ \widehat{\varphi} = \widehat{\varphi} \circ \varphi = [\deg(\varphi)]\), the multiplication-by-\(\deg(\varphi)\) map! In Sage, we can compose isogenies by multiplying them. We can check equality of isogenies by using ==, which compares the kernel of the isogenies.

sage: comp = phi_hat * phi; comp
Composite morphism of degree 4 = 2^2:
  From: Elliptic Curve defined by y^2 = x^3 + x over Rational Field
  To:   Elliptic Curve defined by y^2 = x^3 + x over Rational Field
sage: comp == E0.scalar_multiplication(2)
True

sage: comp_inv = phi * phi_hat; comp_inv
Composite morphism of degree 4 = 2^2:
  From: Elliptic Curve defined by y^2 = x^3 - 4*x over Rational Field
  To:   Elliptic Curve defined by y^2 = x^3 - 4*x over Rational Field
sage: comp_inv == comp_inv.domain().scalar_multiplication(2)
True

Let’s go through isogenies over finite fields quickly as well. First, we want to find the isogeny \(E_0 \to E_1\). If we don’t know the degree of the isogeny, we can also enumerate all isogenies with degree bounded above by a constant, which is by default \(30\):

sage: K = GF(17)
sage: E0 = EllipticCurve(K, [1, 0])
sage: E1 = EllipticCurve(K, [6, 3])
sage: E0.is_isogenous(E1)
True

sage: E0.isogenies_prime_degree()
[Isogeny of degree 2 from Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17 to Elliptic Curve defined by y^2 = x^3 + 13*x over Finite Field of size 17,
 Isogeny of degree 2 from Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17 to Elliptic Curve defined by y^2 = x^3 + 11*x + 5 over Finite Field of size 17,
 Isogeny of degree 2 from Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17 to Elliptic Curve defined by y^2 = x^3 + 11*x + 12 over Finite Field of size 17,
 Isogeny of degree 5 from Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17 to Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17,
 Isogeny of degree 5 from Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17 to Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17,
 Isogeny of degree 13 from Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17 to Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17,
 Isogeny of degree 13 from Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17 to Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17,
 Isogeny of degree 17 from Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17 to Elliptic Curve defined by y^2 = x^3 + 16*x over Finite Field of size 17,
 Isogeny of degree 29 from Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17 to Elliptic Curve defined by y^2 = x^3 + 13*x over Finite Field of size 17,
 Isogeny of degree 29 from Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17 to Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17]

We don’t see our curve \(E_1\) as a codomain of an isogeny. This is because isogenies are only unique up to isomorphism. Thus, we might have to post-compose an isogeny with an isomorphism, i.e. a degree-\(1\) isogeny, to get the desired isogeny \(E_0 \to E_1\).

sage: for phi in E0.isogenies_prime_degree():
....:     if phi.codomain().is_isomorphic(E1):
....:         psi = phi.codomain().isomorphism_to(E1)
....:         phi = psi * phi
....:         print("Found:", phi)
....:         break
Found: Isogeny of degree 2 from Elliptic Curve defined by y^2 = x^3 + x over Finite Field of size 17 to Elliptic Curve defined by y^2 = x^3 + 6*x + 3 over Finite Fiel
d of size 17
sage: phi.rational_maps()
((-4*x^2 - x + 8)/(x - 4), (2*x^2*y + x*y + 2*y)/(x^2 - 8*x - 1))

Furthermore, we can compute the kernel polynomial and the kernel of the isogeny:

sage: f = phi.kernel_polynomial(); f
x + 13
sage: f.roots(multiplicities=False)
[4]
sage: ker = [E0(0)] + [P for x in _ for P in E0.lift_x(x, all=True)]; ker
[(0 : 1 : 0), (4 : 0 : 1)]
sage: all(phi(P) == 0 for P in ker)
True

We can also verify Tate’s isogeny theorem, which states that two elliptic curves over a finite field are isogenous if and only if they have the same number of points:

sage: E0.order() == E1.order()
True

Recall that finite subgroups are in bijection with isogenies, and more specifically, any point \(P \in E\) of finite order defines an isogeny \(\varphi_P : E \to E / \langle P \rangle\). In Sage, we can compute this using the .isogeny(P) method. For example, the isogeny \(E_0 \to E_1\) has kernel \(\langle (4, 0) \rangle\), so we can recover the isogeny (including the isomorphism) as follows:

sage: P = E0(4, 0)
sage: phi_P = E0.isogeny(P, codomain=E1)
sage: phi_P == phi
True

Supersingular Elliptic Curves

I will swap the order of this section and the SIDH section, since it’s more natural to present the protocol at the end. In the first post, I referenced a theorem in Silverman’s AEC, which states that the curve \(E : y^2 = x^3 + x\) is supersingular over \(\mathbb{F}_{p^k}\) when \(p \equiv 3 \pmod{4}\). Let’s check out what special about this curve!

Firstly, supersingularity can be checked in Sage easily:

sage: E = EllipticCurve(GF(13), [1, 0])
sage: E.is_supersingular()
False

sage: p = 19
sage: E = EllipticCurve(GF(p), [1, 0])
sage: E.is_supersingular()
True

sage: E = EllipticCurve(GF((p, 5)), [1, 0])
sage: E.is_supersingular()
True

We can verify the definition we used for supersingularity, i.e. that \(|E(\mathbb{F}_q)| \equiv 1 \pmod{p}\):

sage: E.order()
2476100
sage: E.order() % p
1

In fact, we can verify the “powerful lemma” directly:

sage: p = 19
....: for k in range(1, 21):
....:     E = EllipticCurve(GF((p, k)), [1, 0])
....:     order = E.order()
....:     if k % 2 == 1:
....:         assert order == p^k + 1
....:     elif k % 4 == 0:
....:         assert order == (p^(k // 2) - 1)^2
....:     else:
....:         assert order == (p^(k // 2) + 1)^2

For those who are familiar with the term “trace of Frobenius”, an equivalent formulation of our definition is that a curve is supersingular when the trace of Frobenius is a multiple of \(p\). We can also verify this in Sage:

sage: E = EllipticCurve(GF((p, 5)), [1, 0])
sage: E.trace_of_frobenius() % p
0

Finally, just to give more examples for you to play around, here are two more supersingular curves:

sage: E = EllipticCurve(GF((10^9 + 7, 4)), [608905244, 136566674])
sage: E = EllipticCurve(GF((10^9 + 7, 4)), [249126373, 784186028])

SIDH & Detailed Protocol

Finally, we arrive at the meat of these two posts: implementing SIDH in Sage! We will break this into sub-steps, following the “Detailed protocol” given in the first post.

Public Parameters

Integers \(a = 91, b = 57, f = 1\), so \(p = 2^{91} \cdot 3^{57} - 1\)
Supersingular elliptic curve \(E_0 : y^2 = x^3 + x\) over a fixed A modulus should be chosen. We will use \(x^2 + 1\) here, i.e. \(\mathbb{F}_{p^2} \cong \mathbb{F}_p[X] / (X^2 + 1)\) \(\mathbb{F}_{p^2}\)
basis \(\color{pink} \{P_1, P_2\}\) of \(E_0[2^a]\)
basis \(\color{lime} \{Q_1, Q_2\}\) of \(E_0[3^b]\)

The first two lines are easy to implement in Sage. As hinted, we have to fix a finite field modulus that both parties use (to serialise data etc.). We will use the modulus \(x^2 + 1\), which is irreducible in \(\mathbb{F}_p[x]\).

sage: a, b, f = 91, 57, 1
sage: p = f * 2^a * 3^b - 1
sage: is_prime(p)
True
sage: K.<i> = GF((p, 2), modulus=x^2 + 1)
sage: E0 = EllipticCurve(K, [1, 0])
sage: E0.is_supersingular()
True
sage: E0.order() == (2^a * 3^b)^2
True

Next, we have to compute the basis of \(E_0[2^a]\) and \(E_0[3^b]\). Before we get on to how to compute them, let’s verified that such basis indeed exists over our \(K\). We can use the .abelian_group method mentioned above to do this:

sage: set_random_seed(1337)
sage: G = E0.abelian_group(); G
Additive abelian group isomorphic to Z/3887237936338808897317857226010584482932750535805108224 + Z/3887237936338808897317857226010584482932750535805108224 embedded in Abelian group of points on Elliptic Curve defined by y^2 = x^3 + x over Finite Field in i of size 3887237936338808897317857226010584482932750535805108223^2
sage: G.invariants()
(3887237936338808897317857226010584482932750535805108224,
 3887237936338808897317857226010584482932750535805108224)
sage: G.torsion_subgroup(2^a).invariants() == (2^a, 2^a)
True
sage: G.torsion_subgroup(3^b).invariants() == (3^b, 3^b)
True

To compute the basis, we have three possible approaches. The first one is by sampling random points, clearing out cofactors, and checking linear independence, which can be done via Weil pairing:

sage: set_random_seed(1337)
sage: while True:
....:     P1, P2 = [E0.random_point() * 3^b for _ in range(2)]
....:     if P1.order() == P2.order() == 2^a and P1.weil_pairing(P2, 2^a).multiplicative_order() == 2^a:
....:         break
sage: while True:
....:     Q1, Q2 = [E0.random_point() * 2^a for _ in range(2)]
....:     if Q1.order() == Q2.order() == 3^b and Q1.weil_pairing(Q2, 3^b).multiplicative_order() == 3^b:
....:         break
sage: P1, P2, Q1, Q2
((2177487706377088761481524756747874678745073661077223131*i + 1108955115155556733180492903474041841641462326710538559 : 1791483673166938029372764577113986648425901965918914656*i + 30699774988539859132 68362043012822268244150813011229842 : 1),
 (2739250155562685879061339278671728858838381488025869079*i + 3126094196227962224043822299674404648878032620948909058 : 3381381191889714566762732917934344835001374489965899011*i + 34556042134254054099 60831935147366678573638099906313380 : 1),
 (137090328070107676991555332350414673343031883693734096*i + 1183496910839253242662449756369234335816235741577417009 : 3557424710022084781321522851073720847208486707391198743*i + 153750460360469466647 0234795655106298598747171774972681 : 1),
 (3699069675618859091748960064718480148208746273630917820*i + 856719356620043489334394388698024003020547719259950713 : 1522515716428284276377450141149488374268709345468618716*i + 338163813660978483669 6067878385504842801824816883786531 : 1))

The second approach is to use the generators of the torsion subgroup of \(G\).

sage: P1, P2 = [P.element() for P in G.torsion_subgroup(2^a).gens()]
sage: Q1, Q2 = [P.element() for P in G.torsion_subgroup(2^a).gens()]
sage: P1, P2, Q1, Q2
((577034547453098423246295441600941274636367160973518537*i + 3856700886061726851282173908071324039783471523466286373 : 3741572677216800689171746726408379465941491075911200408*i + 1045538535559434609738118540742109516481511018452719275 : 1),
 (1792801238231187732544887501137930663088705953055645123*i + 643158215428552082131388892999971795151344302067542910 : 3394680284776828087973190048374425699325648809295513339*i + 3266974396952170994613186706071148999637635735807116983 : 1),
 (577034547453098423246295441600941274636367160973518537*i + 3856700886061726851282173908071324039783471523466286373 : 3741572677216800689171746726408379465941491075911200408*i + 1045538535559434609738118540742109516481511018452719275 : 1),
 (1792801238231187732544887501137930663088705953055645123*i + 643158215428552082131388892999971795151344302067542910 : 3394680284776828087973190048374425699325648809295513339*i + 3266974396952170994613186706071148999637635735807116983 : 1))

Finally we can use yx7 / yyyyx4 / Lorenz’s insane amount of work in Sage, and use the .torsion_basis method. There’s also ongoing work on this method:

sage: set_random_seed(1337)
sage: P1, P2 = E0.torsion_basis(2^a)
sage: Q1, Q2 = E0.torsion_basis(3^b)
sage: P1, P2, Q1, Q2
((577034547453098423246295441600941274636367160973518537*i + 3856700886061726851282173908071324039783471523466286373 : 3741572677216800689171746726408379465941491075911200408*i + 1045538535559434609738118540742109516481511018452719275 : 1),
 (1792801238231187732544887501137930663088705953055645123*i + 643158215428552082131388892999971795151344302067542910 : 3394680284776828087973190048374425699325648809295513339*i + 3266974396952170994613186706071148999637635735807116983 : 1),
 (1191174884164132194308790043137117594635624785804862883*i + 733802895733457424111061607997835008793783970107967798 : 1248881161673866277498032215180291816730526556395258723*i + 3840076495717715782127427065540057576594362127145396622 : 1),
 (320841087602825159953017126458292055875254980151601290*i + 769895406929582352581402433656853988805104906090615858 : 1459079840408465712036431191452928440375394532953631833*i + 2016392331831519635604392895272063574514088301502574218 : 1))

Moving forward, we will be using the output of the .torsion_basis.

saviour or something…

Private Keys + Key Exchange: Computing Kernel

Alice: Random integers \(a_1, a_2\)
Bob: Random integers \(b_1, b_2\)
0A. Alice computes \(P = [a_1]P_1 + [a_2]P_2 \in E_0\)
0B. Bob computes \(Q = [b_1]Q_1 + [b_2]Q_2 \in E_0\)

Just sample them randomly.

sage: set_random_seed(1337)
sage: a1, a2 = [randrange(2^a) for _ in range(2)]
sage: b1, b2 = [randrange(3^b) for _ in range(2)]
sage: P = a1 * P1 + a2 * P2
sage: Q = b1 * Q1 + b2 * Q2
sage: P, Q
((1992823321778690380948063358764714540659745500840523635*i + 3537583408445532802804715084327215103301485860402556989 : 742353139296217585692311628634957684803516890167916628*i + 3263560868246409437739750846330397558804247487222510086 : 1),
 (2668499403728360249058277893943034924675249423627051530*i + 3173042250818142638081040193006933616086311464671574899 : 2498749800477666430507377755913063256893874579453495848*i + 1909041617670705834387051814066489298685795205391885923 : 1))

Key Exchange: Isogeny Chain 1

1A. Alice computes \(\varphi_A : E_0 \to E_A\) with \(E_A = E_0 / \langle P \rangle\)
2A. Alice computes \(Q_1’ = \varphi_A(Q_1)\) and \(Q_2’ = \varphi_A(Q_2)\)
1B. Bob computes \(\varphi_B : E_0 \to E_B\) with \(E_B = E_0 / \langle Q \rangle\)
2B. Bob computes \(P_1’ = \varphi_B(P_1)\) and \(P_2’ = \varphi_B(P_2)\)

This part is easy for us… right?

sage: set_random_seed(1337)
sage: phiA = E0.isogeny(P) # hangs

Oh right, we have to implement the isogeny chain method. Just like above, there are two methods. The first one is to implement the isogeny chain decomposition by hand. As explained in the post, the idea is to ``kill’’ the kernel subgroup in multiple iterations, each iteration only killing a small prime order sub-subgroup, and pushing the kernel subgroup through the small isogeny. It can be implemented in a few lines:

sage: def isogeny_chain(P, l):
....:     k = ZZ(log(P.order(), l))
....:     assert is_prime(l) and P.order() == l^k and P.order() % l == 0
....:     E = P.curve()
....:     phi = E.identity_morphism()
....:     for i in range(1, k + 1):
....:         psi = E.isogeny(l^(k - i) * P)
....:         phi = psi * phi
....:         E = psi.codomain()
....:         P = psi(P)
....:     return phi
sage: phiA = isogeny_chain(P, 2)
sage: EA, Q1A, Q2A = phiA.codomain(), phiA(Q1), phiA(Q2)
sage: phiB = isogeny_chain(Q, 3)
sage: EB, P1B, P2B = phiB.codomain(), phiB(P1), phiB(P2)
sage: EA, EB
(Elliptic Curve defined by y^2 = x^3 + (1822530699521883672143422779388649334237212896329791112*i+3568662577460411800587272558215317090731602035089121703)*x + (16844175111449769925773252405268442467211 67502347196824*i+1122161097308449692861069425854966522057551911923384786) over Finite Field in i of size 3887237936338808897317857226010584482932750535805108223^2,
 Elliptic Curve defined by y^2 = x^3 + (2661692422410715315385047943543850875192251068967742147*i+324869182601077271126976073334262797771038412225734649)*x + (327079425764139154762187029117826807006412 4585550544979*i+3565943788469482713346364479211820304421598432274371339) over Finite Field in i of size 3887237936338808897317857226010584482932750535805108223^2)

Alternatively, we can pray to our almighty saviour again, who blesses us with a keyword argument at .isogeny, namely the algorithm= argument. If we pass in algorithm="factored", then the isogeny chain algorithm above is used automatically.

sage: phiA = E0.isogeny(P, algorithm="factored")
sage: EA, Q1A, Q2A = phiA.codomain(), phiA(Q1), phiA(Q2)
sage: phiB = E0.isogeny(Q, algorithm="factored")
sage: EB, P1B, P2B = phiB.codomain(), phiB(P1), phiB(P2)
sage: EA, EB
(Elliptic Curve defined by y^2 = x^3 + (1822530699521883672143422779388649334237212896329791112*i+3568662577460411800587272558215317090731602035089121703)*x + (16844175111449769925773252405268442467211 67502347196824*i+1122161097308449692861069425854966522057551911923384786) over Finite Field in i of size 3887237936338808897317857226010584482932750535805108223^2,
 Elliptic Curve defined by y^2 = x^3 + (2661692422410715315385047943543850875192251068967742147*i+324869182601077271126976073334262797771038412225734649)*x + (327079425764139154762187029117826807006412 4585550544979*i+3565943788469482713346364479211820304421598432274371339) over Finite Field in i of size 3887237936338808897317857226010584482932750535805108223^2)

Hurray!

Key Exchange: Isogeny Chain 2

(Step 3A, 3B omitted for obvious reasons.)

4A. Alice computes \(P’ = [a_1]P_1’ + [a_2]P_2’ \in E_B\)
5A. Alice computes \(\color{gold} E_{AB} = E_B / \langle P’ \rangle\)
4B. Bob computes \(Q’ = [a_1]Q_1’ + [a_2]Q_2’ \in E_A\)
5B. Bob computes \(\color{gold} E_{BA} = E_A / \langle Q’ \rangle\)

Step 4 is easy and step 5 is the same as step 1 above.

sage: PB = a1 * P1B + a2 * P2B
sage: EAB = EB.isogeny(PB, algorithm="factored").codomain()
sage: QA = b1 * Q1A + b2 * Q2A
sage: EBA = EA.isogeny(QA, algorithm="factored").codomain()
sage: EAB, EBA
(Elliptic Curve defined by y^2 = x^3 + (1365367514626378052504473740715639769314121188091153483*i+2651058611096873241886541435045653539489760673562077221)*x + (17394535589046134493683744388964354757527 78617573618578*i+3312769349534295346409154904928871289729503540269712726) over Finite Field in i of size 3887237936338808897317857226010584482932750535805108223^2,
 Elliptic Curve defined by y^2 = x^3 + (1365367514626378052504473740715639769314121188091153483*i+2651058611096873241886541435045653539489760673562077221)*x + (17394535589046134493683744388964354757527 78617573618578*i+3312769349534295346409154904928871289729503540269712726) over Finite Field in i of size 3887237936338808897317857226010584482932750535805108223^2)

Woohoo! We got the same curve.

Key Exchange: Deriving Secret

As mentioned, such isogenies are only defined up to isomorphisms, and there is a chance that EAB and EBA are only isomorphic, not equal. Hence, we need an invariant of elliptic curves that does not change under isomorphisms. The typical choice is the \(j\)-invariant, for which the definition does not matter too much to us, but it is easily computable in Sage:

sage: jAB = EAB.j_invariant(); jAB
3378759616299362693936029133632292064252451552090231734*i + 2483820189230540397711199906295045169174476636348752124
sage: jBA = EBA.j_invariant(); jBA
3378759616299362693936029133632292064252451552090231734*i + 2483820189230540397711199906295045169174476636348752124
sage: jAB == jBA
True

And there you go! This is how you implement the SIDH scheme from scratch in Sage. Hope this blog helps you understand how to do isogeny-based cryptography in Sage!

For no reason, here’s a song I like recently.