ADDED: Reasoning about elliptic curves in library(crypto).

This is useful to establish shared secrets, using ECDH key exchange.

Note that CLP(ℤ) goal expansion is currently disabled due to #445,
and this slows down the computations considerably for the time being.

README.md

@@ -378,7 +378,8 @@ The modules that ship with Scryer Prolog are also called
 * [`sockets`](src/prolog/lib/sockets.pl)
   Predicates for opening and accepting TCP connections as streams.
 * [`crypto`](src/prolog/lib/crypto.pl)
-  Cryptographically secure random numbers and hashes.
+  Cryptographically secure random numbers and hashes, and
+  reasoning about elliptic curves.
 
 To read contents of external files, use `phrase_from_file/2` from
 [`library(pio)`](src/prolog/lib/pio.pl) to apply a DCG to

src/prolog/lib/crypto.pl

@@ -14,13 +14,19 @@
 
 :- module(crypto, [hex_bytes/2,
                    crypto_n_random_bytes/2,
-                   crypto_data_hash/3
+                   crypto_data_hash/3,
+                   crypto_name_curve/2,        % +Name, -Curve
+                   crypto_curve_order/2,       % +Curve, -Order
+                   crypto_curve_generator/2,   % +Curve, -Generator
+                   crypto_curve_scalar_mult/4  % +Curve, +Scalar, +Point, -Result
                    ]).
 
 :- use_module(library(error)).
 :- use_module(library(lists)).
 :- use_module(library(between)).
 :- use_module(library(dcgs)).
+:- use_module(library(clpz)).
+:- use_module(library(arithmetic)).
 
 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    hex_bytes(?Hex, ?Bytes) is det.
@@ -173,3 +179,195 @@ hash_algorithm(sha256).
 hash_algorithm(sha512).
 hash_algorithm(sha384).
 hash_algorithm(sha512_256).
+
+
+/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+   Modular multiplicative inverse.
+
+   Compute Y = X^(-1) mod p, using the extended Euclidean algorithm.
+- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
+
+multiplicative_inverse_modulo_p(X, P, Y) :-
+        eea(X, P, _, _, Y),
+        R #= X*Y mod P,
+        zcompare(C, 1, R),
+        must_be_one(C, X, P, Y).
+
+must_be_one(=, _, _, _).
+must_be_one(>, X, P, Y) :- throw(multiplicative_inverse_modulo_p(X,P,Y)).
+must_be_one(<, X, P, Y) :- throw(multiplicative_inverse_modulo_p(X,P,Y)).
+
+/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+   Extended Euclidean algorithm.
+
+   Computes the GCD and the Bézout coefficients S and T.
+- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
+
+eea(I, J, G, S, T) :-
+        State0 = state(1,0,0,1),
+        eea_loop(I, J, State0, G, S, T).
+
+eea_loop(I, J, State0, G, S, T) :-
+        zcompare(C, 0, J),
+        eea_(C, I, J, State0, G, S, T).
+
+eea_(=, I, _, state(_,_,U,V), I, U, V).
+eea_(<, I0, J0, state(S0,T0,U0,V0), I, U, V) :-
+        Q #= I0 // J0,
+        R #= I0 mod J0,
+        S1 #= U0 - (Q*S0),
+        T1 #= V0 - (Q*T0),
+        eea_loop(J0, R, state(S1,T1,S0,T0), I, U, V).
+
+
+/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+   Operations on Elliptic Curves
+   =============================
+
+   Sample use: Establishing a shared secret S, using ECDH key exchange.
+
+    ?- crypto_name_curve(Name, C),
+       crypto_curve_generator(C, Generator),
+       PrivateKey = 10,
+       crypto_curve_scalar_mult(C, PrivateKey, Generator, PublicKey),
+       Random = 12,
+       crypto_curve_scalar_mult(C, Random, Generator, R),
+       crypto_curve_scalar_mult(C, Random, PublicKey, S),
+       crypto_curve_scalar_mult(C, PrivateKey, R, S).
+
+- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
+
+/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+   An elliptic curve over a prime field F_p is represented as:
+
+   curve(P,A,B,point(X,Y),Order,Cofactor).
+
+   First, we define suitable accessors.
+- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
+
+curve_p(curve(P,_,_,_,_,_), P).
+curve_a(curve(_,A,_,_,_,_), A).
+curve_b(curve(_,_,B,_,_,_), B).
+
+crypto_curve_order(curve(_,_,_,_,Order,_), Order).
+crypto_curve_generator(curve(_,_,_,G,_,_), G).
+
+/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+   Scalar point multiplication.
+
+   R = k*Q.
+
+   The Montgomery ladder method is used to mitigate side-channel
+   attacks such as timing attacks, since the number of multiplications
+   and additions is independent of the private key K. This method does
+   not even reveal the key's Hamming weight (number of 1s).
+- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
+
+crypto_curve_scalar_mult(Curve, K, Q, R) :-
+        msb(K, Upper),
+        scalar_multiplication(Curve, K, Upper, ml(null,Q)-R),
+        must_be_on_curve(Curve, R).
+
+scalar_multiplication(Curve, K, I, R0-R) :-
+        zcompare(C, -1, I),
+        scalar_mult_(C, Curve, K, I, R0-R).
+
+scalar_mult_(=, _, _, _, ml(R,_)-R).
+scalar_mult_(<, Curve, K, I0, ML0-R) :-
+        BitSet #= K /\ (1 << I0),
+        zcompare(C, 0, BitSet),
+        montgomery_step(C, Curve, ML0, ML1),
+        I1 #= I0 - 1,
+        scalar_multiplication(Curve, K, I1, ML1-R).
+
+montgomery_step(=, Curve, ml(R0,S0), ml(R,S)) :-
+        curve_points_addition(Curve, R0, S0, S),
+        curve_point_double(Curve, R0, R).
+montgomery_step(<, Curve, ml(R0,S0), ml(R,S)) :-
+        curve_points_addition(Curve, R0, S0, R),
+        curve_point_double(Curve, S0, S).
+
+/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+   Doubling a point: R = A + A.
+- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
+
+curve_point_double(_, null, null).
+curve_point_double(Curve, point(AX,AY), R) :-
+        curve_p(Curve, P),
+        curve_a(Curve, A),
+        Numerator #= (3*AX^2 + A) mod P,
+        Denom0 #= 2*AY mod P,
+        multiplicative_inverse_modulo_p(Denom0, P, Denom),
+        S #= (Numerator*Denom) mod P,
+        R = point(RX,RY),
+        RX #= (S^2 - 2*AX) mod P,
+        RY #= (S*(AX - RX) - AY) mod P,
+        must_be_on_curve(Curve, R).
+
+/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+   Adding two points.
+- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
+
+curve_points_addition(Curve, P, Q, R) :-
+        curve_points_addition_(P, Curve, Q, R).
+
+curve_points_addition_(null, _, P, P).
+curve_points_addition_(P, _, null, P).
+curve_points_addition_(point(AX,AY), Curve, point(BX,BY), R) :-
+        curve_p(Curve, P),
+        Numerator #= (AY - BY) mod P,
+        Denom0 #= (AX - BX) mod P,
+        multiplicative_inverse_modulo_p(Denom0, P, Denom),
+        S #= (Numerator * Denom) mod P,
+        R = point(RX,RY),
+        RX #= (S^2 - AX - BX) mod P,
+        RY #= (S*(AX - RX) - AY) mod P,
+        must_be_on_curve(Curve, R).
+
+/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+   Validation.
+- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
+
+curve_contains_point(Curve, point(QX,QY)) :-
+        curve_a(Curve, A),
+        curve_b(Curve, B),
+        curve_p(Curve, P),
+        QY^2 mod P #= (QX^3 + A*QX + B) mod P.
+
+must_be_on_curve(Curve, P) :-
+        \+ curve_contains_point(Curve, P),
+        throw(not_on_curve(P)).
+must_be_on_curve(Curve, P) :- curve_contains_point(Curve, P).
+
+/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
+   Predefined curves
+   =================
+
+   List available curves:
+
+   $ openssl ecparam -list_curves
+
+   Show curve parameters for secp256k1:
+
+   $ openssl ecparam -param_enc explicit -conv_form uncompressed \
+                     -text -no_seed -name secp256k1
+
+   You must remove the leading "04:" from the generator.
+- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
+
+crypto_name_curve(secp112r1,
+                  curve(0x00db7c2abf62e35e668076bead208b,
+                        0x00db7c2abf62e35e668076bead2088,
+                        0x659ef8ba043916eede8911702b22,
+                        point(0x09487239995a5ee76b55f9c2f098,
+                              0xa89ce5af8724c0a23e0e0ff77500),
+                        0x00db7c2abf62e35e7628dfac6561c5,
+                        1)).
+crypto_name_curve(secp256k1,
+                  curve(0x00fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffc2f,
+                        0x0,
+                        0x7,
+                        point(0x79be667ef9dcbbac55a06295ce870b07029bfcdb2dce28d959f2815b16f81798,
+                              0x483ada7726a3c4655da4fbfc0e1108a8fd17b448a68554199c47d08ffb10d4b8),
+                        0x00fffffffffffffffffffffffffffffffebaaedce6af48a03bbfd25e8cd0364141,
+                        1)).