1: Elliptic curve crytography¶
In this notebook, we give an introduction to elliptic curve cryptography (ECC). Elliptic curves are everywhere nowadays, e.g.
- Crypto: Many cryptocurrencies, like Bitcoin and Ethereum, utilize ECC as part of their cryptographic algorithms. ECC is used for generating digital signatures to ensure the integrity and authenticity of transactions. It offers strong security with relatively smaller key sizes, which is crucial for efficient blockchain operations. We will explain and build the ECDSA algorithm that they use.
- Email encryption: ECC is used in email encryption and digital signatures, protecting sensitive email communications between individuals and organizations. For example, it is used in in OpenPGP, the most widely used email encryption standard. We will explain and build the ECDH algorithm that they use.
- Web browser and server connection: ECC is used in Transport Layer Security (TLS) and its predecessor Secure Sockets Layer (SSL), which provide encrypted connections between web browsers and servers (both using ECDSA and ECDH).
Hence, elliptic curve cryptography constitute a fundamental concept for a secure and private digital world.
In this tutorial, we will describe the algorithms behind these use-cases:
- Mathematical introduction The mathematical background behind elliptic curves and what they do.
- Elliptic curve cryptography An elliptic curve crytopgrahic system.
- Elliptic-curve Diffie–Hellman (ECDH): a key exchange protocol based on elliptic curve cryptography. It allows two parties to generate a shared secret over an insecure channel, which can be used for secure communication or encryption, without ever transmitting the secret directly.
- Elliptic Curve Digital Signature Algorithm (ECDSA): this is a widely used digital signature algorithm. ECDSA provides a way to generate and verify digital signatures, ensuring the authenticity, integrity, and non-repudiation of messages or data.
Code: you can download the jupyter notebook that created this post here.
Acknowledgement: I took inspiration from other tutorials on the interent,in particular I want to acknowledge Andrea Corbellini's posts and Nick Sullivan's post.
A few tools that we will need later:
from typing import List, Union, Tuple, Literal
import numpy as np
import matplotlib.pyplot as plt
#!pip install bitarray
import bitarray
class Timer:
"""Function to time a piece of code"""
def __init__(self,message):
self.message = message
def __enter__(self):
self.start = time.perf_counter()
return self
def __exit__(self, *args):
self.end = time.perf_counter()
self.interval = self.end - self.start
print(f"Time to execute {self.message} in seconds = {self.interval:.2e}")
def get_binary_representation(number: int):
"""function to get binary representation of number"""
binary_representation = ""
while number:
new_string = str(number % 2)
binary_representation = new_string + binary_representation
number = number // 2
return binary_representation
1. Elliptic Curves¶
Let's start by defining what we mean with elliptic curves. Let $\mathbb{K}$ be a field). In our case, this will be
- $\mathbb{K}=\mathbb{R}$: we use the real number for illustration purposes as it allows us to develop geometric intuition. However, they are generally not used in cryptography.
- $\mathbb{K}=GF(p)=\mathbb{Z}/p\mathbb{Z}$: this field is used for cryptographic purposes for a large primer number $p$.
For two elements $a,b\in\mathbb{K}$, an elliptic curves is a subset $E_{a,b}\subset \mathbb{K}\times \mathbb{K}$ defined by the equation: $$y^2 = x^3+ax+b$$ i.e. $E_{a,b}$ consists of all pairs $(x,y)$ satisfying this equation.
Symmetry of elliptic curves: An immediate property of an elliptic curve is clear: if $(x,y)\in E_{a,b}$, then $(x,-y)\in E_{a,b}$ (as $y^2=(-y)^2$). Therefore, elliptic curves are symmetric around the $x$-axis.
Excluding singular curves: A technnical requirement for elliptic curves is that $4*a^3+27*b^2\neq 0$. The reason for that is that singularities have to be avoided. If you would like to know more about this, we are going to give a full explanation later in the appendix chapter.
1.1. Elliptic Curve Class¶
Let's define an elliptic curve class.
class EllipticCurve:
def __init__(self, a: int, b: int, mod_n=None, allow_singularity: bool = False):
if mod_n is not None: assert isinstance(mod_n,int)
assert ((isinstance(a,int) and isinstance(b,int) and mod_n is not None))|\
((isinstance(a,float) and isinstance(b,float) and mod_n is None)),\
"""K=GF(p) requires ints a,b. K=R requires floats a,b."""
self.a = a
self.b = b
self.mod_n = mod_n
if not allow_singularity:
self._check_singularity()
@property
def determinant(self):
determinant = -4*(self.a**3)-27*(self.b**2)
if self.mod_n is not None:
determinant = determinant % self.mod_n
return determinant
def _check_singularity(self):
assert self.determinant != 0, "Singular curve not allowed."
def __str__(self):
equation = f"y^2 = x^3 + {self.a}*x + {self.b}"
if self.mod_n is not None:
equation = f"{equation} mod {self.mod_n}"
return equation
def __eq__(self, op):
return (self.a == op.a) & (self.b == op.b) & (self.mod_n == op.mod_n)
def mask_points_in_curve(self,x_arr,y_arr):
"""Function that finds in input all points that allow on the curve"""
y_side = y_arr**2
x_side = x_arr**3+self.a*x_arr+self.b
if self.mod_n is not None:
y_side = y_side % self.mod_n
x_side = x_side % self.mod_n
return x_side == y_side
else:
return np.abs(x_side-y_side)<1e-5
def give_point_cloud(self, n_grid: int = 20, min_x: float = -3.0, max_x: float = 3.0):
"""
Returns cloud of points on elliptic curve
"""
if self.mod_n is None:
x_arr = np.linspace(min_x,max_x,n_grid)
y_pos = np.sqrt(x_arr**3+self.a*x_arr+self.b)
x_pos = x_arr
x_neg = x_arr[::-1]
y_neg = - y_pos[::-1]
x = np.concatenate([x_neg,x_pos])
y = np.concatenate([y_neg,y_pos])
mask_nan = np.isnan(y)
x,y = x[~mask_nan], y[~mask_nan]
else:
#Test out all possibilities for x:
x_wanted = np.arange(0,self.mod_n)
y2_wanted = (x_wanted**3+self.a*x_wanted+self.b)%self.mod_n
#Test our all possibilities for y:
y_query = np.arange(0,self.mod_n)
y2_query = (y_query**2)%self.mod_n
#Find combinations that are zero modulo self.mod_n
bool_mask = (((y2_wanted[:,None]-y2_query[None,:])%self.mod_n)==0)
#Extract corresponding x's and y's:
y_query_mat = (bool_mask.astype('int')*y_query[None,:])
x_wanted_mat = (x_wanted[:,None]*bool_mask.astype('int'))
y = y_query_mat[bool_mask]
x = x_wanted_mat[bool_mask]
#Ensure that points are on the curve:
assert self.mask_points_in_curve(x,y).all()
return x,y
def random_points(self,n_samples: int = 1):
"""Give some random points on elliptic curve"""
x,y = self.give_point_cloud()
ids = np.random.randint(len(x),size=n_samples)
return x[ids],y[ids]
def plot(self, ax, label=None,scatter=True,**kwargs):
"""plot the curve"""
if scatter:
return ax.scatter(*self.give_point_cloud(**kwargs),label=label)
else:
return ax.plot(*self.give_point_cloud(**kwargs),label=label)
Let's define some curves:
$\mathbb{K}=GF(p)$
curve_1 = EllipticCurve(10,12,14)
curve_2 = EllipticCurve(10,11,16)
print(curve_1)
print(curve_2)
y^2 = x^3 + 10*x + 12 mod 14 y^2 = x^3 + 10*x + 11 mod 16
$\mathbb{K}=\mathbb{R}$
curve_3 = EllipticCurve(1.0,2.0)
curve_4 = EllipticCurve(1.0,3.0)
print(curve_3)
print(curve_4)
try:
curve_5 = EllipticCurve(-3,2)
except AssertionError:
print("Singular curve. Not allowed...")
y^2 = x^3 + 1.0*x + 2.0 y^2 = x^3 + 1.0*x + 3.0 Singular curve. Not allowed...
1.2. Plot a few curves¶
We plot a few curves for various $a,b$ and for $\mathbb{K}=\mathbb{R},GF(p)$:
fig, axs = plt.subplots(1,3,figsize=(15,5))
for b in np.linspace(-100,100,10):
curve = EllipticCurve(.0,b)
curve.plot(ax=axs[0],min_x=-10.0,max_x=10.0,n_grid=10000,label=f"{curve.b:.0f}")
axs[0].legend(title="b")
axs[0].set_title("a=0| K = R")
for a in np.linspace(-100,100,10):
curve = EllipticCurve(a,0.0)
curve.plot(ax=axs[1],min_x=-10.0,max_x=10.0,n_grid=10000,label=f"{curve.a:.2f}")
axs[1].legend(title="a")
axs[1].set_title("b=1 | K = R")
for a in np.linspace(-100,100,10):
curve = EllipticCurve(a,-1.0)
curve.plot(ax=axs[2],min_x=-10.0,max_x=20.0,n_grid=10000,label=f"{curve.a:.2f}")
axs[2].legend(title="a")
axs[2].set_title("b=-1 | K = R")
fig, axs = plt.subplots(1,3,figsize=(15,5))
curve = EllipticCurve(1,2,23)
curve.plot(ax=axs[0],min_x=-10.0,max_x=10.0,n_grid=10000,label=f"{curve.b:.0f}")
axs[0].set_title("a=1,b=2,K=GF(23)");
curve = EllipticCurve(4,2,231)
curve.plot(ax=axs[1],min_x=-10.0,max_x=10.0,n_grid=10000,label=f"{curve.b:.0f}")
axs[1].set_title("a=4,b=2,K=GF(231)");
curve = EllipticCurve(4,12,2311)
curve.plot(ax=axs[2],min_x=-10.0,max_x=10.0,n_grid=10000,label=f"{curve.b:.0f}")
axs[2].set_title("a=4,b=12,K=GF(2311)");
/var/folders/z8/hy1b8vjj63xfzqjbz_r7pz2nx23r5d/T/ipykernel_11969/817109513.py:54: RuntimeWarning: invalid value encountered in sqrt y_pos = np.sqrt(x_arr**3+self.a*x_arr+self.b)
1.3. Find modular inverse¶
We define the algorithm that finds the modular inverse of an object. As discussed in the RSA tutorial, this is the extended Euclidean algorithm.
def xgcd(a, b):
"""Extended Euclidean Algorithm (Iterative)
Args:
a: int
b: int
NOTES
=====
We can represent gcd(a,b) = a.x + b.y
This function returns gcd(a, b), x, y
REFERENCES
==========
https://anh.cs.luc.edu/331/notes/xgcd.pdf
EXAMPLES
========
>>> xgcd(15, 35)
(5, -2, 1)
>>> xgcd(30, 20)
(10, 1, -1)
"""
assert a > 0 and b > 0
xprev, x = 0, 1
yprev, y = 1, 0
while a:
q = b // a
x, xprev = xprev - q * x, x
y, yprev = yprev - q * y, y
a, b = b % a, a
return b, xprev, yprev
def inverse(a, n):
"""Modular Multiplicative Inverse
Args:
a: integer whose inverse is to be found
n: modular base
NOTES
=====
Returns the modular multiplicative inverse of 'a' under mod n
Return value is always positive
EXAMPLES
========
>>> inverse(2, 5)
3
>>> inverse(17, 39)
23
>>> inverse(2,4)
Traceback (most recent call last):
Exception: Inverse does not exist.
"""
g, x, _ = xgcd(a, n)
if g != 1:
raise Exception("Inverse does not exist.")
return (x % n + n) % n
1.4 Elliptic Curve Algebra¶
The fascinating property of elliptic curves is that we can define an algebraic operation on elliptic curves, i.e. for two points $P,Q\in E_{a,b}$ we can define a new point $P+Q\in E_{a,b}$. The goal of this section is to motivate and to define what "$P+Q$" means.
1.4.1 Geometric intuition¶
Let's start with the geometric idea behind addition in elliptic curves. This intution only holds for $\mathbb{K}=\mathbb{R}$, i.e. for elliptic curves in our known Euclidean geometry. However, it is a good starting point.
The main idea behind elliptic curve addition is that we can add two points $P=(x_P,y_P)$ and $Q=(x_Q,y_Q)$ on a curve by performing the following operation:
- Find the line that intersects $Q$ and $P$ (in green below)
- Find the 3rd point that intersects with the elliptic curve (intersection between orange and green line)
- Mirror this point along the x-axis to reach $P+Q$.
However, there might be a few edge cases when the line $L=\bar{QP}$ does not intersect with a third point on the curve:
- Mirrored points: if $x_P=x_Q$ and $y_P=-x_Q$, then the line connecting both points will have no third intersection as it is just vertical. We define an additional point at infinity which we simply denote as $0$. We set in this case $P+Q=0$. In this case, we also write $Q=-P$.
- Point at infinity: if $P=0$, then set $P+Q=Q$. If $Q=0$, set $P+Q=P$.
- Equal points: if $P=Q$, then one takes the tangential at $P$. In the appendix, we will show that there is exactly one other point where that lines intersects with the curve. Mirror that point to get $P+P$.
- Tangentials: if 1. and 2. are not the case and the line $\bar{PQ}$ does not intersect with a third point, we will show in the appendix that the line $L=\bar{PQ}$ must be tangential either at $P$ or at $Q$. In this case, pick the point where the line is not tangential and mirror that point.
1.4.2 From geometry to algebra¶
Let's make the definition more exact and write it down in equations.
Trivial cases (edge case 1 and 2). If $P=-Q$, $P=0$, or $Q=0$, we already have simple rules on what to do by definition. This defines equations for edge case 1 and 4.
Two distinct points (including edge case 4). Let $R=P+Q=(x_R,y_R)$. One can then show that equations
$$ \begin{aligned} \lambda &= (y_Q - y_P)(x_Q - x_P)^{-1}\\ x_R &= \lambda^2 - x_P - x_Q \\ y_R &= \lambda \cdot (x_P - x_R) - y_P, \end{aligned}$$
define a point on the elliptic curve (by simply inserting the above formula in the equation). One can also see that $-R=(x_R,-y_R)$ is on the line $L(x)=\lambda (x-x_P)+y_P$ connecting $P$ and $Q$. Note that this includes edge case $4$ (Tangentials).
Equal points (edge case 3). 4. If $P = Q$ (i.e., the points are the same), then the sum $R = P + P$ is computed as follows:
$$ \begin{aligned} \lambda &= (3x_P^2 + a)(2y_P)^{-1}\\ x_R &= \lambda^2 - 2x_P \\ y_R &= \lambda \cdot (x_P - x_R) - y_P \end{aligned}$$ One can see that $-R$ lies on the tangential by computing the derivative of the function $f(x)=\sqrt{x^3+ax+b}$ at $x_P$:
$$f'(x_P)=\frac{1}{2}(3x_P^2+a)(x_P^3+ax_P+b)^{-\frac{1}{2}}=(3x_P^2+a)(2y_P)^{-1}=\lambda$$
Therefore, the tangential is given by $L(x)=\lambda (x-x_P)+y_P$ and hence, $-R$ lies on that tangential. To see that the point lies on the elliptic, simply insert the above formula into the elliptic curve equation (you won't learn a lot, so I would not do that and trust that smart people have done it correctly).
Important!: the above rules define an addition on purely algebraic terms. Therefore, we can use the above equations to also define an addition on elliptic curves from other fields such as $GF(p)$ that don't represent our classical Euclidean geometry. This is a completely way of adding points!
Let's implement the addition for an elliptic curve element:
class EllipticCurveElement:
def __init__(self, coordinates: Union[Tuple[float,float],Tuple[int,int],Literal[0]],
curve: EllipticCurve):
if coordinates:
self.x = coordinates[0]
self.y = coordinates[1]
self.is_zero = False
elif coordinates == 0:
self.is_zero = True
else:
raise ValueError("Invalid input for EllipticCurveElement. Must be either tuple or 0")
self.curve = curve
if not self.is_zero:
valid_input_for_reals = ((isinstance(self.x,int) and \
isinstance(self.y,int) and \
curve.mod_n is not None))
valid_input_for_galois = ((isinstance(self.x,float) and \
isinstance(self.y,float) and \
curve.mod_n is None))
assert valid_input_for_reals | valid_input_for_galois,\
"Invalid input. Floats for reals and ints for galois field."
assert self.curve.mask_points_in_curve(np.array([self.x]),np.array([self.y])).all(),\
"Coordinates are not part of that curve."
def __str__(self):
if self.is_zero:
return "0=inf"
else:
return f"({self.x},{self.y})"
def __add__(self, op):
"""Addition on elliptic curves.
op: must be another EllipticCurveElement
Returns a new EllipticCurveElement representing the sum."""
assert self.curve == op.curve, "Elliptic curves are different. Can't add."
#Identity elements:
if self.is_zero:
return op
elif op.is_zero:
return self
else:
#Compute lambda, the slope of the line
#connecting both points:
if self.x!= op.x:
den = (op.x-self.x)
if self.curve.mod_n is not None:
den = den % self.curve.mod_n
inverse_den = inverse(den,self.curve.mod_n)
lambda_ = ((op.y-self.y)*inverse_den) % self.curve.mod_n
else:
lambda_ = (op.y-self.y)/den
else:
#If self=-op, then return 0=inf
sum_ = self.y + op.y
if self.curve.mod_n is not None:
sum_ = sum_ % self.curve.mod_n
if sum_ == 0:
return EllipticCurveElement(0,curve=self.curve)
#Otherwise, they are equal and compute tangential:
nominator = (3*(op.x**2)+self.curve.a)
denominator = (2*op.y)
if self.curve.mod_n is not None:
nominator = nominator % self.curve.mod_n
denominator = denominator % self.curve.mod_n
inverse_den = inverse(denominator,self.curve.mod_n)
lambda_ = (nominator*inverse_den) % self.curve.mod_n
else:
lambda_ = nominator/denominator
x_R = lambda_**2 - self.x - op.x
y_R = lambda_*(op.x-x_R) - op.y
if self.curve.mod_n is not None:
x_R = x_R % self.curve.mod_n
y_R = y_R % self.curve.mod_n
return EllipticCurveElement((x_R,y_R),curve=self.curve)
def __rmul__(self, n: int):
cum_sum = EllipticCurveElement(0,self.curve)
exponent = EllipticCurveElement((self.x,self.y),self.curve)
while n>0:
if n % 2:
cum_sum = cum_sum + exponent
exponent = exponent + exponent
n = n //2
return cum_sum
def __eq__(self, el):
if not isinstance(el, EllipticCurveElement):
return False
if self.curve != el.curve:
return False
if self.is_zero:
return el.is_zero
if el.is_zero:
return False
if (el.x != self.x) | (el.y != self.y):
return False
return True
Let's look at a few additions:
curve = EllipticCurve(10,11,mod_n=31)
point = curve.random_points(1)
el_1 = EllipticCurveElement((point[0].item(),point[1].item()),curve)
print("Element 1: ", el_1)
el_2 = EllipticCurveElement(0,curve)
print("Element 2: ", el_2)
el_3 = EllipticCurveElement(0,curve)
print("Element 3: ", el_3)
print("Sum 1+2: ", el_1+el_2)
try:
el_1 + el_3
except AssertionError:
print("Expected error: elements from different curves can't be combined.")
Element 1: (24,30) Element 2: 0=inf Element 3: 0=inf Sum 1+2: (24,30)
1.4.1 Plottting summation for $\mathbb{K}=\mathbb{R}$¶
curve = EllipticCurve(1.0,10.0)
x_sol,y_sol = curve.give_point_cloud()
n_plots = 5
fig, axs = plt.subplots(1,n_plots,figsize=(6*n_plots,6))
for i in range(n_plots):
ax = axs[i]
rand_idx_1 = np.random.randint(len(x_sol))
rand_idx_2 = np.random.randint(len(x_sol))
el_1 = EllipticCurveElement((x_sol[rand_idx_1],y_sol[rand_idx_1]),curve=curve)
el_2 = EllipticCurveElement((x_sol[rand_idx_2],y_sol[rand_idx_2]),curve=curve)
el_sum = el_1 + el_2
curve.plot(ax=ax,min_x=-10.0,max_x=max(10.0,el_sum.x),n_grid=10000,label=f"{curve.b:.0f}",scatter=False)
m = (el_2.y-el_1.y)/(el_2.x-el_1.x)
x_min,x_max = ax.get_xlim()
x_arr = np.linspace(x_min,x_max,100)
ax.plot(x_arr,m*(x_arr-el_1.x)+el_1.y,color='tab:green',linestyle='--')
ax.plot([el_sum.x,el_sum.x],[el_sum.y,-el_sum.y],color='tab:orange',linestyle='--')
ax.scatter([el_1.x],[el_1.y],marker='X',color='black')
ax.scatter([el_2.x],[el_2.y],marker='X',color='black')
ax.scatter([el_sum.x],[el_sum.y],marker='X',color='red')
#ax.scatter([el_sum.x],[-el_sum.y],marker='X',color='red')
ax.scatter([el_sum.x],[el_sum.y],marker='X',color='red',s=32)
fig, axs = plt.subplots(1,n_plots,figsize=(6*n_plots,6))
curve = EllipticCurve(1.0,10.0)
x_sol,y_sol = curve.give_point_cloud()
for i in range(n_plots):
ax = axs[i]
rand_idx_1 = np.random.randint(len(x_sol))
rand_idx_2 = np.random.randint(len(x_sol))
el_1 = EllipticCurveElement((x_sol[rand_idx_1],y_sol[rand_idx_1]),curve=curve)
el_2 = el_1
el_sum = el_1 + el_2
curve.plot(ax=ax,min_x=-10.0,max_x=max(10.0,el_sum.x),n_grid=10000,label=f"{curve.b:.0f}",scatter=False)
m = (3*el_1.x**2+curve.a)/(2*el_1.y)
x_min,x_max = ax.get_xlim()
x_arr = np.linspace(x_min,x_max,100)
ax.plot(x_arr,m*(x_arr-el_1.x)+el_1.y,color='tab:green',linestyle='--')
ax.plot([el_sum.x,el_sum.x],[el_sum.y,-el_sum.y],color='tab:orange',linestyle='--')
ax.scatter([el_1.x],[el_1.y],marker='X',color='black')
ax.scatter([el_2.x],[el_2.y],marker='X',color='black')
ax.scatter([el_sum.x],[el_sum.y],marker='X',color='red')
#ax.scatter([el_sum.x],[-el_sum.y],marker='X',color='red')
#ax.scatter([el_sum.x],[el_sum.y],marker='X',color='red',s=32)
/var/folders/z8/hy1b8vjj63xfzqjbz_r7pz2nx23r5d/T/ipykernel_11969/817109513.py:54: RuntimeWarning: invalid value encountered in sqrt y_pos = np.sqrt(x_arr**3+self.a*x_arr+self.b)
