Sage and Elliptic Curve Point Operations (NIST P256)Elliptic curves have a base point (\(G\)) and then we perform either point additions, point doubling or point scaling operations. These are simple operations, and allow us to use fast computation of elliptic curve points. In this case we will use \(G+G\), \(2.G\) and double G operations, along with generating a random private key (\(n\)) and then computing the associated public key of \(n.G\). In this case we will use the NIST defined P256 curve (as typically used in TLS). TheoryElliptic curves have a base point (\(G\)) and then we perform either point additions, point doubling or point scaling operations. These are simple operations, and allow us to use fast computation of elliptic curve points. In this case we will use \(G+G\), \(2.G\) and double G operations, along with generating a random private key (\(n\)) and then computing the associated public key of \(n.G\). For the NIST P256, the curve parameters are: \(G\)=(48439561293906451759052585252797914202762949526041747995844080717082404635286, 36134250956749795798585127919587881956611106672985015071877198253568414405109) \(p= 2^{256}-2^{224}+2^{192}+2^{96}-1\) \(n=115792089210356248762697446949407573529996955224135760342422259061068512044369L \) \(y^2=x^3-3x+41058363725152142129326129780047268409114441015993725554835256314039467401291 \pmod p\) and where \(n\) is the order of the curve. This can be coded in Sage as: gx = 48439561293906451759052585252797914202762949526041747995844080717082404635286L gy = 36134250956749795798585127919587881956611106672985015071877198253568414405109L n = 115792089210356248762697446949407573529996955224135760342422259061068512044369L p256 = 2^256-2^224+2^192+2^96-1 a = p256 - 3 b = 41058363725152142129326129780047268409114441015993725554835256314039467401291L F = GF(p256) EC = EllipticCurve([F(a), F(b)]) EC.set_order(n) G = EC(F(gx), F(gy)) Point addition is then coded as: def point_add(P, Q): x1 = P[0] y1 = P[1] x2 = Q[0] y2 = Q[1] if P == [0, 0]: return Q if Q == [0, 0]: return P if x1 == x2 and y1 + y2 == 0: return [0, 0] if P == Q: return point_double(P) slope = F((y2 - y1)/(x2 - x1)) x3 = F(slope**2 - x1 - x2) y3 = F(slope * (x1 - x3) - y1) return [x3, y3] Point doubling: def point_double(P): if P == [0, 0]: return [0, 0] x1 = P[0] y1 = P[1] slope = F((3 * x1**2 + a)/(2 * y1)) x3 = F(slope**2 - 2 * x1) y3 = F(slope * (x1 - x3) - y1) return [x3, y3] And for point multiplication we use the Montgomery method of multiply a point by a scalar value: def point_mul(k, P): if k >= n: k = k % n if k == 0: return [0, 0] kbits = [] while k > 0: kbits.append(k % 2) k = k // 2 knbits = len(kbits) Q = [0, 0] for i in range(knbits): Q = point_double(Q) if (kbits[knbits - i - 1] == 1): Q = point_add(Q, P) return Q CodingIn the following we will compute \(2.G\) with \(G+G\), \(2.G\) and double \(G\), and should even up with the same point. Overall, in P256, the point \(2.G\) is at (56515219790691171413109057904011688695424810155802929973526481321309856242040, 337703184371225825922371145149145259808867551975154856711245809463549758356). You can check [here] gx = 48439561293906451759052585252797914202762949526041747995844080717082404635286L gy = 36134250956749795798585127919587881956611106672985015071877198253568414405109L n = 115792089210356248762697446949407573529996955224135760342422259061068512044369L p256 = 2^256-2^224+2^192+2^96-1 a = p256 - 3 b = 41058363725152142129326129780047268409114441015993725554835256314039467401291L F = GF(p256) EC = EllipticCurve([F(a), F(b)]) EC.set_order(n) G = EC(F(gx), F(gy)) # Ref: https://github.com/asanso/CryptoWithSageMath/blob/master/ec_arithm.sage def point_add(P, Q): x1 = P[0] y1 = P[1] x2 = Q[0] y2 = Q[1] if P == [0, 0]: return Q if Q == [0, 0]: return P if x1 == x2 and y1 + y2 == 0: return [0, 0] if P == Q: return point_double(P) slope = F((y2 - y1)/(x2 - x1)) x3 = F(slope**2 - x1 - x2) y3 = F(slope * (x1 - x3) - y1) return [x3, y3] def point_double(P): if P == [0, 0]: return [0, 0] x1 = P[0] y1 = P[1] slope = F((3 * x1**2 + a)/(2 * y1)) x3 = F(slope**2 - 2 * x1) y3 = F(slope * (x1 - x3) - y1) return [x3, y3] def point_mul(k, P): if k >= n: k = k % n if k == 0: return [0, 0] kbits = [] while k > 0: kbits.append(k % 2) k = k // 2 knbits = len(kbits) Q = [0, 0] for i in range(knbits): Q = point_double(Q) if (kbits[knbits - i - 1] == 1): Q = point_add(Q, P) return Q G_2 = point_add(G, G) print (f"2G (G + G) = {G_2}") G_2 = point_mul(2, G) print (f"2G (2xG) = {G_2}") G_2 = point_double(G) print (f"2G (double G) = {G_2}") neg_G = [G[0], -G[1]] print (f"-G = {neg_G}") P = point_add(G, neg_G) print (f"G-G = {P}") priv = randint(1, n-1) pub = point_mul(priv, G) print (f"Private key: {priv}") print (f"Public key: {pub}") A sample run is: 2G (G + G) = [56515219790691171413109057904011688695424810155802929973526481321309856242040, 3377031843712258259223711451491452598088675519751548567112458094635497583569] 2G (2xG) = [56515219790691171413109057904011688695424810155802929973526481321309856242040, 3377031843712258259223711451491452598088675519751548567112458094635497583569] 2G (double G) = [56515219790691171413109057904011688695424810155802929973526481321309856242040, 3377031843712258259223711451491452598088675519751548567112458094635497583569] -G = [48439561293906451759052585252797914202762949526041747995844080717082404635286, 79657838253606452964112319029819691573475036742305299123656433055298683448842] G-G = [0, 0] Private key: 90055993546656703951884589685384927043650043173802613969260635715884780544629 Public key: [37682788206461863749596019382352165160658362606507883608751739356126844745643, 40943739187827880135926596988467882876876679195407584258200054591496014081818] |