Elliptic butterflies

Jean-Marc Couveignes, Reynald Lercier

► To cite this version:

Jean-Marc Couveignes, Reynald Lercier. Elliptic butterflies. 2025. hal-05342789

HAL Id: hal-05342789

https://hal.science/hal-05342789v1

Preprint submitted on 3 Nov 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
ELLIPTIC BUTTERFLIES

JEAN-MARC COUVEIGNES AND REYNALD LERCIER

*In memory of Tony Ezome*

**ABSTRACT.** We study natural evaluation and interpolation problems for elliptic functions and prove that they allow a recursive treatment using a variant of classical butterflies first introduced by Gauss. We deduce the existence of straight-line programs with complexity in $d \log(d)$ for these problems and present applications to finite field arithmetic, coding theory and cryptography.

<figure><img src="image_2.png" /><figcaption>Butterflyfish (freely adapted from [Sha24])</figcaption></figure>

CONTENTS

<table><tr><td>1. Introduction</td><td>2</td></tr><tr><td>2. Summary and organization of the article</td><td>4</td></tr><tr><td>3. Algorithms and data structures</td><td>5</td></tr><tr><td>4. Elliptic functions</td><td>6</td></tr><tr><td>5. Bases for elliptic function field extensions</td><td>7</td></tr><tr><td>6. A Vélu isogeny of degree two</td><td>8</td></tr><tr><td>7. Evaluation using elliptic butterflies</td><td>8</td></tr><tr><td>8. Interpolation</td><td>12</td></tr><tr><td>9. Cyclic bidagonal linear systems</td><td>14</td></tr><tr><td>10. Reduction</td><td>15</td></tr><tr><td>11. Multiplication in the residue ring of a fiber</td><td>19</td></tr><tr><td>12. Implementation and experimental results</td><td>22</td></tr><tr><td>13. Elliptic bases for finite field extensions</td><td>22</td></tr><tr><td>14. Elliptic Goppa codes</td><td>24</td></tr><tr><td>15. Elliptic LWE cryptography</td><td>25</td></tr><tr><td>References</td><td>27</td></tr></table>
# 1. INTRODUCTION

Let $\mathbf{K}$ be a field and let $X$ be a smooth, projective, absolutely integral curve over $\mathbf{K}$. Let $D$ be a divisor on $X$. We denote by $\mathcal{L}(D)$ the associated linear space and call $d$ its dimension. Let $b_0, b_1, \dots, b_{d-1}$ be $\mathbf{K}$-points on $X$ and not in the support of $D$. We consider the evaluation map

$$ (1) \qquad \begin{aligned} \text{ev}: \mathcal{L}(D) &\longrightarrow \mathbf{K}^d, \\ f &\longmapsto f(b_l)_{0 \le l \le d-1}. \end{aligned} $$

which we assume to be invertible. We are interested in the complexity of evaluation (evaluating the map ev) and interpolation (evaluating the reciprocal map ev⁻¹). Elements in $\mathbf{K}^d$ are represented by their entries. Elements in $\mathcal{L}(D)$ are represented by their coordinates in a basis $u = (u_0, u_1, \dots, u_{d-1})$ that must be specified. It is traditional to use bases satisfying one or several among the following properties:

(1) the sum of the degrees of the $u_l$ is minimal among all bases,

(2) the $u_l$ are eigenvalues of some automorphism of $X$,

(3) the $u_l$ are permuted by some automorphism of $X$.

One looks for a straight-line program that solves the evaluation problem using additions and multiplications by constants in $\mathbf{K}$. The complexity of evaluation is then the minimum number of operations in such a straight-line program. One defines the complexity of interpolation in the same way.

The classical and most important case is when $X = \mathbf{P}^1$ and $D = (d-1)[\infty]$ and $\mathcal{L}(D)$ is the space of polynomials of degree $\le d-1$ in one variable $x$. In this situation, and assuming that $\mathbf{K}$ contains a primitive $d$-th root of unity, the monomial basis $(1, x, x^2, \dots, x^{d-1})$ satisfies the first two conditions above. If further $d$ is a power of two, which implies that the characteristic of $\mathbf{K}$ is not two, then the complexities of evaluation and interpolation are bounded by

$$ (2) \qquad (3/2)d \log(d) + 1. $$

This is achieved using the well-known fast Fourier transform or FFT, first discovered by Gauss. See [BCS97, Theorem 2.6] and the historical account in [HJB85].

Bostan and Schost give in [BS05] a range of nice upper bounds for the complexity of evaluation and interpolation when the basis is the monomial basis or the Newton basis, and the evaluation set is an arithmetic sequence, a geometric sequence, or any set. They make no assumption on $\mathbf{K}$. For example, without making any restriction on $\mathbf{K}$ nor on the evaluation set, the complexity of evaluation and interpolation in the monomial basis is bounded by an absolute constant times $M(d) \log(d)$ where $M(d)$ is the complexity of multiplying two polynomials of degree $\le d$ with coefficients in $\mathbf{K}$. Since the latter is bounded by a constant times $d \log(d) \log(\log(d))$ operations in $\mathbf{K}$, according to a result of Schönhage and Strassen [BCS97, Theorem 2.13], we end up with a general bound of a constant times

$$ (3) \qquad d \log^2(d) \log(\log(d)) $$

for evaluation and interpolation in the monomial basis, at a general set of points.

Comparing Equations (3) and (2) we see that we save a factor $\log(d) \log(\log(d))$ when $d$ is a power of two, $\mathbf{K}$ contains a primitive $d$-th root of unity, $\omega_d$, and the evaluation set is the set of all $d$-th roots of 1. This gain is not marginal as is demonstrated by the many implementations and applications of FFT. Indeed general evaluation and interpolation methods often use specific ones as subroutines.

The key idea behind FFT is to use the involution $x \mapsto -x$ to decompose the polynomial to be evaluated $P(x) = \sum_{0 \le l \le d-1} p_l x^l$ as a sum

$$ P(x) = P^{+}(x) + P^{-}(x) $$
where

$$P^{+}(x) = \sum_{0 \le 2l \le d-1} p_{2l} x^{2l} \quad \text{and} \quad P^{-}(x) = \sum_{0 \le 2l+1 \le d-1} p_{2l+1} x^{2l+1} = x \sum_{0 \le 2l+1 \le d-1} p_{2l+1} x^{2l}$$

are the even and odd parts of $P(x)$. This decomposition reduces the evaluation of $P(x)$ to two similar problems of halved size, allowing a recursive approach that is classically illustrated using diagrams called butterflies (see Figure 1).

<figure><img src="image_3.png" /><figcaption>FIGURE 1. An 8-point FFT butterfly</figcaption></figure>

In this work we define and study a family of evaluation and interpolation problems on an elliptic curve $E$ having a $K$-rational point $t$ of order $d$, a power of two. We show that there are natural analogues of butterflies in this context. We replace roots of unity by multiples of the point $t$ and the associated translations, while polynomials are replaced by functions in the space $\mathcal{L}(\langle t \rangle)$ equipped with a basis $u$ of functions of degree $\le 2$. In particular the role of the involution $x \mapsto -x$ is now played by the translation by $T = (d/2)t$, a point of order two. This results in complexity bounds $Qd \log(d)$ for the evaluation and interpolation problems<fnref n="1" />.

To be complete we should say that evaluation and interpolation are not quite enough for our purposes. We need another linear map on functions that we call reduction and which we now define in a general context. We now have three divisors $D, E$ and $B$ on the curve $X$. We denote by $d$ the dimension of $\mathcal{L}(D)$. We assume that $B$ is effective of degree $d$. We assume that $D$ and $B$ have disjoint support. We denote by $\mathbf{L}$ the residue ring at $B$. This is a $K$-algebra of dimension $d$. We assume that the residue map

$$\mathrm{ev}_D^B: \mathcal{L}(D) \longrightarrow \mathbf{L}, \\ f \longmapsto f \bmod B.$$

is invertible. We also assume that $E$ and $B$ have disjoint support and define

$$\mathrm{ev}_E^B: \mathcal{L}(E) \longrightarrow \mathbf{L}, \\ f \longmapsto f \bmod B.$$

the residue map. For every function $F$ in $\mathcal{L}(E)$ there exists a unique function $f$ in $\mathcal{L}(D)$ such that

$$F \equiv f \pmod{B} \quad \text{that is} \quad \mathrm{ev}_E^B(F) = \mathrm{ev}_D^B(f) \in \mathbf{L}.$$

<footnote n="1">¹Note that everywhere in this article, the notation $\mathcal{Q}$ stands for a positive absolute constant. Any sentence containing this symbol becomes true if the symbol is replaced in every occurrence by some large enough real number.</footnote>
We choose a basis for $\mathcal{L}(D)$ and a basis for $\mathcal{L}(E)$. The reduction problem takes as input the coordinates of $F$ and returns the coordinates of $f$. We will explain how to efficiently solve a family of reduction problems on elliptic curves. In the situation we are interested in, the divisor $E$ is $2D$ because we need to reduce the product of two functions in $\mathcal{L}(D)$.

We will present three applications of these fast evaluation, interpolation, and reduction algorithms. The first one is the construction of normal bases for extensions of finite fields of degree a power of two, allowing particularly fast multiplication. The second one is the description of $[d, d/2, d/2 + 1]$-error correcting codes that can be encoded and checked in time $Qd \log(d)$. The third one is a discussion on the relevance of using residue rings on elliptic curves in the context of LWE cryptography.

To conclude this introduction we recall that there is a long tradition of interpolating on algebraic curves especially over finite fields. Among the major results, one may mention linear upper bounds on the bilinear complexity of multiplication in finite field extensions by Chudnovsky and Chudnovsky [CC88], Shparlinski, Tsfasman, Vlăduț [STV92], Shokrollahi [Sho92], Ballet and Rolland [BR04; Bal99], Chaumine [Cha08], Randriambolona [Ran12] and others. Another major achievement was the construction of excellent codes by Goppa, Tsfasman, Vlăduț, Zink, Ihara, García, Stichtenoth [GS95; Gop81; Gop82; Iha81; TVZ82]. Chudnovsky and Chudnovsky express in [CC89, Section 6] the intuition that elliptic curves with smooth order could be of some use to accelerate some specific interpolation problems in genus 0. Using similar ideas, Ben-Sasson, Carmon, Kopparty, and Levit show in [Ben+23] how to quickly compute the values of a polynomial at a set $S'$ from its value at another set $S$ when $S$ and $S'$ are specific subsets of finite fields i.e. $x$-coordinates of well chosen points on a well chosen elliptic curve. Starting from Chudnovsky's intuition we prefer to stay in the elliptic world where all these ideas and objects originate from, and develop an explicit analogue of butterflies in this context. We also prefer to invoke evaluation and interpolation rather than Fourier transform, because in our humble opinion, a Fourier transform is something different. See e.g. [BCS97, Section 13.5] or [CG23, Section 6.1].

## 2. SUMMARY AND ORGANIZATION OF THE ARTICLE

In this article we consider specific evaluation and interpolation problems in genus one. We consider an elliptic curve $E$ over a field $K$ of odd characteristic. We let $d$ be a power of two and we assume that there is a point $t$ of order $d$ in $E(K)$. We let

$$D = \sum_{0 \le l \le d-1} [lt]$$

be the degree $d$ divisor associated to the group $\langle t \rangle$ generated by $t$. We take $b$ to be a $K$-rational point on $E$ such that

$$db \neq 0.$$

We set

$$b_l = b + lt \quad \text{for } l \in [0, d-1]$$

and we check that none of these $d$ points is in $\langle t \rangle$. The evaluation map $ev$ in Equation (1) is well defined and invertible.

To complete the description of the interpolation problem we are interested in, we need interesting bases for $\mathcal{L}(D)$. These bases are made of functions of degree $\le 2$. In Section 4 we recall nice properties of degree two functions on an elliptic curve. We use these functions in Section 5 to define two nice bases for $\mathcal{L}(D)$. One of these bases satisfies the first condition stated in the introduction. The sum of the degrees of its elements is $2d - 2$. The other basis satisfies the third condition. It is invariant by translation by $t$. And the sum of the degrees of its elements is $2d$. We show that base
change between these two bases requires no more than $Qd$ operations in $K$. This means that these two bases are essentially the same basis from the point of view of complexity theory.

The automorphism $x \mapsto -x$ is crucial to define classical butterflies. Its elliptic counterpart is translation by $T$ where $T = (d/2)t$ is the unique point of order two in $\langle t \rangle$. We study this translation and the corresponding quotient $E' = E/\langle T \rangle$ in Section 6. This enables us to describe the recursive evaluation algorithm in Section 7. An important point is that we use a degree two function $\theta$ that is odd for translation by $T$. We decompose $f \in \mathcal{L}(D)$ as

$$f = f^+ + f^-$$

where $f^+$ is even and $f^-$ is odd. We then set

$$f^0 = f^+ \quad \text{and} \quad f^1 = \theta f^-$$

and we notice that both are functions on the quotient $E'$. In the classical setting the role of $\theta$ is played by the variable $x$ itself.

The interpolation algorithm is already a bit more tricky than evaluation. It is described in Section 8. A difficulty is to solve a linear system of a very special and classical form: a cyclic bid diagonal system. We recall in Section 9 what is needed to solve efficiently such a system in our context.

If we multiply two functions $f_1$ and $f_2$ in $\mathcal{L}(D)$ we obtain a function $F$ in $\mathcal{L}(2D)$. Given such a function $F$, there is a unique function $f$ in $\mathcal{L}(D)$ such that $f(b_l) = F(b_l)$ for every $l \in [0, d-1]$. Finding $f$ once given $F$ is a reduction problem. In Section 10 we define and solve a natural generalization of this problem. The resulting efficient multiplication algorithm in the residue ring at $B$ is presented in Section 11. Section 12 documents a public implementation of the algorithms presented in this article and provides experimental results.

A first application, presented in Section 13, is the construction of normal bases for finite field extensions of degree a power of two, allowing multiplication in time $Qd \log(d)$. We build on the construction in [CL09] of an equivariant version of Chudnovsky's interpolation method [CC88], and we use our complexity estimates for evaluation, interpolation and reduction in this context.

In Section 14, we construct Goppa codes in genus one obtained by restricting the evaluation map $ev$ to the subspace of functions that are invariant by the elliptic involution $P \mapsto -P$ on $E$. The resulting code is a genus zero Goppa code. The advantage of seeing it as a genus one code is that we benefit from the acceleration due to the presence of the automorphism group $\langle t \rangle$. This produces MDS $[d, d/2, d/2+1]$-codes that can be encoded and checked in time $Qd \log(d)$ and decoded up to half the minimum distance at the expense of $Qd \log^2(d) \log(\log(d))$ operations and comparisons in $K$. The only constraint for the existence of such codes is that the size of $K$ be large enough with respect to $d$. This extends Lacan and Soro's result [SL10] for Fermat fields.

The purpose of the final Section 15 is to discuss the advantages and limitations of using elliptic residue rings in the context of LWE-cryptography.

All the straight-line programs presented in this work are summarized in pseudo-code right after their algorithmic description. An implementation of these algorithms and the underlying data structures are introduced in Section 3 and Section 12.

### 3. ALGORITHMS AND DATA STRUCTURES

The algorithms presented in this article are designed to manipulate elements in vector spaces. Vectors are typically coordinate vectors of functions in well-chosen bases of Riemann-Roch spaces, or evaluation vectors of functions at well chosen points. We therefore use the classical $d$-dimensional vector notation $\vec{a}$ and $(a_l)_{l \in \mathbb{Z}/d\mathbb{Z}}$ for inputs, outputs, variables and constants.
For ease of reading, our pseudo-code descriptions of these algorithms are given in a recursive form. However, like the classical butterfly formulation of FFT algorithms à la Cooley–Tukey, they can be very easily unrolled to yield iterative functions. In fact, our algorithms are all straight-line programs: apart from precomputations of constants, they contain no loops or conditionals, and consist of a sequence of operations, each applied to previously computed elements.

Beyond the usual operations of addition and scalar multiplication, denoted by “+” and “:”, we introduce the following notations for other operations on vectors.

* We denote by $\langle \vec{a}, \vec{b} \rangle$ the scalar product of $\vec{a}$ and $\vec{b}$.

* We denote by $\sigma(\vec{a})$ the cyclic rotation, i.e. the vector $(a_{\sigma(l)})_{l \in \mathbb{Z}/d\mathbb{Z}}$, where $\sigma(l) = l - 1$.

* We denote by $e + \vec{a}$ the addition by a constant, i.e. the vector $(e + a_l)_{l \in \mathbb{Z}/d\mathbb{Z}}$.

* We denote by $\vec{a} \oplus \vec{b}$ the direct sum, i.e. the concatenation of $\vec{a}$ and $\vec{b}$.

* We denote by $\vec{a} \diamond \vec{b}$ the component-wise product of $\vec{a}$ and $\vec{b}$.

* We denote by $\vec{a}^{\text{inv}} = (a_l^{-1})_{l \in \mathbb{Z}/d\mathbb{Z}}$ the inverse vector of $\vec{a}$ for the product $\diamond$.

## 4. ELLIPTIC FUNCTIONS

We recall a few formulae from [CL09] regarding low-degree functions on elliptic curves. We let $\mathbf{K}$ be a field and $E$ the elliptic curve over $\mathbf{K}$ defined by the Weierstrass equation

$$Y^2Z + a_1XYZ + a_3YZ^2 = X^3 + a_2X^2Z + a_4XZ^2 + a_6Z^3.$$

Call $O = (0, 1, 0)$ the point at infinity. We set $x = X/Z$, $y = Y/Z$ and $z = -x/y = -X/Y$. The Taylor expansions of $x$ and $y$ at $O$ in the local parameter $z$ are

$$x = \frac{1}{z^2} - \frac{a_1}{z} - a_2 - a_3z + O(z^2),$$

$$y = -\frac{1}{z^3} + \frac{a_1}{z^2} + \frac{a_2}{z} + a_3 + O(z).$$

If $A$ is a point on $E$, we denote by $\tau_A$ the translation by $A$. We denote by

$$z_A = z \circ \tau_{-A}$$

the composition of $z$ with the translation by $-A$. This is a local parameter at $A$. We define $x_A$ and $y_A$ in a similar way. If $A$ and $B$ are two distinct points on $E$, we denote by $u_{A,B}$ the function on $E$ defined as

$$u_{A,B} = \frac{y_A - y(A - B)}{x_A - x(A - B)}.$$

It has polar divisor $[A] + [B]$. It is invariant by the involution exchanging $A$ and $B$,

$$u_{A,B}(A + B - P) = u_{A,B}(P).$$

If $A$, $B$ and $C$ are three points on $E$ we denote by $\Gamma(A, B, C)$ the slope of the secant (resp. tangent) to $E$ going through $C - A$ and $A - B$. As a function on $E^3$ it is well defined for any three points $A$, $B$, $C$ such that $\#\{A, B, C\} \ge 2$. The Taylor expansions of $u_{A,B}$ at $A$ and $B$ are

$$ (4) \qquad \begin{aligned} u_{A,B} &= -\frac{1}{z_A} - x_A(B)z_A + (y_A(B) + a_3)z_A^2 + O(z_A^3) \\ &= \frac{1}{z_B} - a_1 + x_A(B)z_B + (y_A(B) + a_1x_A(B))z_B^2 + O(z_B^3). \end{aligned} $$
We deduce

$$u_{B,A} = -u_{A,B} - a_1,$$

$$ (5) \qquad u_{A,B} + u_{B,C} = u_{A,C} + \Gamma(A, B, C), $$

$$ \Gamma(A, B, C) = u_{B,C}(A) = u_{C,A}(B) = u_{A,B}(C) = -u_{B,A}(C) - a_1, $$

$$ u_{B,C} = u_{B,C}(A) - (x_A(C) - x_A(B))z_A + (y_A(C) - y_A(B))z_A^2 + O(z_A^3), $$

$$ u_{A,B}u_{A,C} = x_A + \Gamma(A, B, C)u_{A,C} + \Gamma(A, C, B)u_{A,B} + a_2 + x_A(B) + x_A(C), $$

$$ u_{A,B}^2 = x_A + x_B - a_1 u_{A,B} + x_A(B) + a_2. $$

## 5. BASES FOR ELLIPTIC FUNCTION FIELD EXTENSIONS

Let $\mathbf{K}$ be a field with odd characteristic and let $E$ be an elliptic curve over $\mathbf{K}$. Let $t$ be a point of order $d = 2^\delta$ in $E(\mathbf{K})$. We assume that $\delta \ge 1$. The group $\langle t \rangle$ generated by $t$ can be seen as an effective divisor of degree $d$ on $E$. We let $\mathcal{L}(\langle t \rangle)$ be the associated $\mathbf{K}$-linear space. We construct two bases for $\mathcal{L}(\langle t \rangle)$ consisting of functions of degree $\le 2$. We start with the following lemma from [CL09, Lemma 4].

**Lemma 1.** The sum $\sum_{l \in \mathbb{Z}/d\mathbb{Z}} u_{lt,(l+1)t}$ is a constant $\mathfrak{a} \in \mathbf{K}$.

*Proof.* The sum $\sum_{l \in \mathbb{Z}/d\mathbb{Z}} u_{lt,(l+1)t}$ is invariant by translations in $\langle t \rangle$. So it can be seen as a function on $E/\langle t \rangle$. As such, it has no more than one pole. Therefore it is constant. Repeated use of Equation (5) shows that

$$ \mathfrak{a} = -a_1 + \Gamma(O, t, 2t) + \Gamma(O, 2t, 3t) + \Gamma(O, 3t, 4t) + \dots + \Gamma(O, (d-2)t, (d-1)t). $$

<figure><img src="image_4.png" /></figure>

For $l$ in $\mathbb{Z}/d\mathbb{Z}$ we set

$$ u_l = u_{lt,(l+1)t} + (1 - \mathfrak{a})/d $$

and check that

$$ (6) \qquad \sum_{l \in \mathbb{Z}/d\mathbb{Z}} u_l = 1 \in \mathbf{K}. $$

We let $(v_l)_{l \in \mathbb{Z}/d\mathbb{Z}}$ be the system defined by

$$ v_0 = 1 \text{ and } v_l = u_{O,lt} \text{ for } l \neq 0 \bmod d. $$

Examination of poles shows that both $(u_l)_{l \in \mathbb{Z}/d\mathbb{Z}}$ and $(v_l)_{l \in \mathbb{Z}/d\mathbb{Z}}$ are bases of $\mathcal{L}(\langle t \rangle)$. For every integer $l$ such that $0 \le l \le d-1$ there exists a constant $\mathfrak{a}_l$ such that

$$ (7) \qquad v_l = \sum_{m=0}^{l-1} u_m + \mathfrak{a}_l. $$

This is proved using Equation (5) repeatedly. We compute $\mathfrak{a}_0 = 1$, $\mathfrak{a}_1 = (\mathfrak{a}-1)/d$ and for $2 \le l \le d-1$

$$ (8) \qquad \mathfrak{a}_l = l(\mathfrak{a}-1)/d - \Gamma(O, t, 2t) - \Gamma(O, 2t, 3t) - \Gamma(O, 3t, 4t) - \dots - \Gamma(O, (l-1)t, lt). $$

We deduce from Equations (6) and (7) that base change between $(u_l)_{l \in \mathbb{Z}/d\mathbb{Z}}$ and $(v_l)_{l \in \mathbb{Z}/d\mathbb{Z}}$ requires no more than $Qd$ additions and multiplications by a constant in $\mathbf{K}$.
## 6. A VÉLU ISOGENY OF DEGREE TWO

Under the hypotheses and notation at the beginning of Section 5 we call

$$T = 2^{\delta-1}t$$

the point of order two in the group $\langle t \rangle$. Following Vélu [Vél71; Vél78] we write $T = (x(T), y(T))$ and set

$$b_2 = a_1^2 + 4a_2,$$

$$b_4 = a_1a_3 + 2a_4,$$

$$b_6 = a_3^2 + 4a_6,$$

$$w_4 = 3x(T)^2 + 2a_2x(T) + a_4 - a_1y(T),$$

$$w_6 = 7x(T)^3 + (2a_2 + b_2)x(T)^2 + (a_4 + (-y(T)a_1 + 2b_4))x(T) + b_6,$$

$$A_4 = a_4 - 5w_4,$$

$$A_6 = a_6 - b_2w_4 - 7w_6$$

and

$$x' = x + x_T - x(T),$$

(9)

$$y' = y + y_T - y(T),$$

$$z' = -x'/y'$$

We denote by

$$\varphi: E \longrightarrow E', \\ P \longmapsto (x'(P), y'(P)).$$

the quotient by $\langle T \rangle$ isogeny. This is a degree two separable isogeny. The quotient curve $E'$ has Weierstrass equation

$$y'^2 + a_1x'y' + a_3y' = x'^3 + a_2x'^2 + A_4x' + A_6.$$

The meaning of the lemma below is that the local parameters $z$ and $z'$ are very close in the neighborhood of the origin.

**Lemma 2.** In the neighborhood of $O$ we have $z' = z + O(z^4)$ and $\frac{1}{z'} = \frac{1}{z} + O(z^2)$.

*Proof.* According to Vélu's formulae (9)

$$z' = -\frac{x'}{y'} = z \frac{1 + \frac{x_T - x(T)}{x}}{1 + \frac{y_T - y(T)}{y}} = z(1 + O(z^3)).$$

□

## 7. EVALUATION USING ELLIPTIC BUTTERFLIES

Under the hypotheses and notation of Sections 5 and 6, we let $b$ be a point in $E(\mathbf{K})$ such that

$$db \neq 0.$$

We set

$$b' = \varphi(b) \quad \text{and} \quad d' = 2^{\delta-1} = d/2$$

and we notice that

$$\hat{\varphi}(d'b') \neq 0$$
where $\hat{\varphi}$ is the dual isogeny to $\varphi$. Given $f$ in $\mathcal{L}(\langle t \rangle)$ by its coordinates in basis $u$, we want to evaluate $f$ at every point in the coset $b + \langle t \rangle$. We reduce this problem to two similar problems on $E'$ each of halved size. We show that this reduction is achieved at the expense of $Qd$ additions and multiplications by a constant in $\mathbf{K}$. Using this recursion, the evaluation of $f$ at $b + \langle t \rangle$ is achieved at the expense of $Qd \log(d)$ additions and multiplications by a constant in $\mathbf{K}$. The reduction step is the elliptic analogue of the butterfly diagram appearing in the standard FFT.

We call $O' = (0, 1, 0)$ the origin on $E'$. We denote by

$$t' = \varphi(t)$$

a point of order $d'$ on $E'$. The group $\langle t' \rangle$ can be seen as an effective divisor on $E'$ with degree $d'$. We consider the linear space $\mathcal{L}(\langle t' \rangle)$ associated with $\langle t' \rangle$. In case $d' \neq 1$ we denote by

$$u' = (u'_{l})_{l \in \mathbb{Z}/d'\mathbb{Z}} \quad \text{and} \quad v' = (v'_{l})_{l \in \mathbb{Z}/d'\mathbb{Z}}$$

the two bases of $\mathcal{L}(\langle t' \rangle)$ as constructed in Section 5. For $l$ in $\mathbb{Z}/d'\mathbb{Z}$

$$u'_{l} = u_{lt', (l+1)t'} + (1 - \mathfrak{a}')/d'$$

where

$$\mathfrak{a}' = -a_1 + \Gamma(O, t', 2t') + \Gamma(O, 2t', 3t') + \Gamma(O, 3t', 4t') + \dots + \Gamma(O, (d' - 2)t', (d' - 1)t').$$

In case $d' = 1$ we set $u'_0 = v'_0 = 1$ by convention. We say that a function $f$ on $E$ is

$$T\text{-even if } f \circ \tau_T = f \quad \text{and} \quad T\text{-odd if } f \circ \tau_T = -f.$$

We denote by

$$(10) \qquad f^+ = \frac{f + f \circ \tau_T}{2} \quad \text{and} \quad f^- = \frac{f - f \circ \tau_T}{2}$$

the $T$-even and $T$-odd parts of $f$. We shall need a small degree $T$-odd function on $E$. We let

$$\theta = u_{0,T} + \frac{a_1}{2}$$

and check that $\theta$ is $T$-odd. We deduce that the divisor of $\theta$ is

$$(\theta) = [U] + [U + T] - [O] - [T]$$

where $U$ is any 2-torsion point on $E$ not in $\langle T \rangle$. As a consequence, the product of $\theta$ and its translate $\theta \circ \tau_U$ is a non-zero constant. We write

$$f^0 = f^+ \quad \text{and} \quad f^1 = \theta f^-$$

and we check that

$$f = f^0 + \theta^{-1} f^1$$

where both $f^0$ and $f^1$ are $T$-even and thus can be seen as functions on $E'$. Indeed

$$f^0 \in \mathcal{L}(\langle t' \rangle) \quad \text{and} \quad f^1 \in \mathcal{L}(\langle t' \rangle + O' - U')$$

where

$$U' = \varphi(U).$$

Assume that $f$ is given in the basis $u = (u_l)_{l \in \mathbb{Z}/d\mathbb{Z}}$ so

$$f = \sum_{l \in \mathbb{Z}/d\mathbb{Z}} f_l u_l$$
and

$$f^0 = f^+ = \frac{1}{2} \sum_{l \in \mathbb{Z}/d\mathbb{Z}} (f_l + f_{l+d'}) u_l = \frac{1}{2} \sum_{0 \le l \le d' - 1} (f_l + f_{l+d'})(u_l + u_{l+d'}).$$

From Equation (4) and Lemma 2 we deduce that there is a constant $u$ such that for every integer $l$ in $[0, d' - 1]$

$$(11) \qquad u_l + u_{l+d'} = u'_l + u.$$

Summing out over $l$ we see that $u = 0$. So

$$(12) \qquad f^0 = f^+ = \frac{1}{2} \sum_{0 \le l \le d' - 1} (f_l + f_{l+d'}) u'_l,$$

and we can compute the coordinates of $f^0$ in the basis $u'$ at the expense of $Qd$ additions and multiplications by a constant in $K$. Similarly

$$(13) \qquad f^1 = \theta f^- = \frac{1}{2} \sum_{l \in \mathbb{Z}/d\mathbb{Z}} (f_l - f_{l+d'}) \cdot \theta \cdot u_l = \frac{1}{2} \sum_{0 \le l \le d' - 1} (f_l - f_{l+d'}) \cdot \theta \cdot (u_l - u_{l+d'}).$$

For every integer $l$ in $[1, d' - 2]$ there exist constants $b_l$ and $c_{l+1}$ in $K$ such that

$$(14) \qquad \theta \cdot (u_l - u_{l+d'}) = b_l \cdot (v'_l - v'_l(U')) + c_{l+1} \cdot (v'_{l+1} - v'_{l+1}(U')).$$

In the special case $l = 0$ we find constants $b_*$ and $c_1$, such that

$$(15) \qquad \theta \cdot (u_0 - u_{d'}) = b_* \cdot (x' - x'(U')) + c_1 \cdot (v'_1 - v'_1(U')).$$

In the special case $l = d' - 1$ we find two constants $b_{d'-1}$ and $c_*$, such that

$$(16) \qquad \theta \cdot (u_{d'-1} - u_{d-1}) = b_{d'-1} \cdot (v'_{d'-1} - v'_{d'-1}(U')) + c_* \cdot (x' - x'(U')).$$

When $d' = 1$, Equations (14), (15) and (16) simplify to

$$(17) \qquad \theta \cdot (u_0 - u_1) = (b_* + c_*)(x' - x'(U')).$$

We deduce that $f^1 = f^2 + f^3$ where

$$(18) \qquad f^2 = \frac{1}{2} \sum_{1 \le l \le d' - 1} (f_l - f_{l+d'}) \cdot b_l \cdot (v'_l - v'_l(U')) \\ \qquad \qquad + \frac{1}{2} \sum_{0 \le l \le d' - 2} (f_l - f_{l+d'}) \cdot c_{l+1} \cdot (v'_{l+1} - v'_{l+1}(U'))$$

and

$$(19) \qquad f^3 = \frac{1}{2} ((f_0 - f_{d'})b_* + (f_{d'-1} - f_{d-1})c_*) (x' - x'(U')).$$

The function $f^2$ belongs to $\mathcal{L}(\langle t' \rangle - U')$ and according to Equation (18) one can compute its coordinates in the basis $v'$ at the expense of $Qd$ additions and multiplications by a constant in $K$. As explained in Section 5 we can deduce the coordinates of $f^2$ in the basis $u'$ at the expense of $Qd$ more such operations in $K$. When $d' \ne 1$ we evaluate $f^0$ and $f^2$ recursively as functions on $E'$. In the special case when $d' = 1$ these two functions are given constants. Using Equation (19) we evaluate $f^3$ at the expense of 3 additions and $d' + 2$ multiplications provided we have precomputed the values of $x' - x'(U')$ at $b' + \langle t' \rangle$. We deduce the values of $f^1 = f^2 + f^3$ at the expense of $d'$ additions. We deduce the values of $f = f^0 + \theta^{-1}f^1$ at the expense of $d$ multiplications and $d$ additions. As for the
constants appearing in Equations (18) and (19) we deduce from the examination of Taylor expansions
that

$$
b_{*} = c_{*} = 1,
$$

$$
(20) \qquad b_l = -\Gamma(O, T, lt) - a_1/2 = -\theta(lt) \quad \text{for } l \in [1, d' - 1],
$$

$$
c_l = \Gamma(O, T, lt) + a_1/2 = \theta(lt) \quad \text{for } l \in [1, d' - 1].
$$

function BUTTERFLY_EVALUATE ($\vec{f}$)

if $d = 1$ then return $\vec{f}$; end if

▷ Recursion end

$$
\vec{f}^{\text{low}} \leftarrow \vec{f}_{[0, d'-1]}; \vec{f}^{\text{high}} \leftarrow \vec{f}_{[d', d-1]}
$$

▷ Split $\vec{f}$ in half

$$
\vec{f}^{+} \leftarrow \left(\frac{1}{2}\right) \cdot \left(\vec{f}^{\text{low}} + \vec{f}^{\text{high}}\right); \vec{f}^{-} \leftarrow \left(\frac{1}{2}\right) \cdot \left(\vec{f}^{\text{low}} - \vec{f}^{\text{high}}\right)
$$

▷ Symmetrization

$$
\vec{\alpha}^{+} \leftarrow \text{BUTTERFLY_EVALUATE}(\vec{f}^{+})
$$

▷ Recursive call

$$
r \leftarrow b_0 f_0^{-} + c_0 f_{d'-1}^{-}; t \leftarrow 0
$$

▷ Constants

if $d > 2$ then

$$
t \leftarrow \mathfrak{m}'_0 f_0^{-} + (\mathfrak{n}'_{d'-1} - \mathfrak{m}'_{d'-1}) f_{d'-1}^{-} - \langle \vec{\mathfrak{n}}', \vec{f}^{-} \rangle
$$

end if

$$
\vec{f}^{-} \leftarrow \text{BIDIAGONAL_EVALUATE}(\vec{\mathfrak{b}}, \vec{\mathfrak{c}}, \vec{f}^{-}); f_{0}^{-} \leftarrow t
$$

▷ Bi-diagonal formulas

$$
\vec{f}^{-} \leftarrow \text{VTOU_BASECHANGE}(\vec{\alpha}', \vec{f}^{-})
$$

▷ Base change from $(v'_l)$ to $(u'_l)$

$$
\vec{\alpha}^{-} \leftarrow \text{BUTTERFLY_EVALUATE}(\vec{f}^{-})
$$

▷ Recursive call

$$
\vec{\alpha}^{-} \leftarrow (\vec{\alpha}^{-} + r \cdot \vec{\xi}') \diamond \vec{t}^{\text{inv}}
$$

▷ θ-multiplication

$$
\textbf{return} \ (\vec{\alpha}^{+} + \vec{\alpha}^{-}) \oplus (\vec{\alpha}^{+} - \vec{\alpha}^{-})
$$

end function

Description

• The constants $\vec{\alpha}', \vec{b}$ and $\vec{c}$ are defined by $b_0 = b_*, c_0 = c_*$ and by Eq. (8) and (20).

• The constant $\vec{\xi}'$ denotes $(x'(b' + lt') - x'(U'))_{l \in \mathbb{Z}/d'\mathbb{Z}}$.

• The constant $\vec{v}'$ denotes $(v'_l(U'))_{l \in \mathbb{Z}/d'\mathbb{Z}}$.

• The constant $\vec{t}$ denotes $(\theta(b + l t))_{l=0, \dots, d'-1}$.

• The constant $\vec{m}'$ denotes $\vec{b} \diamond \vec{v}'$ and $\vec{n}'$ denotes $\vec{m}' + \sigma^{-1}(\vec{c} \diamond \vec{v}')$.

Parameters

<caption>FIGURE 2. Butterfly Evaluation</caption>

**Proposition 3** (Fast elliptic evaluation). There exists a constant $Q$ such that the following is true. Let $\mathbf{K}$ be a field with odd characteristic. Let $E$ be an elliptic curve over $\mathbf{K}$. Let $\delta \ge 1$ be an integer. Let $t$ be a point of order $d = 2^\delta$ in $E(\mathbf{K})$. Let $b$ be a point in $E(\mathbf{K})$ such that $\mathrm{db} \ne 0$. There exists a straight-line program that on input a function $f$ in $\mathcal{L}(\langle t \rangle)$ given by its coordinates in either basis $u$ or $v$ defined in Section 5, computes the values $f(b + lt)$ for $l \in [0, d-1]$ at the expense of $Qd \log(d)$ operations in $\mathbf{K}$.

*Proof.* This results from the discussion above. For clarity, the evaluation procedure is presented as concise pseudo-code in Figure 2.
8. INTERPOLATION

In the context of the beginning of Section 7 and using a similar recursion we explain how to recover the coefficients of a function $f \in \mathcal{L}(\langle t \rangle)$ in the basis $u$ from its values at $b + \langle t \rangle$. We assume that we are given $d$ scalars $(\alpha_l)_{0 \le l \le d-1}$ in $\mathbf{K}$. We want to compute the coordinates $(f_l)_{l \in \mathbf{Z}/d\mathbf{Z}}$ in the basis $u$ of the unique function $f$ in $\mathcal{L}(\langle t \rangle)$ such that

$$f(b + lt) = \alpha_l \quad \text{for every } l \in [0, d-1].$$

We use a recursion again. We first compute

$$\alpha_l^+ = \frac{\alpha_l + \alpha_{l+d'}}{2} \quad \text{for } 0 \le l \le d' - 1$$

the values of $f^+$ at $b' + \langle t' \rangle$. If $d' \neq 1$ we deduce, by recursion, the coordinates $(f_l^+)^{0 \le l \le d'-1}$ of $f^+$ in the basis $u'$. In the special case when $d' = 1$ we simply have $f_0^+ = \alpha_0^+$. From Equations (12) we deduce

$$\begin{aligned} f^+ &= \sum_{0 \le l \le d' - 1} f_l^+ \cdot u'_l \\ &= \frac{1}{2} \sum_{0 \le l \le d' - 1} (f_l + f_{l+d'}) u'_l. \end{aligned}$$

So for $0 \le l \le d' - 1$

$$f_l^+ = \frac{f_l + f_{l+d'}}{2}.$$

This gives half the information we need to quickly compute the $(f_l)_{l \in \mathbf{Z}/d\mathbf{Z}}$. We will be done when we compute the

$$\frac{f_l - f_{l+d'}}{2} \quad \text{for } 0 \le l \le d' - 1.$$

According to Equation (13) these $d'$ scalars are the coordinates of $f^1 = \theta \cdot f^-$ in the basis of $\mathcal{L}(\langle t' \rangle + O' - U')$ made of the $\theta \cdot (u_l - u_{l+d'})$ for $l \in [0, d' - 1]$. We will use an auxiliary function $\xi_b$ on $E'$ having degree $d' + 1$, polar divisor $\langle t' \rangle + O'$, and vanishing at $b' + lt'$ for every $l \in [0, d' - 1]$. Such a function is unique up to a multiplicative constant, and it has no pole and no zero at $U'$ because $db \neq 0$. We therefore assume that $\xi_b(U') = 1$. There exists $d' + 1$ scalars $\mathfrak{d}_0, \mathfrak{d}_1, \dots, \mathfrak{d}_{d'-1}, \mathfrak{d}_*$ such that

$$(21) \qquad \xi_b = \mathfrak{d}_* \cdot x' + \sum_{0 \le l \le d' - 1} \mathfrak{d}_l \cdot v_l' = 1 + \mathfrak{d}_* \cdot (x' - x'(U')) + \sum_{1 \le l \le d' - 1} \mathfrak{d}_l \cdot (v_l' - v_l'(U')).$$

We assume that we have precomputed these $d' + 1$ scalars. We now compute

$$\alpha_l^- = \frac{\alpha_l - \alpha_{l+d'}}{2} \quad \text{and} \quad \alpha_l^1 = \theta(b + lt) \cdot \alpha_l^- \quad \text{for } 0 \le l \le d' - 1$$

at the expense of $Qd$ operations in $\mathbf{K}$, assuming we have precomputed the values of $\theta$ at $b + \langle t \rangle$. The $\alpha_l^1$ are the values of $f^1$ at $b' + \langle t' \rangle$. We let $f^*$ be the function in $\mathcal{L}(\langle t' \rangle)$ that takes value $\alpha_l^1$ at $b' + lt$. We obtain by recursion the coordinates of $f^*$ in the basis $u'$. We base change to $v'$ and obtain $d'$ scalars $f_l^*$ such that

$$f^* = \sum_{0 \le l \le d' - 1} f_l^* v_l'$$
We are not quite done. The function $f^*$ is not the one we are looking for because it does not vanish at $U'$. So we compute

$$f^*(U') = \sum_{0 \le l \le d' - 1} f_l^* \cdot v_l'(U')$$

at the expense of $d'$ multiplications and $d' - 1$ additions, provided we have precomputed the $v_l'(U')$. And we write

$$f^* = f^*(U') + \sum_{1 \le l \le d' - 1} f_l^* \cdot (v_l' - v_l'(U'))$$

We deduce

$$ \begin{aligned} f^1 &= f^* - f^*(U') \cdot \xi_b \\ &= -f^*(U') \cdot \mathfrak{d}_* \cdot (x' - x'(U')) + \sum_{1 \le l \le d' - 1} (f_l^* - \mathfrak{d}_l \cdot f^*(U')) \cdot (v_l' - v_l'(U')). \end{aligned} $$

This gives the coordinates

$(s_*, s_1, s_2, \dots, s_{d'-1})$

of $f^1$ in the basis of $\mathcal{L}(\langle t' \rangle + O' - U')$ made of $x' - x'(U')$ and the $v_l' - v_l'(U')$ for $l \in [1, d' - 1]$. Indeed

$$s_l = f_l^* - \mathfrak{d}_l \cdot f^*(U') \quad \text{for } 1 \le l \le d' - 1, \text{ and } s_* = -f^*(U') \cdot \mathfrak{d}_*.$$

We want to deduce the coordinates

$(t_0, t_1, t_2, \dots, t_{d-1})$

of $f^1$ in the basis of $\mathcal{L}(\langle t' \rangle + O' - U')$ made of the $\theta \cdot (u_l - u_{l+d'})$ for $l \in [0, d' - 1]$. According to Equations (14), (15), and (16) this amounts to solving the following linear system.

$$ \begin{pmatrix} \mathfrak{b}_* & 0 & \cdots & 0 & \mathfrak{c}_* \\ \mathfrak{c}_1 & \mathfrak{b}_1 & 0 & \cdots & 0 \\ 0 & \mathfrak{c}_2 & \mathfrak{b}_2 & \cdots & \vdots \\ \vdots & \cdots & \cdots & \cdots & 0 \\ 0 & \cdots & 0 & \mathfrak{c}_{d'-1} & \mathfrak{b}_{d'-1} \end{pmatrix} \begin{pmatrix} t_0 \\ t_1 \\ \vdots \\ t_{d-2} \\ t_{d'-1} \end{pmatrix} = \begin{pmatrix} s_* \\ s_1 \\ \vdots \\ s_{d'-2} \\ s_{d'-1} \end{pmatrix} $$

This system has a very special form. One says it is cyclic bidiagonal. Solving it requires $Qd$ operations in $K$ as recalled in section 9. In the special case when $d' = 1$ we simply invoke Equation (17).

**Proposition 4 (Fast elliptic interpolation).** There exists a constant $Q$ such that the following is true. Let $K$ be a field with odd characteristic. Let $E$ be an elliptic curve over $K$. Let $\delta \ge 1$ be an integer. Let $t$ be a point of order $d = 2^\delta$ in $E(K)$. Let $b$ be a point in $E(K)$ such that $db \ne 0$. There exists a straight-line program that on input $d$ scalars $(\alpha_l)_{0 \le l \le d-1}$, computes a function $f$ in $\mathcal{L}(\langle t \rangle)$, given by its coordinates in either basis $u$ or $v$ defined in Section 5, such that $f(b+lt) = \alpha_l$ for every $l \in [0, d-1]$, at the expense of $Qd \log(d)$ operations in $K$.

*Proof.* This results from the discussion above. To aid understanding, the interpolation procedure is shown in compact pseudo-code in Figure 3.

□
**function BUTTERFLY_INTERPOLATE(α̃)**

if $d = 1$ then return $\vec{\alpha}$; end if

$\vec{\alpha}^{\text{low}} \leftarrow \vec{\alpha}_{[0, d'-1]}$; $\vec{\alpha}^{\text{high}} \leftarrow \vec{\alpha}_{[d', d-1]}$

$\vec{\alpha}^{+} \leftarrow (1/2) \cdot (\vec{\alpha}^{\text{low}} + \vec{\alpha}^{\text{high}})$; $\vec{\alpha}^{-} \leftarrow (1/2) \cdot (\vec{\alpha}^{\text{low}} - \vec{\alpha}^{\text{high}})$

$\vec{f}^{+} \leftarrow \text{BUTTERFLY_INTERPOLATE}(\vec{\alpha}^{+})$

$\vec{f}^{-} \leftarrow \text{BUTTERFLY_INTERPOLATE}(\vec{t} \diamond \vec{\alpha}^{-})$

$\vec{f}^{-} \leftarrow \text{UTOV_BASECHANGE}(\vec{\alpha}', \vec{f}^{-})$

$f_{*} \leftarrow \langle \vec{v}', \vec{f}^{-} \rangle$; $\vec{f}^{-} \leftarrow \vec{f}^{-} - f_{*} \cdot \vec{v}$; $f_{0}^{-} \leftarrow -\vec{v}_{0} f_{*}$

$\vec{f}^{-} \leftarrow \text{BIDIAGONAL_INVERT}(\vec{b}, \vec{c}, \vec{f}^{-})$

return $(\vec{f}^{+} + \vec{f}^{-}) \oplus (\vec{f}^{+} - \vec{f}^{-})$

end function

Description

▷ Recursion end

▷ Split $\vec{\alpha}$ in half

▷ Symmetrization

▷ Recursive call

▷ Recursive call

▷ Base change from $(u'_{l})$ to $(v'_{l})$

▷ $\zeta_b$-normalization

▷ Bi-diagonal inversion

* The constants $\vec{\alpha}', \vec{b}$ and $\vec{c}$ are defined by $b_0 = b_*, c_0 = c_*$ and by Eq. (8) and (20).

* The constant $\vec{v}'$ is defined by Eq. (21).

* The constant $\vec{v}'$ denotes $(v'_{l}(U'))_{l \in \mathbb{Z}/d'\mathbb{Z}}$.

* The constant $\vec{t}$ denotes $(\theta(b+l t))_{l=0, \dots, d'-1}$.

Parameters

<caption>FIGURE 3. Butterfly Interpolation</caption>

## 9. CYCLIC BIDIAGONAL LINEAR SYSTEMS

Let $d \ge 2$ be an integer. Let $\mathbf{R}$ be a commutative local ring. Let $b_0, b_1, \dots, b_{d-1}$ and $c_0, c_1, \dots, c_{d-1}$ be scalars in $\mathbf{R}$. Let $\mathbf{M}$ be the matrix

$$ M = \begin{pmatrix} b_0 & 0 & \cdots & 0 & c_0 \\ c_1 & b_1 & 0 & \cdots & 0 \\ 0 & c_2 & b_2 & \cdots & \vdots \\ \vdots & \ddots & \ddots & \ddots & 0 \\ 0 & \cdots & 0 & c_{d-1} & b_{d-1} \end{pmatrix}. $$

The determinant of $\mathbf{M}$ is

$$ \det(\mathbf{M}) = \prod_{i=0}^{d-1} b_i - (-1)^d \prod_{i=0}^{d-1} c_i. $$

We assume that it is invertible in $\mathbf{R}$. The ring $\mathbf{R}$ being local, this implies that either every $b_i$ is invertible or every $c_i$ is invertible. Let $s_0, s_1, \dots, s_{d-1}$ be scalars in $\mathbf{R}$. We solve the system

$$ \begin{pmatrix} b_0 & 0 & \cdots & 0 & c_0 \\ c_1 & b_1 & 0 & \cdots & 0 \\ 0 & c_2 & b_2 & \cdots & \vdots \\ \vdots & \ddots & \ddots & \ddots & 0 \\ 0 & \cdots & 0 & c_{d-1} & b_{d-1} \end{pmatrix} \begin{pmatrix} t_0 \\ t_1 \\ \vdots \\ t_{d-2} \\ t_{d-1} \end{pmatrix} = \begin{pmatrix} s_0 \\ s_1 \\ \vdots \\ s_{d-2} \\ s_{d-1} \end{pmatrix} $$
in the $d$ unknowns $t_0, t_1, \dots, t_{d-1}$. We eliminate $t_0$ from the first equation using the second equation. We then eliminate $t_1$ from the resulting equation using the third equation. And so on. We find

$$ \det(M) \cdot t_{d-1} = - \sum_{0 \le k \le d-1} (-1)^k \left( \prod_{0 \le l \le k-1} b_l \right) \cdot s_k \cdot \left( \prod_{k+1 \le l \le d-1} c_l \right). $$

Assume that all the $b_i$ are invertible. We compute $t_0 = (s_0 - c_0t_{d-1})/b_0$, $t_1 = (s_1 - c_1t_0)/b_1$, ..., $t_{d-2} = (s_{d-2} - c_{d-2}t_{d-3})/b_{d-2}$. Otherwise we know that all the $c_i$ are invertible. We then compute $t_{d-2} = (s_{d-1} - b_{d-1}t_{d-1})/c_{d-1}$, $t_{d-3} = (s_{d-2} - b_{d-2}t_{d-2})/c_{d-2}$, ..., $t_1 = (s_2 - b_2t_2)/c_2$, $t_0 = (s_1 - b_1t_1)/c_1$.

## 10. REDUCTION

Under the hypotheses at the beginning of Section 5 and assuming that the base field $K$ has odd characteristic, we let $b$ be a point on $E$ defined over a separable closure $K_s$ of $K$. We assume that $db \neq 0$. We call $B$ the coset of $b$ under the action of $\langle t \rangle$. So

$$ B = b + \langle t \rangle = \{b, b+t, b+2t, \dots, b+(d-1)t\}. $$

We assume that $B$ is left invariant by the absolute Galois group of $K$. For every $l$ in $Z/dZ$ we denote by $x_l$ the function defined by

$$ (22) \qquad x_l = x_{lt} = x \circ \tau_{-lt}. $$

We assume that we are given scalars $(F_l)_{l \in Z/dZ}$ in $K$ and we consider the function

$$ F = \sum_{l \in Z/dZ} F_l x_l. $$

We are interested in such functions because they form a $t$-invariant supplementary subspace to $\mathcal{L}(\langle t \rangle)$ in $\mathcal{L}(2\langle t \rangle)$. Given such a function $F$ there exists a unique function

$$ f = \sum_{l \in Z/dZ} f_l u_l $$

such that $f$ and $F$ agree on $B$ that is

$$ f(b+lt) = F(b+lt) \quad \text{for every } l \in Z/dZ. $$

We write this condition $f \equiv F \pmod{B}$. The purpose of this section is to compute the $(f_l)_{l \in Z/dZ}$ once given the $(F_l)_{l \in Z/dZ}$. We shall achieve this at the expense of $Qd \log(d)$ additions and scalar multiplications thanks to a recursion as in the previous sections. We define $f^+$ and $f^-$ as in Equation (10) and similarly

$$ F^+ = \frac{F + F \circ \tau_T}{2} \quad \text{and} \quad F^- = \frac{F - F \circ \tau_T}{2}. $$

We set

$$ B' = \varphi(B) = b' + \langle t' \rangle = \{b', b' + t', b' + 2t', \dots, b' + (d' - 1)t'\} $$
and notice that $F^+ \equiv f^+$ (mod $B$) and $F^- \equiv f^-$ (mod $B$). From Equation (9) we deduce

$$ (23) \qquad \begin{aligned} F^+ &= \sum_{l \in \mathbb{Z}/d\mathbb{Z}} \frac{F_l + F_{l+d'}}{2} x_l \\ &= \sum_{0 \le l \le d' - 1} \frac{F_l + F_{l+d'}}{2} (x_l + x_{l+d'}) \\ &= \sum_{0 \le l \le d' - 1} \frac{F_l + F_{l+d'}}{2} (x'_l + x(T)) \\ &= F_*^+ + \sum_{0 \le l \le d' - 1} F_l^+ \cdot x'_l \end{aligned} $$

where

$$ F_*^+ = (x(T)/2) \sum_{0 \le l \le d-1} F_l \quad \text{and} \quad F_l^+ = \frac{F_l + F_{l+d'}}{2} \quad \text{for } 0 \le l \le d' - 1. $$

When $d' \neq 1$, we obtain by recursion $d'$ scalars $(f_l^+)$$_{0 \le l \le d' - 1}$ such that

$$ \sum_{0 \le l \le d' - 1} F_l^+ \cdot x'_l \equiv \sum_{0 \le l \le d' - 1} f_l^+ \cdot u'_l \mod B'. $$

Using Equation (11) we deduce

$$ \begin{aligned} F^+ &= F_*^+ + \sum_{0 \le l \le d' - 1} f_l^+ \cdot u'_l \mod B \\ &= F_*^+ + \sum_{0 \le l \le d' - 1} f_l^+ \cdot (u_l + u_{l+d'}) \\ &= F_*^+ \sum_{0 \le l \le d - 1} u_l + \sum_{0 \le l \le d' - 1} f_l^+ \cdot (u_l + u_{l+d'}) \\ &= \sum_{0 \le l \le d' - 1} (f_l^+ + F_*^+)(u_l + u_{l+d'}). \end{aligned} $$

So we have reduced $F^+$ in that case. In the special case when $d' = 1$ we notice that

$$ (24) \qquad x'_0 = x' \equiv x'(b') \pmod{B'} $$

and we deduce from Equation (23) that

$$ F^+ \equiv (F_*^+ + F_0^+ \cdot x'(b')) \cdot (u_0 + u_1) \pmod{B}. $$
We now proceed to reducing $F^-$. We define

$$ 
\begin{align} 
F^1 &= \theta \cdot F^- \notag \\ 
&= \sum_{l \in \mathbb{Z}/d\mathbb{Z}} \frac{F_l - F_{l+d'}}{2} \cdot \theta \cdot x_l \tag{25} \\ 
&= \sum_{0 \le l \le d' - 1} \frac{F_l - F_{l+d'}}{2} \cdot \theta \cdot (x_l - x_{l+d'}). \notag 
\end{align} 
$$

For every integer $l$ in $[1, d' - 1]$ the product $\theta \cdot (x_l - x_{l+d'})$ is a function on $E'$. It belongs to the linear space $\mathcal{L}(-2[lt'] - [O'] + [U'])$. So there exist constants $\mathfrak{e}_l$ and $\mathfrak{f}_l$ in $\mathbf{K}$ such that

$$ (26) \qquad \theta \cdot (x_l - x_{l+d'}) = \mathfrak{f}_l \cdot (x'_l - x'_l(U')) + \mathfrak{e}_l \cdot (v'_l - v'_l(U')). $$

In the special case $l = 0$ the product $\theta \cdot (x_0 - x_{d'})$ belongs to the linear space $\mathcal{L}(-3[O'] + [U'])$. So there exist constants $\mathfrak{f}_0$ and $\mathfrak{g}$, such that

$$ (27) \qquad \theta \cdot (x_0 - x_{d'}) = \mathfrak{f}_0 \cdot (x' - x'(U')) + \mathfrak{g} \cdot (y' - y'(U')). $$

Using Equations (25), (26), and (27) we write

$$ 
\begin{align} 
F^1 &= F^2 + F^3 + F^4 \quad \text{with} \notag \\ 
F^2 &= \frac{1}{2} \sum_{1 \le l \le d' - 1} (F_l - F_{l+d'}) \cdot \mathfrak{e}_l \cdot v'_l \notag \\ 
&\quad - \frac{1}{2} \sum_{1 \le l \le d' - 1} (F_l - F_{l+d'}) \cdot (\mathfrak{f}_l \cdot x'_l(U') + \mathfrak{e}_l \cdot v'_l(U')) \tag{28} \\ 
&\quad - \frac{1}{2} (F_0 - F_{d'}) \cdot (\mathfrak{f}_0 \cdot x'_0(U') + \mathfrak{g} \cdot y'(U')), \notag \\ 
F^3 &= \frac{1}{2} \sum_{0 \le l \le d' - 1} (F_l - F_{l+d'}) \cdot \mathfrak{f}_l \cdot x'_l, \notag \\ 
F^4 &= \frac{1}{2} (F_0 - F_{d'}) \cdot \mathfrak{g} \cdot y'. \notag 
\end{align} 
$$

The above expression for $F^2$ provides its coordinates $(F_l^2)_{0 \le l \le d' - 1}$ in the basis $v'$ of $\mathcal{L}(\langle t' \rangle)$. There exist scalars $(F_l^3)_{0 \le l \le d' - 1}$ such that

$$ F^3 \equiv \sum_{0 \le l \le d' - 1} F_l^3 \cdot v'_l \quad (\text{mod } B'). $$

To obtain these scalars we use the expression for $F^3$ in Equation (28). If $d' \ne 1$ we apply reduction recursively, then base change from $u'$ to $v'$. In the special case when $d' = 1$ we deduce from Equation (24) that

$$ F^3 \equiv \frac{1}{2}(F_0 - F_1) \cdot \mathfrak{f}_0 \cdot x'(b') \cdot v'_0 \quad (\text{mod } B'). $$
As for $F^4$, we assume that we have precomputed scalars $(\mathfrak{h}_l)_{0\le l\le d'-1}$ such that

$$ (29) \qquad y' \equiv \sum_{0 \le l \le d'-1} \mathfrak{h}_l \cdot v'_l \pmod{B'}. $$

We compute

$$ F_l^4 = (1/2)(F_0 - F_{d'}) \cdot \mathfrak{g} \cdot \mathfrak{h}_l \quad \text{for } l \in [0, d' - 1] $$

and we have

$$ F^4 \equiv \sum_{0 \le l \le d' - 1} F_l^4 \cdot v'_l \pmod{B'}. $$

Finally we compute

$$ f_l^* = F_l^2 + F_l^3 + F_l^4. $$

The reduction of $F^1$ modulo $B'$ is

$$ f^* = \sum_{0 \le l \le d' - 1} f_l^* \cdot v_l' \equiv F^1 \pmod{B'}. $$

Since $F^1 = \theta \cdot F^-$ we would now divide $f^*$ by $\theta$. But $\theta$ vanishes at $U$ and $f^*$ does not. Again we overcome this difficulty thanks to the function $\xi_b$ introduced before Equation (21). We compute

$$ f^*(U') = \sum_{0 \le l \le d' - 1} f_l^* \cdot v_l'(U') $$

and check that the function

$$ f^1 = f^* - f^*(U') \cdot \xi_b \in \mathcal{L}(\langle t' \rangle + O' - U') $$

is congruent to $F^1$ modulo $B'$. So

$$ \begin{align} f^{-} &= f^{1} \cdot \theta^{-1} \in \mathcal{L}(\langle t \rangle) \nonumber \\ (30) \qquad \begin{aligned} &= f^{*} \cdot \theta^{-1} - f^{*}(U') \cdot \xi_{b} \cdot \theta^{-1} \\ &= f^{*}(U') \cdot (1 - \xi_{b}) \cdot \theta^{-1} + \sum_{1 \le l \le d' - 1} f_{l}^{*} \cdot (v_{l}' - v_{l}'(U')) \cdot \theta^{-1} \in \mathcal{L}(\langle t \rangle). \end{aligned} \end{align} $$

Indeed this function is congruent to $F^-$ and it belongs to $\mathcal{L}(\langle t \rangle)$. We are not quite done because $f^-$ is given as a linear combination of $(1 - \xi_b) \cdot \theta^{-1}$ and the $(v'_l - v'_l(U')) \cdot \theta^{-1}$. Fortunately, for every $l$ in $[1, d' - 1]$ there exist constants $i_l$ and $j_l$ such that

$$ (31) \qquad (v'_l - v'_l(U')) \cdot \theta^{-1} = i_l \cdot (v_l - v_{l+d'}) + j_l. $$

And there exist constants $(\mathfrak{l})_{0 \le l \le d-1}$ such that

$$ (32) \qquad (1 - \xi_b) \cdot \theta^{-1} = \sum_{0 \le l \le d-1} \mathfrak{l}_l \cdot v_l. $$

Substituting Equations (31) and (32) in (30) we terminate the reduction of $F^-$. We finally compute $f = f^+ + f^-$.
To complete this discussion we give simple expressions for the constants involved in Equations (26), (27), (31). Thus $g = 1$, $f_0 = a_1/2$, and for $l \in [1, d' - 1]$

$$ 
\begin{align} 
\mathfrak{e}_l &= x(lt) - x(lt + T), \notag \\ 
f_l &= \theta(lt), \tag{33} \\ 
i_l &= \theta^{-1}(lt), \notag \\ 
j_l &= 1. \notag 
\end{align} 
$$

**Proposition 5 (Fast elliptic reduction).** There exists a constant $Q$ such that the following is true. Let $K$ be a field with odd characteristic. Let $E$ be an elliptic curve over $K$. Let $\delta \ge 1$ be an integer. Let $t$ be a point of order $d = 2^\delta$ in $E(K)$. Let $b$ be a point on $E$ defined over a separable extension of $K$. Assume $db \ne 0$. Let

$$ B = b + \langle t \rangle = \{b, b+t, b+2t, \dots, b+(d-1)t\}. $$

Assume $B$ is left invariant by the absolute Galois group of $K$. So $B$ is a reduced $K$-subscheme of $E$ having dimension 0 and degree $d$. Let $(x_l)_{0 \le l \le d-1}$ be the functions on $E$ defined in Equation 22. There exists a straight-line program that on input $d$ scalars $(F_l)_{0 \le l \le d-1}$, computes $d$ scalars $(f_l)_{0 \le l \le d-1}$, such that the functions

$$ F = \sum_{l \in \mathbb{Z}/d\mathbb{Z}} F_l x_l \quad \text{and} \quad f = \sum_{l \in \mathbb{Z}/d\mathbb{Z}} f_l u_l $$

are congruent modulo $B$, meaning

$$ f(b+lt) = F(b+lt) \quad \text{for every integer } l \in [0, d-1]. $$

This straight-line program consists of no more than $Qd \log(d)$ operations in $K$.

*Proof.* This results from the discussion above. To make the procedure clearer, we present the reduction routine as brief pseudo-code in Figure 4.

<figure><img src="image_5.png" /></figure>

A consequence of Propositions 5 and 3 is that, in case when $b$ is a $K$-rational point, we can evaluate $F = \sum_{l \in \mathbb{Z}/d\mathbb{Z}} F_l x_l$ at all $b + lt$ by first computing the reduction $f$ of $F$ modulo $B$, then evaluating $f$ at the $b + lt$. This takes time $Qd \log(d)$.

## 11. MULTIPLICATION IN THE RESIDUE RING OF A FIBER

Under the hypotheses and notation of Section 10 we explain how to multiply in the residue ring $L$ at $B$. We assume that there exists a $K$-rational point $R$ on $E$ such that $dR \ne 0$. The residue ring $L$ at $B$ is an algebra of dimension $d$ over $K$. We let $(u_l)_{l \in \mathbb{Z}/d\mathbb{Z}}$ be the functions on $E$ defined in Section 5. For every $l$ in $\mathbb{Z}/d\mathbb{Z}$ we let

$$ (34) \qquad \theta_l = u_l \mod B $$

be the image of $u_l$ in $L$. The $\theta_l$ form a basis $\Theta$ of $L$ over $K$. We want to multiply two elements

$$ f = \sum_{l \in \mathbb{Z}/d\mathbb{Z}} f_l \theta_l \quad \text{and} \quad g = \sum_{l \in \mathbb{Z}/d\mathbb{Z}} g_l \theta_l $$
**function BUTTERFLY_REDUCE($\vec{F}$)**

if $d = 1$ then return $x'(b') \cdot \vec{F}$; end if

$\vec{F}^{\text{low}} \leftarrow \vec{F}_{[0, d'-1]}$; $\vec{F}^{\text{high}} \leftarrow \vec{F}_{[d', d-1]}$

$\vec{F}^{+} \leftarrow (1/2) \cdot (\vec{F}^{\text{low}} + \vec{F}^{\text{high}})$; $\vec{F}^{-} \leftarrow (1/2) \cdot (\vec{F}^{\text{low}} - \vec{F}^{\text{high}})$

$F_{*}^{+} \leftarrow x(T) \sum_{l=0}^{d'-1} F_{l}^{+}$

$\vec{f}^{+} \leftarrow F_{*}^{+} + \text{BUTTERFLY_REDUCE}(\vec{F}^{+})$

$\vec{f}^{1} \leftarrow \text{BUTTERFLY_REDUCE}(\vec{f} \diamond \vec{F}^{-})$

$\vec{f}^{1} \leftarrow \text{UTOV_BASECHANGE}(\vec{a}', \vec{f}^{1})$

$\vec{f}^{1} \leftarrow \vec{f}^{1} + \vec{e} \diamond \vec{F}^{-} + F_{0}^{-} \cdot \vec{h}$

$f_{0}^{1} \leftarrow f_{0}^{1} + \mathfrak{p}_{*} F_{0}^{-} - \langle \vec{p}, \vec{F}^{-} \rangle$

$f_{*}^{-} \leftarrow \langle \vec{v}', \vec{f}^{1} \rangle$

$\vec{f}^{-} \leftarrow -\vec{f}^{1} \diamond \vec{i} - f_{*}^{-} \cdot \vec{l}$; $f_{0}^{-} \leftarrow l_{*} f_{*}^{-}$

$\vec{f}^{-} \leftarrow \text{VTOU_BASECHANGE}(\vec{a}^{1}, \vec{f}^{-})$

$\vec{f}^{-} \leftarrow f_{*}^{-} \cdot (l_{0} - l_{*}) + \sum_{l=0}^{d'} f_{l}^{1} + \vec{f}^{-}$

**return** $(\vec{f}^{+} + \vec{f}^{-}) \oplus (\vec{f}^{+} - \vec{f}^{-})$

**end function**

▷ Recursion end

▷ Split $\vec{F}$ in half

▷ Symmetrization

▷ Star constant

▷ Recursive call and shift by $F_{*}^{+}$

▷ Recursive call

▷ Base change from $(u'_{l})$ to $(v'_{l})$

▷ Normalization

▷ Recover $f^{-}$

▷ Back base change

▷ Shift vector entries

Description

• The constants $\vec{a}'$ is defined by Equation (8) and $\vec{a}^1$ denotes $\vec{a}_{[d', d-1]} - \vec{a}_{[0, d'-1]}$.

• The constants $\vec{e}, \vec{f}$ and $\vec{i}$ are defined by Equation (33).

• The constant $\vec{h}$ is defined by Equation (29).

• The constant $\vec{l}$ is defined by Equation (32), restricted to $l \in [0, d' - 1]$ and $l_{*}$ denotes the constant $l_{d'}$.

• The constant $\vec{v}'$ denotes $(v'_{l}(U'))_{l \in \mathbb{Z}/d'\mathbb{Z}}$.

• The constant $\vec{p}$ denotes the vector $\vec{e} \diamond \vec{v}' + \vec{f} \diamond (x'(U' + l t'))_{l \in \mathbb{Z}/d'\mathbb{Z}}$, and $\mathfrak{p}_{*}$ denotes the constant $\mathfrak{p}_{0} - \mathfrak{f}_{0} x'(U') - y'(U')$.

Parameters

<caption>FIGURE 4. Butterfly Reduction</caption>

in $\mathbf{L}$, given by their coordinates in the basis $\Theta$. We follow [CL09, Section 4.3.4.]. We first define two functions

$$ \mathcal{F} = \sum_{l \in \mathbb{Z}/d\mathbb{Z}} f_l u_l \quad \text{and} \quad \mathcal{G} = \sum_{l \in \mathbb{Z}/d\mathbb{Z}} g_l u_l $$

in $\mathcal{L}(\langle t \rangle)$ such that $f = \mathcal{F} \mod B$ and $g = \mathcal{G} \mod B$. As a consequence of Equation (5) we can decompose the product $\mathcal{FG}$ as a sum $\mathcal{C} + \mathcal{D}$ where

$$ \mathcal{C} = \sum_{l \in \mathbb{Z}/d\mathbb{Z}} (f_l - f_{l-1})(g_l - g_{l-1}) x_l \quad \text{and} \quad \mathcal{D} = \mathcal{F}\mathcal{G} - \mathcal{C} \in \mathcal{L}(\langle t \rangle). $$

Using the method presented in Section 10 we find scalars $C_l$ such that

$$ \mathcal{C} \equiv \sum_{l \in \mathbb{Z}/d\mathbb{Z}} C_l u_l \pmod{B}. $$
It remains to compute the coordinates of $\mathcal{D}$ in the basis $(u_l)_{l \in \mathbb{Z}/d\mathbb{Z}}$. To this end we evaluate $\mathcal{F}$ and $\mathcal{G}$ at the points $R+lt$ for $l \in [0, d-1]$ as explained in Section 7. We also evaluate $\mathcal{C}$ at the points $R+lt$. This is achieved in two steps. We first reduce $\mathcal{C}$ modulo $R+\langle t \rangle$ using the method from Section 10. We then evaluate the resulting function at the $R+lt$ using the method in Section 7. We can thus compute

$$ \mathcal{D}(R+lt) = \mathcal{F}(R+lt) \cdot \mathcal{G}(R+lt) - \mathcal{C}(R+lt) \quad \text{for } l \in [0, d-1]. $$

The interpolation method from Section 8 finds scalars $D_l$ such that

$$ \mathcal{D} = \sum_{l \in \mathbb{Z}/d\mathbb{Z}} D_l u_l. $$

Finally the product of $f$ and $g$ is

$$ fg = \sum_{l \in \mathbb{Z}/d\mathbb{Z}} (C_l + D_l)\theta_l. $$

Every step in this calculation is achieved in time $Qd \log(d)$.

This finishes the proof of Proposition 6 below. One can summarize this proof by saying that the equi-variant version of Chudnovsky's multiplication algorithm summarized in [CL09, Lemma 6] involves convolutions by constant vectors corresponding to evaluation, interpolation and reduction. And the main consequence of the recursive approach presented in Sections 7, 8, and 10 is that these convolution products are achieved in time $Od \log(d)$ when $d$ is a power of two.

The multiplication is depicted in streamlined pseudo-code in Figure 5.

**function** NORMALBASIS_MULTIPLY($\vec{f}, \vec{g}$)

$\vec{\alpha} \leftarrow BUTTERFLY_EVALUATE(\vec{f})$

$\vec{\beta} \leftarrow BUTTERFLY_EVALUATE(\vec{g})$

$\vec{H} \leftarrow (\vec{f} - \sigma(\vec{f})) \diamond (\vec{g} - \sigma(\vec{g}))$

$\vec{h} \leftarrow BUTTERFLY_REDUCE_R(\vec{H})$

▷ Reduction modulo $R + \langle t \rangle$

$\vec{\gamma} \leftarrow BUTTERFLY_EVALUATE(\vec{h})$

$\vec{\delta} \leftarrow \vec{\alpha} \diamond \vec{\beta} - \vec{\gamma}$

$\vec{k} \leftarrow BUTTERFLY_INTERPOLATE(\vec{\delta})$

$\vec{h} \leftarrow BUTTERFLY_REDUCE_B(\vec{H})$

▷ Reduction modulo $b + \langle t \rangle$

return $\vec{h} + \vec{k}$

end function

Description

* The function `BUTTERFLY_REDUCE_R` (resp. `BUTTERFLY_REDUCE_B`) denotes the reduction routine defined in Figure 4, where the constants $h$ and $l$ depend on the point $R \in E(\mathbb{K})$ (resp. $b \in E(\mathbb{K}_s)$) as defined in Section 13.

Parameters

<caption>FIGURE 5. Butterfly Multiplication</caption>
**Proposition 6 (Multiplication).** There exists a constant $Q$ such that the following is true. Let $\mathbf{K}$ be a field with odd characteristic. Let $E$ be an elliptic curve over $\mathbf{K}$. Let $\delta \ge 1$ be an integer. Let $t$ be a point of order $d = 2^\delta$ in $E(\mathbf{K})$. We assume that there exists a point $R$ in $E(\mathbf{K})$ such that $dR \ne 0$. Let $b$ be a point on $E$ defined over a separable extension of $\mathbf{K}$. Assume $db \ne 0$. Let

$$B = b + \langle t \rangle = \{b, b+t, b+2t, \dots, b+(d-1)t\}.$$

Assume $B$ is left invariant by the absolute Galois group of $\mathbf{K}$. Let $\mathbf{L}$ be the residue ring at $B$. Let $\Theta = (\theta_l)_{l \in \mathbb{Z}/l\mathbb{Z}}$ be the basis of $\mathbf{L}$ defined in Equation 34. There exists a straight-line program that on input $2d$ scalars $(f_l)_{0 \le l \le d-1}$ and $(g_l)_{0 \le l \le d-1}$ in $\mathbf{K}$ computes $d$ scalars $(h_l)_{0 \le l \le d-1}$, such that

$$\sum_{l \in \mathbb{Z}/d\mathbb{Z}} h_l \theta_l = \left( \sum_{l \in \mathbb{Z}/d\mathbb{Z}} f_l \theta_l \right) \times \left( \sum_{l \in \mathbb{Z}/d\mathbb{Z}} g_l \theta_l \right) \in \mathbf{L}.$$

This straight-line program consists of no more than $Qd \log(d)$ operations in $\mathbf{K}$.

We denote by $\otimes$ the multiplication law on $\mathbf{K}^d$ defined in Proposition 6 above. We call it an elliptic multiplication law.

**Remark 7.** In the special case where $b = R$, the pseudo-code in Figure 5 simplifies significantly. It consists of two evaluations,

$$\vec{\alpha} \leftarrow BUTTERFLY \_EVALUATE(\vec{f}) \quad \text{and} \quad \vec{\beta} \leftarrow BUTTERFLY \_EVALUATE(\vec{g}),$$

followed by a single interpolation, $BUTTERFLY \_INTERPOLATE(\vec{\alpha} \diamond \vec{\beta})$, which yields the result.

## 12. IMPLEMENTATION AND EXPERIMENTAL RESULTS

A public implementation of the algorithms presented in this work, developed in the computational algebra system Magma [BCP97], is available online [CL25]. In particular, it includes the classical Fast Fourier Transform (FFT) with the Cooley–Tukey algorithm, as well as our routines, both requiring only $O(d \log(d))$ operations in $\mathbb{F}_q$ for evaluations or interpolations at $d$ points.

We present in Table 1 the timings measured for our Magma routines on a standard laptop, working modulo a 64-bit prime. These results are provided for comparison purposes only, given the likely overhead of the Magma interpreter in these highly recursive routines. The relative differences in timings, namely, a constant factor of about 5 between Cooley-Tukey and the elliptic butterflies, are nonetheless significant in practice.

## 13. ELLIPTIC BASES FOR FINITE FIELD EXTENSIONS

We have shown in [CL09, Theorem 2] how evaluation and interpolation of elliptic functions provide normal bases for finite field extensions having quasi-linear time multiplication. We consider in this section extensions of degree $d = 2^\delta$, a power of two. We construct normal bases allowing multiplication in time $Qd \log(d)$ for these extensions. We prove the theorem below.

**Theorem 8 (Fast normal bases).** There exists an absolute constant $Q$ such that the following is true. Let $\mathbf{K}$ be a finite field with $q$ elements and odd characteristic. Let $\delta \ge 1$ be an integer and let $d = 2^\delta$. Assume that $4d^4 \le q$. Let $\mathbf{L}/\mathbf{K}$ be a field extension of degree $d$. There exists a normal $\mathbf{K}$-basis $\Theta$ of $\mathbf{L}$ and a straight-line program that takes as input the coordinates in $\Theta$ of two elements in $\mathbf{L}$ and returns the coordinates of their product, at the expense of $Qd \log(d)$ operations in $\mathbf{K}$.
<table><caption>TABLE 1. Fast evaluation and interpolation timings</caption><thead><tr><th rowspan="2">log<sub>2</sub>(d)</th><th colspan="2">Cooley-Tukey [Wik]</th><th colspan="2">Elliptic Butterflies</th></tr><tr><th>Evaluate</th><th>Interpolate</th><th>Evaluate (Fig. 2)</th><th>Interpolate (Fig. 3)</th></tr></thead><tbody><tr><td>8</td><td>0.00 s</td><td>0.01 s</td><td>0.00 s</td><td>0.01 s</td></tr><tr><td>9</td><td>0.01 s</td><td>0.01 s</td><td>0.02 s</td><td>0.01 s</td></tr><tr><td>10</td><td>0.02 s</td><td>0.03 s</td><td>0.03 s</td><td>0.03 s</td></tr><tr><td>11</td><td>0.04 s</td><td>0.03 s</td><td>0.06 s</td><td>0.07 s</td></tr><tr><td>12</td><td>0.03 s</td><td>0.03 s</td><td>0.13 s</td><td>0.15 s</td></tr><tr><td>13</td><td>0.05 s</td><td>0.06 s</td><td>0.28 s</td><td>0.31 s</td></tr><tr><td>14</td><td>0.12 s</td><td>0.12 s</td><td>0.56 s</td><td>0.66 s</td></tr><tr><td>15</td><td>0.25 s</td><td>0.25 s</td><td>1.17 s</td><td>1.38 s</td></tr><tr><td>16</td><td>0.54 s</td><td>0.55 s</td><td>2.45 s</td><td>2.72 s</td></tr></tbody></table>

*Proof.* Let $K_s$ be a separable closure of $K$. Let $\nu$ be the 2-valuation of $q - 1$.

We first consider the case $\nu \ge \delta + 1$. This is a classical and favorable case because we have enough roots of unity. Let $b$ be an element of order $2^{\delta+\nu}$ in $K_s^*$. Let

$$t = b^q/b = b^{q-1}.$$

This is an element of multiplicative order $d$ in $K^*$. So $K(b)$ is a degree $d$ extension of $K$. We call it $\mathbf{L}$. We set

$$a = b^d.$$

This is an element of order $2^\nu$ in $K^*$. So $\mathbf{L}$ is a Kummer extension of $K$ and $(1, b, b^2, \dots, b^{d-1})$ is a basis for $\mathbf{L}$ over $K$. We call it $\Pi$. Multiplication of two elements of $\mathbf{L}$ given by their coordinates in the basis $\Pi$ requires $Qd \log(d)$ operations in $K$ using FFT because we have an element of order $2d$ in $K^*$. Of course $\Pi$ is not a normal basis. But the element

$$\theta = \sum_{l=0}^{d-1} b^k$$

is a normal element in $\mathbf{L}/K$. We denote by $\Theta$ the associated normal basis and notice that passing from $\Theta$ to $\Pi$ is a Fourier transform of order $d$. This requires no more than $Qd \log(d)$ operations in $K$ because we have $d$-th roots of unity in $K^*$ and $d$ is a power of two.

We now consider the case $\nu \le \delta$ and use the constructions in [CL09]. The integer $d_q$ as defined in [CL09, Definition 1] is either $d^2$ or $2d^2$. Since $4d^4 \le q$ we have $d_q \le \sqrt{q}$ and we can apply [CL09, Lemma 9]. There exists an elliptic curve $E$ over $K$, a point $t$ of order $d$ in $E(K)$, a point $R$ in $E(K)$ such that $dR \ne 0$, and a point $b$ in $E(K_s)$ such that $db \ne 0$ and the conjugate of $b$ by the Frobenius of $E/K$ is $b+t$. We set $B = b + \langle t \rangle$. This is a degree $d$ divisor on $E$. It is irreducible over $K$. The residue field at $B$ is a degree $d$ extension field of $K$ which we denote by $\mathbf{L}$. The $\theta_l$ introduced in Equation 34 form a normal basis $\Theta$ of $\mathbf{L}$. According to Proposition 6 there exists a straight-line program that multiplies to elements in $\mathbf{L}$ given by their coordinates in the basis $\Theta$ at the expense of $Qd \log(d)$ operations in $K$.

□
## 14. ELLIPTIC GOPPA CODES

In this section we study a family of Goppa codes in genus one allowing particularly fast encoding. We let $p$ be an odd prime and $q$ a power of $p$. Let $d = 2^\delta$ with $\delta \ge 1$. We set $d' = d/2$ and assume that

$$ (35) \qquad q \ge \max\left(\frac{d^4}{4}, (2d+1)^2+1\right). $$

The length of the Hasse interval is $4\sqrt{q}$. So there are two consecutive multiples of $d^2$ in it. At least one of them is not congruent to 1 modulo $p$. We deduce that there exists an elliptic curve $E$ over $\mathbb{K}$ such that $E(\mathbb{K})$ contains a point $t$ of order $d$. We set

$$ T = d't. $$

Since $\#E(\mathbb{K}) \ge (\sqrt{q}-1)^2 > 4d^2$, there exists a point $Q$ in $E(\mathbb{K})$ such that

$$ 2dQ \neq 0. $$

We denote by

$$ \iota : E \to E $$

the involution $P \mapsto -P$ and by

$$ x : E \to E/\langle\iota\rangle $$

the quotient map. We call

$$ \mathcal{L} = \mathcal{L}(\langle t \rangle)^\iota $$

the subspace of $\mathcal{L}(\langle t \rangle)$ fixed by $\iota$. We can see $\mathcal{L}$ as the linear space on the genus zero curve $E/\langle\iota\rangle$ associated to the divisor

$$ \sum_{l \in [1, d'-1]} [x(lt)]. $$

In particular functions in $\mathcal{L}$ have no pole at $\{O, T\}$. And $\mathcal{L}$ has dimension $d'$. The functions

$$ \ell_0 = 1 \quad \text{and} \quad \ell_l = v_l - v_{-l} \quad \text{for } l \in [1, d' - 1] $$

form a basis $\ell$ of $\mathcal{L}$. We let

$$ \text{ev} : \mathcal{L} \to \mathbb{K}^d $$

be the evaluation map at the $Q + lt$ for $0 \le l \le d-1$. Every function $f$ in $\mathcal{L}$ has degree $\le d-2$ as a function on $E$ and degree $\le d' - 1$ as a function on $E/\langle\iota\rangle$. If $f$ vanishes at $d'$ points in $Q + \langle t \rangle$ then it has $d'$ zeros on $E/\langle\iota\rangle$ also because $\iota(Q + \langle t \rangle)$ and $Q + \langle t \rangle$ do not intersect. So $f$ must be zero. We deduce that $\text{ev}$ is injective and its image $C \subset \mathbb{K}^d$ is a linear code of length $d$, dimension $d'$ and minimum distance $d' + 1$. We can see $C$ as a subcode of a genus one Goppa code or as a genus zero Goppa code.

**Encoding** amounts to evaluating the map $\text{ev}$ at a function $f$ in $\mathcal{L}$ given by its coordinates $(m_l)_{0 \le l \le d' - 1}$ in the basis $\ell$. To this end we first compute the coordinates $(n_l)_{l \in \mathbb{Z}/\mathbb{Z}}$ of $f$ in the basis $v$ of $\mathcal{L}(\langle t \rangle)$. We have

$$ n_0 = m_0, n_{d'} = 0, \text{ and } n_l = m_l, \text{ and } n_{-l} = -m_l \text{ for } 1 \le l \le d' - 1. $$

We deduce the coordinates of $f$ in the basis $u$ using Equations (6) and (7). We finally evaluate $f$ at $Q + \langle t \rangle$ as explained in Section 7. So encoding takes time $Qd \log(d)$.

**Checking** a received word is achieved in time $Qd \log(d)$ also by first interpolating the received values as explained in Section 8. We obtain the coordinates $(n_l)_{l \in \mathbb{Z}}$ in the basis $v$ of a function $f \in \mathcal{L}(\langle t \rangle)$ taking the received values. We check that $f \in \mathcal{L}$ that is

$$ n_{d'} = 0 \text{ and } n_l + n_{-l} = 0 \text{ for } 1 \le l \le d' - 1. $$
The initial message is given by the coordinates $(m_l)_{0 \le l \le d' - 1}$ of $f$ in the basis $\ell$ of $\mathcal{L}$. Namely

$$m_0 = n_0 \text{ and } m_l = n_l \text{ for } 1 \le l \le d' - 1.$$

If the check is failed we deduce that there are errors. We then remind that $C$ is a genus 0 code under its clothes of genus 1 code. In this context, **decoding** up to half the minimum distance is achieved thanks to an half-gcd computation as explained in [Sug+75]. The time complexity is $Q \log(d)$ times the complexity of multiplying two polynomials of degree $d$ and coefficients in $K$. See [GG13, Chapter 11]. The following theorem summarizes the content of this section

**Theorem 9 (Fast MDS codes).** *There exists a constant $Q$ such that the following is true. Let $d \ge 2$ be a power of two and let $q$ be an odd prime power such that inequality (35) holds true. Let $K$ be a field with $q$ elements. There exists a $[d, d/2, d/2 + 1]$ linear code $C$ over $K$, a straight-line program that encodes $C$ at the expense of $Qd \log(d)$ operations in $K$, a straight-line program that checks $C$ at the expense of $Qd \log(d)$ operations in $K$, and a computation tree that corrects $C$ up to $d/4$ errors at the expense of $Qd \log^2(d) \log(\log(d))$ operations and comparisons in $K$.*

## 15. ELLIPTIC LWE CRYPTOGRAPHY

Our third application concerns the construction of secure cryptographic schemes within quantum computational models, and in particular in the presence of quantum algorithms such as those of Shor [Sho97] and Grover [Gro96]. Among the cryptosystems most extensively studied in this context are those based on computational problems over Euclidean lattices.

The Learning With Errors (LWE) problem, introduced by Regev, constitutes a foundational hardness assumption underlying these constructions.

**Definition 10 (Short-LWE assumptions [Reg05], [Reg09, Lemma 4.4]).** Let $d, q \ge 2$ be integers, let $\chi$ be an error distribution over $\mathbb{Z}$ (typically a discrete Gaussian), and let $\vec{s} \leftarrow \chi^m$ be a secret vector. Given pairs $(A, A\vec{s} + \vec{e})$, where $A \leftarrow (\mathbb{Z}/q\mathbb{Z})^{m \times d}$ is uniformly random and $\vec{e} \leftarrow \chi^m$,

* the short-search-LWE problem is to recover $\vec{s}$,

* the short-decisional-LWE problem is to distinguish, with non-negligible advantage, such pairs from uniformly random pairs $(A, \vec{v})$.

Regev establishes that solving the decisional Learning With Errors (LWE) problem on average is at least as hard as solving worst-case approximation problems on Euclidean lattices using a quantum algorithm. This remarkable reduction renders the LWE problem particularly attractive for cryptographic applications.

We wonder if replacing matrix products with the elliptic multiplication law defined in Proposition 6 could lead to cryptosystems with improved practical characteristics. This question motivates the definition of an *Elliptic-LWE assumption*, which transposes LWE into the setting of elliptic multiplications.

**Definition 11 (Short-Elliptic-LWE assumptions).** Let $E$ be an elliptic curve defined modulo a prime $q$, with a point $t \in E(\mathbb{Z}/q\mathbb{Z})$ of order $d = 2^\delta$, and let $b$ be another point in $E(\mathbb{Z}/q\mathbb{Z})$ such that $db \ne 0$. Let $\otimes$ denote the multiplication law modulo $b + \langle t \rangle$ defined in Proposition 6. Finally, let $\chi$ be an error distribution over $\mathbb{Z}$ (typically a discrete Gaussian), and let $\vec{s} \leftarrow \chi^d$ be a secret vector.

Given pairs $(\vec{a}, \vec{v} = \vec{a} \otimes \vec{s} + \vec{e})$, where $\vec{a} \in (\mathbb{Z}/q\mathbb{Z})^d$ is uniformly random and $\vec{e} \leftarrow \chi^d$:

* the short-search Elliptic-LWE problem is to recover $\vec{s}$,

* the short-decisional Elliptic-LWE problem is to distinguish such pairs from uniformly random pairs $(\vec{a}, \vec{v})$ with non-negligible advantage.

Transposing the Regev construction in this setting yields a CPA-secure encryption scheme.
**Theorem 12.** *Under the Short-Elliptic-LWE assumption, the encryption scheme defined in Figure 6 is CPA-secure.*

The correctness follows from the equation $p = \langle \vec{r}, \vec{e} \rangle - \langle \vec{e}_1, \vec{s} \rangle + e_2 + \mu \lfloor q/2 \rfloor$. Since the error term remains small for appropriate parameter choices, decryption succeeds with high probability. Security in the chosen-plaintext attack (CPA) model follows directly from the decisional LWE assumption. The vector $\vec{w}$ is indistinguishable from a random vector, which makes $\vec{c}_1$ and $\langle \vec{r}, \vec{w} \rangle + e_2$ two LWE instances, each indistinguishable from uniform, thereby perfectly hiding the message.

In terms of complexity, applying $\phi_{\vec{a}}$ requires $Q d \log(d)$ operations in $\mathbb{Z}/q\mathbb{Z}$, using the algorithms described in Section 10. By Tellegen's principle, applying the transpose map $^t\phi_{\vec{a}}$ has the same asymptotic complexity [BLS03; KKB88]. More generally, all steps can be carried out within this time bound. The sizes of the keys and ciphertexts are $O(d)$ bits.

Let $E$ be an elliptic curve defined modulo a prime $q$ with a point $t \in E(\mathbb{Z}/q\mathbb{Z})$ of order $d = 2^\delta$, and let $b$ be another rational point in $E(\mathbb{Z}/q\mathbb{Z})$ such that $db \neq 0$. Let $\otimes$ denote the multiplication law modulo $b + \langle t \rangle$ defined by Proposition 6. Let $\phi_{\vec{a}}$ (and $^t\phi_{\vec{a}}$) denote the linear endomorphism $\vec{x} \mapsto \vec{a} \otimes \vec{x}$ in the coordinate system defined by the basis $(u_l)_{l \in \mathbb{Z}/d\mathbb{Z}}$ (and its transposed for the canonical scalar product $\langle \cdot, \cdot \rangle$ in this basis).

Parameters

Sample a uniformly random $d$-dimensional vector $\vec{a}$, a secret key $\vec{s} \leftarrow \chi^d$ and an error vector $\vec{e} \leftarrow \chi^d$. Compute $\vec{w} = \phi_{\vec{a}}(\vec{s}) + \vec{e}$. The public key is $(\vec{a}, \vec{w})$.

Key Generation

To encrypt a message $\mu \in \{0, 1\}$, sample random vectors $\vec{r} \leftarrow \chi^d$, $\vec{e}_1 \leftarrow \chi^d$ and $e_2 \leftarrow \chi$. Then compute $\vec{c}_1 = {}^t\phi_{\vec{a}}(\vec{r}) + \vec{e}_1$ and $c_2 = \langle \vec{r}, \vec{w} \rangle + e_2 + \lfloor q/2 \rfloor \mu$. The ciphertext is the pair $(\vec{c}_1, c_2)$.

Encryption

To decrypt a ciphertext $(\vec{c}_1, c_2)$, compute $p = c_2 - \langle \vec{c}_1, \vec{s} \rangle$. The plaintext bit $\mu$ is then recovered by comparing $p$ to $\lfloor q/2 \rfloor$: if $p$ is closer to 0 than to $\lfloor q/2 \rfloor$ modulo $q$, output 0; otherwise, output 1.

Decryption

<caption>**FIGURE 6.** Elliptic-LWE encryption scheme variant of Regev's construction</caption>

We note that the elliptic-decisional LWE assumption bears similarities to the decisional Ring-LWE assumption studied by Lyubashevsky, Peikert, and Regev [LPR13], in which matrix products are replaced by polynomial multiplications in the quotient ring $\mathbb{Z}_q[x]/\langle\Phi(x)\rangle$ with $\Phi(x)$ a cyclotomic polynomial. When the degree $d$ is a power of two, this construction operates in the cyclotomic field defined by $\Phi_{2d}(x) = x^d + 1$. Since polynomial multiplication modulo such a polynomial can be implemented using Fast-Fourier-Transform algorithms provided that $2d$ divides $q-1$, this approach achieves quasi-linear time complexity in $d$ for encryption and decryption operations, while also significantly reducing key sizes. This substantial improvement forms the core of cryptosystems currently undergoing standardization [Ava+24; Bai+24].
Ideally, one would like in our case to encrypt similarly *d* bits at once, as is possible with Ring-LWE. However, the elliptic multiplication law is not compatible with the norms on the operands, and it appears difficult to control the error in the computations; as a result, the decryption algorithm simply would not work. The reason for this incompatibility of the elliptic multiplication law with the norms on the operands is that the coefficients of the multiplication tensor ⊗ have no reason to be small as is the case in the cyclotomic context.

On the other hand, the difficulty of conveniently lifting the Elliptic-LWE problem to characteristic zero may be seen as an argument for the Short-Elliptic-LWE assumption that the Ring-LWE assumption lacks. Also it would be desirable to design more efficient schemes relying on the Short-Elliptic-LWE assumption.

Nevertheless, it is still possible to make the algorithm in Figure 6 practical using standard techniques from LWE-based cryptography (see, for instance, the Frodo construction [Bos+16]). More precisely, let $\beta$ and $\ell$ be two integer parameters. We can encrypt a $\beta\ell^2$-bit message, denoted $(\mu_{i,j})$, by applying the following modifications to the algorithm in Figure 6.

* Set the modulus $q$ so that $\beta$ bits of message can be recovered, instead of a single one, in the decryption equation;

* Generate secret keys consisting of $\ell$ vectors $\vec{s}_0, \dots, \vec{s}_{\ell-1}$, associated with public keys $(\vec{a}, \vec{w}_0, \dots, \vec{w}_{\ell-1})$;

* Replace the ciphertext vector $\vec{c}_1$ by $\ell$ vectors computed using $\ell$ ephemeral secret vectors $\vec{r}_0, \dots, \vec{r}_{\ell-1}$;

* Compute $\ell \times \ell$ elements $c_2$, one for each inner product $\langle \vec{r}_i, \vec{w}_j \rangle_{i,j}$ and addition of $e_{i,j} + \lfloor q/2^{\lfloor 1+\log_2 q \rfloor - \beta} \rfloor \mu_{i,j}$.

The overall time complexity is now $O(\ell d \log d + \ell^2 d)$ operations in $\mathbb{Z}/q\mathbb{Z}$. The first term arises from the application of the endomorphisms $\phi_{\vec{a}}$ and ${}^t\phi_{\vec{a}}$ to the $\ell$ secret vectors $\vec{s}_i$ and $\vec{r}_j$ during key generation and encryption. The second term corresponds to the $\ell^2$ inner products in the encryption or decryption. The size of the exchanged data is $\ell d \log q$ bits.

Using the same parameters as those employed in standards, for instance $\beta = 4$ and $\ell = 8$ to encrypt 256-bit data [Bos+16], we observe that this results in timings and data sizes that are not fundamentally larger than in the single-bit version.

More generally, most lattice-based cryptographic schemes roughly follow the following asymptotic guideline, expressed in terms of a security parameter $\lambda$ tending to infinity: $d = O(\lambda \log \lambda)$, $q = O(d)$, and $\chi$ a discrete gaussian distribution with standard deviation $\sigma = O(\sqrt{d})$ [Pei16, Chapter 4]. We can then extract $\beta = O(\log d)$ bits per inner product. Setting furthermore $\ell \approx O(\sqrt{d}/\log d)$ so that we can encrypt on the order of $\lambda$ bits. The total time complexity becomes $O(d^2/\log d)$ bit-operations, while the data sizes become $O(d^{3/2})$ bits. These complexities are significantly better than the $O(d^{5/2})$ bit-operations of [Bos+16] for same data size.

REFERENCES

[Ava+24] R. Avanzi et al. CRYSTALS-Kyber Algorithm Specifications and Supporting Documentation. Tech. rep. NIST PQC Standard. National Institute of Standards and Technology, 2024.

[Bai+24] S. Bai et al. CRYSTALS-Dilithium Algorithm Specifications and Supporting Documentation. Tech. rep. NIST PQC Standard. National Institute of Standards and Technology, 2024.
[BR04] S. Ballet and R. Rolland. "Multiplication algorithm in a finite field and tensor rank of the multiplication". In: *J. Algebra* 272.1 (2004), pp. 173–185.

[Bal99] S. Ballet. "Curves with many points and multiplication complexity in any extension of $F_q$". In: *Finite Fields Appl.* 5.4 (1999), pp. 364–377.

[Ben+23] E. Ben-Sasson, D. Carmon, S. Kopparty, and D. Levit. "Elliptic curve fast Fourier transform (ECFFT). I: Low-degree extension in time $O(n \log n)$ over all finite fields". English. In: *Proceedings of the 34th annual ACM-SIAM symposium on discrete algorithms, SODA 2023*. New York, NY: Association for Computing Machinery (ACM), 2023, pp. 700–737.

[Bos+16] J. Bos, C. Costello, L. Ducas, I. Mironov, M. Naehrig, V. Nikolaenko, A. Raghunathan, and D. Stebila. "Frodo: Take off the Ring! Practical, Quantum-Secure Key Exchange from LWE". In: *Proceedings of the 2016 ACM SIGSAC Conference on Computer and Communications Security*. Vienna, Austria: Association for Computing Machinery, 2016, pp. 1006–1018. Available at https://eprint.iacr.org/2016/659.

[BCP97] W. Bosma, J. Cannon, and C. Playoust. "The Magma algebra system. I. The user language". In: *J. Symbolic Comput.* 24.3-4 (1997). Computational algebra and number theory (London, 1993), pp. 235–265. Available at https://doi.org/10.1006/jsco.1996.0125.

[BLS03] A. Bostan, G. Lecerf, and É. Schost. "Tellegen's principle into practice". In: *Proceedings of the 2003 International Symposium on Symbolic and Algebraic Computation*. Philadelphia, PA, USA: Association for Computing Machinery, 2003, pp. 37–44.

[BS05] A. Bostan and É. Schost. "Polynomial evaluation and interpolation on special sets of points". English. In: *J. Complexity* 21.4 (2005), pp. 420–446.

[BCS97] P. Bürgisser, M. Clausen, and M. A. Shokrollahi. *Algebraic complexity theory*. Vol. 315. With the collaboration of Thomas Lickteig. Springer-Verlag, Berlin, 1997, pp. xxiv+618.

[Cha08] J. Chaumine. "Multiplication in small finite fields using elliptic curves". In: *Algebraic geometry and its applications*. Vol. 5. World Sci. Publ., Hackensack, NJ, 2008, pp. 343–350.

[CC88] D. V. Chudnovsky and G. V. Chudnovsky. "Algebraic complexities and algebraic curves over finite fields". English. In: *J. Complexity* 4.4 (1988), pp. 285–316.

[CC89] D. V. Chudnovsky and G. V. Chudnovsky. *Computational problems in arithmetic of linear differential equations. Some diophantine applications*. English. Number theory, Semin., New York/NY 1985-88, Lect. Notes Math. 1383, 12-49 (1989). 1989.

[CG23] J.-M. Couveignes and J. Gasnier. *Explicit Riemann-Roch spaces in the Hilbert class field*. 2023. arXiv: 2309.06754 [math.NT]. Available at https://arxiv.org/abs/2309.06754.

[CL09] J.-M. Couveignes and R. Lercier. "Elliptic periods for finite fields". In: *Finite Fields Their Appl.* 15.1 (2009), pp. 1–22.

[CL25] J.-M. Couveignes and R. Lercier. *Elliptic Butterflies*. Version 1.0.0. Computer Software. Oct. 2025. Available at https://github.com/rlercier/Elliptic-Butterflies.

[GS95] A. García and H. Stichtenoth. "A tower of Artin-Schreier extensions of function fields attaining the Drinfeld-Vladut bound". In: *Invent. Math.* 121.1 (1995), pp. 211–222.

[GG13] J. von zur Gathen and J. Gerhard. *Modern Computer Algebra* (3. ed.) Cambridge University Press, 2013.

[Gop81] V. D. Goppa. "Codes on algebraic curves". In: *Dokl. Akad. Nauk SSSR* 259.6 (1981), pp. 1289–1290.

[Gop82] V. D. Goppa. "Algebraic-geometric codes". In: *Izv. Akad. Nauk SSSR Ser. Mat.* 46.4 (1982), pp. 762–781, 896.
[Gro96] L. K. Grover. "A fast quantum mechanical algorithm for database search". In: *Proceedings of the Twenty-Eighth Annual ACM Symposium on Theory of Computing*. Philadelphia, Pennsylvania, USA: Association for Computing Machinery, 1996, pp. 212–219.

[HJB85] M. T. Heideman, D. H. Johnson, and C. S. Burrus. "Gauss and the history of the fast Fourier transform". English. In: *Arch. Hist. Exact Sci.* 34 (1985), pp. 265–277.

[Iha81] Y. Ihara. "Some remarks on the number of rational points of algebraic curves over finite fields". In: *J. Fac. Sci. Univ. Tokyo Sect. IA Math.* 28.3 (1981), 721–724 (1982).

[KKB88] M. Kaminski, D. G. Kirkpatrick, and N. H. Bshouty. "Addition requirements for matrix and transposed matrix products". In: *J. Algorithms* 9.3 (Sept. 1988), pp. 354–364.

[LPR13] V. Lyubashevsky, C. Peikert, and O. Regev. "On Ideal Lattices and Learning with Errors over Rings". In: *Journal of the ACM* 60.6 (2013), 43:1–43:35.

[Pei16] C. Peikert. "A Decade of Lattice Cryptography". In: *Foundations and Trends*® in Theoretical Computer Science 10.4 (2016), pp. 283–424.

[Ran12] H. Randriambololona. "Bilinear complexity of algebras and the Chudnovsky-Chudnovsky interpolation method". In: *J. Complexity* 28.4 (2012), pp. 489–517.

[Reg05] O. Regev. "On Lattices, Learning with Errors, Random Linear Codes, and Cryptography". In: *Proceedings of the 37th Annual ACM Symposium on Theory of Computing*. 2005, pp. 84–93.

[Reg09] O. Regev. "On lattices, learning with errors, random linear codes, and cryptography". In: *Proceedings of the 41st annual ACM symposium on Theory of computing (STOC)*. ACM. 2009, pp. 84–93.

[Sha24] C. J. Sharp. *Pacific double-saddle butterflyfish (Chaetodon ulietensis) and other Chaetodon, Moorea — Wikimedia Commons*. 2024. Available at https://commons.wikimedia.org/wiki/File:Pacific_double-saddle_butterflyfish_(Chaetodon_ulietensis)_and_other_Chaetodon_Moorea.jpg.

[Sho92] M. A. Shokrollahi. "Optimal algorithms for multiplication in certain finite fields using elliptic curves". In: *SIAM J. Comput.* 21.6 (1992), pp. 1193–1198.

[Sho97] P. W. Shor. "Polynomial-Time Algorithms for Prime Factorization and Discrete Logarithms on a Quantum Computer". In: *SIAM Journal on Computing* 26.5 (1997), pp. 1484–1509.

[STV92] I. E. Shparlinski, M. A. Tsfasman, and S. G. Vladut. "Curves with many points and multiplication in finite fields". In: *Coding theory and algebraic geometry (Luminy, 1991)*. Vol. 1518. Springer, Berlin, 1992, pp. 145–169.

[SL10] A. Soro and J. Lacan. "FNT-Based Reed-Solomon Erasure Codes". In: *7th IEEE Consumer Communications and Networking Conference, CCNC 2010, Las Vegas, NV, USA, January 9-12, 2010*. IEEE, 2010, pp. 1–5.

[Sug+75] Y. Sugiyama, M. Kasahara, S. Hirasawa, and T. Namekawa. "A Method for Solving Key Equation for Decoding Goppa Codes". In: *Inf. Control.* 27.1 (1975), pp. 87–99. Available at https://doi.org/10.1016/S0019-9958(75)90090-X.

[TVZ82] M. A. Tsfasman, S. G. Vlăduț, and T. Zink. "Modular curves, Shimura curves, and Goppa codes, better than Varshamov-Gilbert bound". In: *Math. Nachr.* 109 (1982), pp. 21–28.

[Vél71] J. Vélu. "Isogénies entre courbes elliptiques". In: *C. R. Acad. Sci. Paris Sér. A-B* 273 (1971), A238–A241.

[Vél78] J. Vélu. "Courbes elliptiques munies d'un sous-groupe $\mathbb{Z}/n\mathbb{Z} \times \mu_n$". French. In: *Bull. Soc. Math. Fr., Suppl., Mém.* 57 (1978), p. 152. Available at https://eudml.org/doc/94779.

[Wik] Wikipedia. *Cooley–Tukey FFT algorithm*. Online; accessed 2025. Available at https://en.wikipedia.org/wiki/Cooley–Tukey_FFT_algorithm.
JEAN-MARC COUVEIGNES, UNIV. BORDEAUX, CNRS, INRIA, BORDEAUX-INP, IMB, UMR 5251, F-33400 TALENCE, FRANCE.

Email address: jean-marc.couveignes@u-bordeaux.fr

REYNALD LERCIER, DGA & UNIV. RENNES, CNRS, IRMAR - UMR 6625, F-35000 RENNES, FRANCE.

Email address: reynald.lercier@univ-rennes.fr