How do I convert the definition of E-521 into a curve definition a la Bouncy Castle?

Question

I am currently trying to create an ECCCurve for E-521. Unfortunately, it is not currently a named curve in the library I am using, so I will have to define it manually.

I am using the definition of the E-521 curve provided here:

E521 is a curve over GF(2^521-1), formula x^2+y^2=1-376014x^2y^2,
basepoint
(1571054894184995387535939749894317568645297350402905821437625
18115230499438118852963259119606760410077267392791511426719338990
5003276673749012051148356041324, 12), order 2^519 -
3375547632585017057891076304187826360719049612140512266186351500
85779108655765, cofactor 4.

The two curve definitions I can choose from are FpCurve and F2mCurve. I'm not sure which one to choose for this. I have looked here. Unfortunately, I don't know how to extrapolate which one to use based just on the knowledge that E-521 is an Edwards curve.

Is an Edwards curve not compatible with the base classes that are defined in the library? Do I need to wait for further development there, or can I risk attempting to implement it myself? The latter is obviously a question I will likely have to answer myself, but some extra information on the risks/other implementations I could examine would be helpful.

Answer 1

Bouncy Castle does not support curves of the Edwards form (such as the one E-521 has). See related question: Elliptic Curves of different forms

Bouncy Castle's FpCurve supports curves of the Weierstrass simple form. You can convert Edwards curves to this form by a series of (complex) calculations. SafeCurves shows you the necessary equations. First convert the Edwards curve to Montgomery form, then to Weierstrass form.

Unfortunately doing this will result in absolutely huge numbers, which - on top of the fact that Edwards curves are computed faster than Weierstrass curves - results in poor performance.

---

Appendix: Curve forms

• Weierstrass: $y^2 = x^3 + ax + b$
• Montgomery: $By^2 = x^3 + Ax^2 + x$
• Edwards: $x^2 + y^2 = 1 + dx^2y^2$

Answer 2

Edwards Curves can be converted to Montgomery or Weierstrass curves.

The only tricky part is that the point (0, 1) goes to the point at infinity. Here is some haskell that will do it. (Requires arithmoi package for extended gcd)

{-# LANGUAGE NoImplicitPrelude #-}
module E521 (p521, IntModp521, toIntModp521, fromIntModp521,
             d, a_montgomery, b_montgomery, a_weierstrass, b_weierstrass,
             bp, bp_montgomery, bp_weierstrass)
  where

import Data.Ratio (numerator, denominator)
import Math.NumberTheory.GCD (extendedGCD)
import Prelude ( Eq(..), Fractional(..), Num(..), Ord(..), Show(..),
                 Integer,
                 (^), ($), mod, otherwise)

p521 :: Integer
p521 = 2^521 - 1

newtype IntModp521 = IntModp521 Integer deriving (Show, Eq)
toIntModp521 :: Integer -> IntModp521
toIntModp521 x | x < 0     = IntModp521 $ (x `mod` p521) + p521
                   | otherwise = IntModp521 $ (x `mod` p521)
fromIntModp521 (IntModp521 x) = x

instance Num IntModp521 where
  (IntModp521 x) + (IntModp521 y) = IntModp521 $ (x + y) `mod` p521
      (IntModp521 x) * (IntModp521 y) = IntModp521 $ (x * y) `mod` p521
  abs (IntModp521 x) = IntModp521 $ x
  negate (IntModp521 x) = toIntModp521 (-x)
  fromInteger = toIntModp521
  signum (IntModp521 0) = 0
  signum _              = 1

instance Fractional IntModp521 where
  recip (IntModp521 x) = IntModp521 $ middle $ extendedGCD x p521
    where
      middle (a, b, c) = b
  fromRational a = (toIntModp521 (numerator a)) / (toIntModp521 (denominator a))

d :: IntModp521
d = -376014

bp :: (IntModp521, IntModp521)
bp = (x, y)
  where
    x = 1571054894184995387535939749894317568645297350402905821437625181152304994381188529632591196067604100772673927915114267193389905003276673749012051148356041324
    y = 12

a_montgomery :: IntModp521
a_montgomery = 2 * (1 + d) / (1 - d)

b_montgomery :: IntModp521
b_montgomery = 4 / (1 - d)

a_weierstrass :: IntModp521
a_weierstrass = (3 - a_montgomery ^ 2) / (3 * b_montgomery ^ 2)

b_weierstrass :: IntModp521
b_weierstrass = (2 * a_montgomery ^ 3 - 9 * a_montgomery) / (27 * b_montgomery ^ 3)

bp_montgomery = (u, v)
  where
    u = (1 + y) / (1 - y)
    v = u / x
    (x, y) = bp

bp_weierstrass = (u, v)
  where
    u = (x + a_montgomery / 3) / b_montgomery
    v = y / b_montgomery
    (x, y) = bp_montgomery