timoreo/olm
1
1
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity

Squashed 'lib/curve25519-donna/' content from commit 28772f3

Browse source 
git-subtree-dir: lib/curve25519-donna
git-subtree-split: 28772f37a4b8a57ab9439b9e79b19f9abee686da
This commit is contained in:
Mark Haines 2015-02-26 16:40:56 +00:00
commit e50ac70731
23 changed files with 2353 additions and 0 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

12
.gitignore vendored Normal file
View file

@ -0,0 +1,12 @@
/curve25519-donna-c64.a
/curve25519-donna.a
/test-curve25519-donna
/speed-curve25519-donna
/test-curve25519-donna-c64
/speed-curve25519-donna-c64
/test-sc-curve25519-donna-c64
/build
*.o
*.pyc
/dist
/MANIFEST

46
LICENSE.md Normal file
View file

@ -0,0 +1,46 @@
Copyright 2008, Google Inc.
All rights reserved.
Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are
met:
* Redistributions of source code must retain the above copyright
notice, this list of conditions and the following disclaimer.
* Redistributions in binary form must reproduce the above
copyright notice, this list of conditions and the following disclaimer
in the documentation and/or other materials provided with the
distribution.
* Neither the name of Google Inc. nor the names of its
contributors may be used to endorse or promote products derived from
this software without specific prior written permission.
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
curve25519-donna: Curve25519 elliptic curve, public key function
http://code.google.com/p/curve25519-donna/
Adam Langley <agl@imperialviolet.org>
Derived from public domain C code by Daniel J. Bernstein <djb@cr.yp.to>
More information about curve25519 can be found here
http://cr.yp.to/ecdh.html
djb's sample implementation of curve25519 is written in a special assembly
language called qhasm and uses the floating point registers.
This is, almost, a clean room reimplementation from the curve25519 paper. It
uses many of the tricks described therein. Only the crecip function is taken
from the sample implementation.

56
Makefile Normal file
View file

@ -0,0 +1,56 @@
CFLAGS=-Wmissing-prototypes -Wdeclaration-after-statement -O2 -Wall
CFLAGS_32=-m32
targets: curve25519-donna.a curve25519-donna-c64.a
test: test-donna test-donna-c64
clean:
rm -f *.o *.a *.pp test-curve25519-donna test-curve25519-donna-c64 speed-curve25519-donna speed-curve25519-donna-c64 test-noncanon-curve25519-donna test-noncanon-curve25519-donna-c64
curve25519-donna.a: curve25519-donna.o
ar -rc curve25519-donna.a curve25519-donna.o
ranlib curve25519-donna.a
curve25519-donna.o: curve25519-donna.c
gcc -c curve25519-donna.c $(CFLAGS) $(CFLAGS_32)
curve25519-donna-c64.a: curve25519-donna-c64.o
ar -rc curve25519-donna-c64.a curve25519-donna-c64.o
ranlib curve25519-donna-c64.a
curve25519-donna-c64.o: curve25519-donna-c64.c
gcc -c curve25519-donna-c64.c $(CFLAGS)
test-donna: test-curve25519-donna
./test-curve25519-donna | head -123456 | tail -1
test-donna-c64: test-curve25519-donna-c64
./test-curve25519-donna-c64 | head -123456 | tail -1
test-curve25519-donna: test-curve25519.c curve25519-donna.a
gcc -o test-curve25519-donna test-curve25519.c curve25519-donna.a $(CFLAGS) $(CFLAGS_32)
test-curve25519-donna-c64: test-curve25519.c curve25519-donna-c64.a
gcc -o test-curve25519-donna-c64 test-curve25519.c curve25519-donna-c64.a $(CFLAGS)
speed-curve25519-donna: speed-curve25519.c curve25519-donna.a
gcc -o speed-curve25519-donna speed-curve25519.c curve25519-donna.a $(CFLAGS) $(CFLAGS_32)
speed-curve25519-donna-c64: speed-curve25519.c curve25519-donna-c64.a
gcc -o speed-curve25519-donna-c64 speed-curve25519.c curve25519-donna-c64.a $(CFLAGS)
test-sc-curve25519-donna-c64: test-sc-curve25519.c curve25519-donna-c64.a
gcc -o test-sc-curve25519-donna-c64 -O test-sc-curve25519.c curve25519-donna-c64.a test-sc-curve25519.s $(CFLAGS)
test-noncanon-donna: test-noncanon-curve25519-donna
./test-noncanon-curve25519-donna
test-noncanon-donna-c64: test-noncanon-curve25519-donna-c64
./test-noncanon-curve25519-donna-c64
test-noncanon-curve25519-donna: test-noncanon.c curve25519-donna.a
gcc -o test-noncanon-curve25519-donna test-noncanon.c curve25519-donna.a $(CFLAGS) $(CFLAGS_32)
test-noncanon-curve25519-donna-c64: test-noncanon.c curve25519-donna-c64.a
gcc -o test-noncanon-curve25519-donna-c64 test-noncanon.c curve25519-donna-c64.a $(CFLAGS)

40
README Normal file
View file

@ -0,0 +1,40 @@
See http://code.google.com/p/curve25519-donna/ for details.
BUILDING:
If you run `make`, two .a archives will be built, similar to djb's curve25519
code. Alternatively, read on:
The C implementation is contained within curve25519-donna.c. It has no external
dependancies and is BSD licenced. You can copy/include/link it directly in with
your program. Recommended C flags: -O2
The x86-64 bit implementation is contained within curve25519-donna-x86-64.c and
curve25519-donna-x86-64.s. Build like this:
% cpp curve25519-donna-x86-64.s > curve25519-donna-x86-64.s.pp
% as -o curve25519-donna-x86-64.s.o curve25519-donna-x86-64.s.pp
% gcc -O2 -c curve25519-donna-x86-64.c
Then the two .o files can be linked in
USAGE:
The usage is exactly the same as djb's code (as described at
http://cr.yp.to/ecdh.html) expect that the function is called curve25519_donna.
In short,
To generate a private key just generate 32 random bytes.
To generate the public key, just do:
static const uint8_t basepoint[32] = {9};
curve25519_donna(mypublic, mysecret, basepoint);
To generate an agreed key do:
uint8_t shared_key[32];
curve25519_donna(shared_key, mysecret, theirpublic);
And hash the shared_key with a cryptographic hash function before using.

118
contrib/Curve25519Donna.c Normal file
View file

@ -0,0 +1,118 @@
/*
James Robson
Public domain.
*/
#include "Curve25519Donna.h"
#include <stdio.h>
#include <stdlib.h>
extern void curve25519_donna(unsigned char *output, const unsigned char *a,
const unsigned char *b);
unsigned char*
as_unsigned_char_array(JNIEnv* env, jbyteArray array, int* len);
jbyteArray as_byte_array(JNIEnv* env, unsigned char* buf, int len);
jbyteArray as_byte_array(JNIEnv* env, unsigned char* buf, int len) {
jbyteArray array = (*env)->NewByteArray(env, len);
(*env)->SetByteArrayRegion(env, array, 0, len, (jbyte*)buf);
//int i;
//for (i = 0;i < len;++i) printf("%02x",(unsigned int) buf[i]); printf(" ");
//printf("\n");
return array;
}
unsigned char*
as_unsigned_char_array(JNIEnv* env, jbyteArray array, int* len) {
*len = (*env)->GetArrayLength(env, array);
unsigned char* buf = (unsigned char*)calloc(*len+1, sizeof(char));
(*env)->GetByteArrayRegion (env, array, 0, *len, (jbyte*)buf);
return buf;
}
JNIEXPORT jbyteArray JNICALL Java_Curve25519Donna_curve25519Donna
(JNIEnv *env, jobject obj, jbyteArray a, jbyteArray b) {
unsigned char o[32] = {0};
int l1, l2;
unsigned char* a1 = as_unsigned_char_array(env, a, &l1);
unsigned char* b1 = as_unsigned_char_array(env, b, &l2);
if ( !(l1 == 32 && l2 == 32) ) {
fprintf(stderr, "Error, must be length 32");
return NULL;
}
curve25519_donna(o, (const unsigned char*)a1, (const unsigned char*)b1);
free(a1);
free(b1);
return as_byte_array(env, (unsigned char*)o, 32);
}
JNIEXPORT jbyteArray JNICALL Java_Curve25519Donna_makePrivate
(JNIEnv *env, jobject obj, jbyteArray secret) {
int len;
unsigned char* k = as_unsigned_char_array(env, secret, &len);
if (len != 32) {
fprintf(stderr, "Error, must be length 32");
return NULL;
}
k[0] &= 248;
k[31] &= 127;
k[31] |= 64;
return as_byte_array(env, k, 32);
}
JNIEXPORT jbyteArray JNICALL Java_Curve25519Donna_getPublic
(JNIEnv *env, jobject obj, jbyteArray privkey) {
int len;
unsigned char* private = as_unsigned_char_array(env, privkey, &len);
if (len != 32) {
fprintf(stderr, "Error, must be length 32");
return NULL;
}
unsigned char pubkey[32];
unsigned char basepoint[32] = {9};
curve25519_donna(pubkey, private, basepoint);
return as_byte_array(env, (unsigned char*)pubkey, 32);
}
JNIEXPORT jbyteArray JNICALL Java_Curve25519Donna_makeSharedSecret
(JNIEnv *env, jobject obj, jbyteArray privkey, jbyteArray their_pubkey) {
unsigned char shared_secret[32];
int l1, l2;
unsigned char* private = as_unsigned_char_array(env, privkey, &l1);
unsigned char* pubkey = as_unsigned_char_array(env, their_pubkey, &l2);
if ( !(l1 == 32 && l2 == 32) ) {
fprintf(stderr, "Error, must be length 32");
return NULL;
}
curve25519_donna(shared_secret, private, pubkey);
return as_byte_array(env, (unsigned char*)shared_secret, 32);
}
JNIEXPORT void JNICALL Java_Curve25519Donna_helowrld
(JNIEnv *env, jobject obj) {
printf("helowrld\n");
}

53
contrib/Curve25519Donna.h Normal file
View file

@ -0,0 +1,53 @@
/* DO NOT EDIT THIS FILE - it is machine generated */
#include <jni.h>
/* Header for class Curve25519Donna */
#ifndef _Included_Curve25519Donna
#define _Included_Curve25519Donna
#ifdef __cplusplus
extern "C" {
#endif
/*
* Class: Curve25519Donna
* Method: curve25519Donna
* Signature: ([B[B)[B
*/
JNIEXPORT jbyteArray JNICALL Java_Curve25519Donna_curve25519Donna
(JNIEnv *, jobject, jbyteArray, jbyteArray);
/*
* Class: Curve25519Donna
* Method: makePrivate
* Signature: ([B)[B
*/
JNIEXPORT jbyteArray JNICALL Java_Curve25519Donna_makePrivate
(JNIEnv *, jobject, jbyteArray);
/*
* Class: Curve25519Donna
* Method: getPublic
* Signature: ([B)[B
*/
JNIEXPORT jbyteArray JNICALL Java_Curve25519Donna_getPublic
(JNIEnv *, jobject, jbyteArray);
/*
* Class: Curve25519Donna
* Method: makeSharedSecret
* Signature: ([B[B)[B
*/
JNIEXPORT jbyteArray JNICALL Java_Curve25519Donna_makeSharedSecret
(JNIEnv *, jobject, jbyteArray, jbyteArray);
/*
* Class: Curve25519Donna
* Method: helowrld
* Signature: ()V
*/
JNIEXPORT void JNICALL Java_Curve25519Donna_helowrld
(JNIEnv *, jobject);
#ifdef __cplusplus
}
#endif
#endif

77
contrib/Curve25519Donna.java Normal file
View file

@ -0,0 +1,77 @@
/*
James Robson
Public domain.
*/
public class Curve25519Donna {
final protected static char[] hexArray = "0123456789ABCDEF".toCharArray();
public static String bytesToHex(byte[] bytes) {
char[] hexChars = new char[bytes.length * 2];
int v;
for ( int j = 0; j < bytes.length; j++ ) {
v = bytes[j] & 0xFF;
hexChars[j * 2] = hexArray[v >>> 4];
hexChars[j * 2 + 1] = hexArray[v & 0x0F];
}
return new String(hexChars);
}
public native byte[] curve25519Donna(byte[] a, byte[] b);
public native byte[] makePrivate(byte[] secret);
public native byte[] getPublic(byte[] privkey);
public native byte[] makeSharedSecret(byte[] privkey, byte[] theirPubKey);
public native void helowrld();
// Uncomment if your Java is 32-bit:
//static { System.loadLibrary("Curve25519Donna"); }
// Otherwise, load this 64-bit .jnilib:
static { System.loadLibrary("Curve25519Donna_64"); }
/*
To give the old tires a kick (OSX):
java -cp `pwd` Curve25519Donna
*/
public static void main (String[] args) {
Curve25519Donna c = new Curve25519Donna();
// These should be 32 bytes long
byte[] user1Secret = "abcdefghijklmnopqrstuvwxyz123456".getBytes();
byte[] user2Secret = "654321zyxwvutsrqponmlkjihgfedcba".getBytes();
// You can use the curve function directly...
//byte[] o = c.curve25519Donna(a, b);
//System.out.println("o = " + bytesToHex(o));
// ... but it's not really necessary. Just use the following
// convenience methods:
byte[] privKey = c.makePrivate(user1Secret);
byte[] pubKey = c.getPublic(privKey);
byte[] privKey2 = c.makePrivate(user2Secret);
byte[] pubKey2 = c.getPublic(privKey2);
System.out.println("'user1' privKey = " + bytesToHex(privKey));
System.out.println("'user1' pubKey = " + bytesToHex(pubKey));
System.out.println("===================================================");
System.out.println("'user2' privKey = " + bytesToHex(privKey2));
System.out.println("'user2' pubKey = " + bytesToHex(pubKey2));
System.out.println("===================================================");
byte[] ss1 = c.makeSharedSecret(privKey, pubKey2);
System.out.println("'user1' computes shared secret: " + bytesToHex(ss1));
byte[] ss2 = c.makeSharedSecret(privKey2, pubKey);
System.out.println("'user2' computes shared secret: " + bytesToHex(ss2));
}
}

68
contrib/make-snippets Normal file
View file

@ -0,0 +1,68 @@
CFLAGS=-Wmissing-prototypes -Wdeclaration-after-statement -O2 -Wall
CC=clang
targets: curve25519-donna.a curve25519-donna-c64.a
test: test-donna test-donna-c64
clean:
rm -f java-src/*.class java-src/*.jnilib *.dylib *.o *.a *.pp test-curve25519-donna test-curve25519-donna-c64 speed-curve25519-donna speed-curve25519-donna-c64
curve25519-donna.a: curve25519-donna.o
ar -rc curve25519-donna.a curve25519-donna.o
ranlib curve25519-donna.a
##### OSX dynamic library (32- & 64-bit)
curve25519donna.dylib: curve25519-donna.a curve25519-donna-c64.a
$(CC) -m32 -fpic -shared -Wl,-all_load curve25519-donna.a -Wl,-all_load -o libcurve25519donna.dylib
$(CC) -fpic -shared -Wl,-all_load curve25519-donna-c64.a -Wl,-all_load -o libcurve25519donna_64.dylib
##### OSX/Java section hence
# Java JNI - compiled for OSX (32- & 64-bit)
Curve25519Donna.class:
cd java-src; javah -jni Curve25519Donna; cd ..
cd java-src; javac Curve25519Donna.java; cd ..
Curve25519Donna.jnilib: curve25519-donna.a curve25519-donna-c64.a Curve25519Donna.class
@echo "Building 32-bit..."
clang -o java-src/libCurve25519Donna.jnilib $(CFLAGS) -lc -shared -m32 -I /System/Library/Frameworks/JavaVM.framework/Headers curve25519-donna.o java-src/Curve25519Donna.c
@echo "Building 64-bit..."
clang -o java-src/libCurve25519Donna_64.jnilib $(CFLAGS) -lc -shared -I /System/Library/Frameworks/JavaVM.framework/Headers curve25519-donna-c64.o java-src/Curve25519Donna.c
##### OSX/Java section end
curve25519-donna.o: curve25519-donna.c
$(CC) -c curve25519-donna.c $(CFLAGS) -m32
curve25519-donna-c64.a: curve25519-donna-c64.o
ar -rc curve25519-donna-c64.a curve25519-donna-c64.o
ranlib curve25519-donna-c64.a
curve25519-donna-c64.o: curve25519-donna-c64.c
$(CC) -c curve25519-donna-c64.c $(CFLAGS)
test-donna: test-curve25519-donna
./test-curve25519-donna | head -123456 | tail -1
test-donna-c64: test-curve25519-donna-c64
./test-curve25519-donna-c64 | head -123456 | tail -1
test-curve25519-donna: test-curve25519.c curve25519-donna.a
$(CC) -o test-curve25519-donna test-curve25519.c curve25519-donna.a $(CFLAGS) -m32
test-curve25519-donna-c64: test-curve25519.c curve25519-donna-c64.a
$(CC) -o test-curve25519-donna-c64 test-curve25519.c curve25519-donna-c64.a $(CFLAGS)
speed-curve25519-donna: speed-curve25519.c curve25519-donna.a
$(CC) -o speed-curve25519-donna speed-curve25519.c curve25519-donna.a $(CFLAGS) -m32
speed-curve25519-donna-c64: speed-curve25519.c curve25519-donna-c64.a
$(CC) -o speed-curve25519-donna-c64 speed-curve25519.c curve25519-donna-c64.a $(CFLAGS)
test-sc-curve25519-donna-c64: test-sc-curve25519.c curve25519-donna-c64.a
$(CC) -o test-sc-curve25519-donna-c64 -O test-sc-curve25519.c curve25519-donna-c64.a test-sc-curve25519.s $(CFLAGS)

449
curve25519-donna-c64.c Normal file
View file

@ -0,0 +1,449 @@
/* Copyright 2008, Google Inc.
* All rights reserved.
*
* Code released into the public domain.
*
* curve25519-donna: Curve25519 elliptic curve, public key function
*
* http://code.google.com/p/curve25519-donna/
*
* Adam Langley <agl@imperialviolet.org>
*
* Derived from public domain C code by Daniel J. Bernstein <djb@cr.yp.to>
*
* More information about curve25519 can be found here
* http://cr.yp.to/ecdh.html
*
* djb's sample implementation of curve25519 is written in a special assembly
* language called qhasm and uses the floating point registers.
*
* This is, almost, a clean room reimplementation from the curve25519 paper. It
* uses many of the tricks described therein. Only the crecip function is taken
* from the sample implementation.
*/
#include <string.h>
#include <stdint.h>
typedef uint8_t u8;
typedef uint64_t limb;
typedef limb felem[5];
// This is a special gcc mode for 128-bit integers. It's implemented on 64-bit
// platforms only as far as I know.
typedef unsigned uint128_t __attribute__((mode(TI)));
#undef force_inline
#define force_inline __attribute__((always_inline))
/* Sum two numbers: output += in */
static inline void force_inline
fsum(limb *output, const limb *in) {
output[0] += in[0];
output[1] += in[1];
output[2] += in[2];
output[3] += in[3];
output[4] += in[4];
}
/* Find the difference of two numbers: output = in - output
* (note the order of the arguments!)
*
* Assumes that out[i] < 2**52
* On return, out[i] < 2**55
*/
static inline void force_inline
fdifference_backwards(felem out, const felem in) {
/* 152 is 19 << 3 */
static const limb two54m152 = (((limb)1) << 54) - 152;
static const limb two54m8 = (((limb)1) << 54) - 8;
out[0] = in[0] + two54m152 - out[0];
out[1] = in[1] + two54m8 - out[1];
out[2] = in[2] + two54m8 - out[2];
out[3] = in[3] + two54m8 - out[3];
out[4] = in[4] + two54m8 - out[4];
}
/* Multiply a number by a scalar: output = in * scalar */
static inline void force_inline
fscalar_product(felem output, const felem in, const limb scalar) {
uint128_t a;
a = ((uint128_t) in[0]) * scalar;
output[0] = ((limb)a) & 0x7ffffffffffff;
a = ((uint128_t) in[1]) * scalar + ((limb) (a >> 51));
output[1] = ((limb)a) & 0x7ffffffffffff;
a = ((uint128_t) in[2]) * scalar + ((limb) (a >> 51));
output[2] = ((limb)a) & 0x7ffffffffffff;
a = ((uint128_t) in[3]) * scalar + ((limb) (a >> 51));
output[3] = ((limb)a) & 0x7ffffffffffff;
a = ((uint128_t) in[4]) * scalar + ((limb) (a >> 51));
output[4] = ((limb)a) & 0x7ffffffffffff;
output[0] += (a >> 51) * 19;
}
/* Multiply two numbers: output = in2 * in
*
* output must be distinct to both inputs. The inputs are reduced coefficient
* form, the output is not.
*
* Assumes that in[i] < 2**55 and likewise for in2.
* On return, output[i] < 2**52
*/
static inline void force_inline
fmul(felem output, const felem in2, const felem in) {
uint128_t t[5];
limb r0,r1,r2,r3,r4,s0,s1,s2,s3,s4,c;
r0 = in[0];
r1 = in[1];
r2 = in[2];
r3 = in[3];
r4 = in[4];
s0 = in2[0];
s1 = in2[1];
s2 = in2[2];
s3 = in2[3];
s4 = in2[4];
t[0] = ((uint128_t) r0) * s0;
t[1] = ((uint128_t) r0) * s1 + ((uint128_t) r1) * s0;
t[2] = ((uint128_t) r0) * s2 + ((uint128_t) r2) * s0 + ((uint128_t) r1) * s1;
t[3] = ((uint128_t) r0) * s3 + ((uint128_t) r3) * s0 + ((uint128_t) r1) * s2 + ((uint128_t) r2) * s1;
t[4] = ((uint128_t) r0) * s4 + ((uint128_t) r4) * s0 + ((uint128_t) r3) * s1 + ((uint128_t) r1) * s3 + ((uint128_t) r2) * s2;
r4 *= 19;
r1 *= 19;
r2 *= 19;
r3 *= 19;
t[0] += ((uint128_t) r4) * s1 + ((uint128_t) r1) * s4 + ((uint128_t) r2) * s3 + ((uint128_t) r3) * s2;
t[1] += ((uint128_t) r4) * s2 + ((uint128_t) r2) * s4 + ((uint128_t) r3) * s3;
t[2] += ((uint128_t) r4) * s3 + ((uint128_t) r3) * s4;
t[3] += ((uint128_t) r4) * s4;
r0 = (limb)t[0] & 0x7ffffffffffff; c = (limb)(t[0] >> 51);
t[1] += c; r1 = (limb)t[1] & 0x7ffffffffffff; c = (limb)(t[1] >> 51);
t[2] += c; r2 = (limb)t[2] & 0x7ffffffffffff; c = (limb)(t[2] >> 51);
t[3] += c; r3 = (limb)t[3] & 0x7ffffffffffff; c = (limb)(t[3] >> 51);
t[4] += c; r4 = (limb)t[4] & 0x7ffffffffffff; c = (limb)(t[4] >> 51);
r0 += c * 19; c = r0 >> 51; r0 = r0 & 0x7ffffffffffff;
r1 += c; c = r1 >> 51; r1 = r1 & 0x7ffffffffffff;
r2 += c;
output[0] = r0;
output[1] = r1;
output[2] = r2;
output[3] = r3;
output[4] = r4;
}
static inline void force_inline
fsquare_times(felem output, const felem in, limb count) {
uint128_t t[5];
limb r0,r1,r2,r3,r4,c;
limb d0,d1,d2,d4,d419;
r0 = in[0];
r1 = in[1];
r2 = in[2];
r3 = in[3];
r4 = in[4];
do {
d0 = r0 * 2;
d1 = r1 * 2;
d2 = r2 * 2 * 19;
d419 = r4 * 19;
d4 = d419 * 2;
t[0] = ((uint128_t) r0) * r0 + ((uint128_t) d4) * r1 + (((uint128_t) d2) * (r3 ));
t[1] = ((uint128_t) d0) * r1 + ((uint128_t) d4) * r2 + (((uint128_t) r3) * (r3 * 19));
t[2] = ((uint128_t) d0) * r2 + ((uint128_t) r1) * r1 + (((uint128_t) d4) * (r3 ));
t[3] = ((uint128_t) d0) * r3 + ((uint128_t) d1) * r2 + (((uint128_t) r4) * (d419 ));
t[4] = ((uint128_t) d0) * r4 + ((uint128_t) d1) * r3 + (((uint128_t) r2) * (r2 ));
r0 = (limb)t[0] & 0x7ffffffffffff; c = (limb)(t[0] >> 51);
t[1] += c; r1 = (limb)t[1] & 0x7ffffffffffff; c = (limb)(t[1] >> 51);
t[2] += c; r2 = (limb)t[2] & 0x7ffffffffffff; c = (limb)(t[2] >> 51);
t[3] += c; r3 = (limb)t[3] & 0x7ffffffffffff; c = (limb)(t[3] >> 51);
t[4] += c; r4 = (limb)t[4] & 0x7ffffffffffff; c = (limb)(t[4] >> 51);
r0 += c * 19; c = r0 >> 51; r0 = r0 & 0x7ffffffffffff;
r1 += c; c = r1 >> 51; r1 = r1 & 0x7ffffffffffff;
r2 += c;
} while(--count);
output[0] = r0;
output[1] = r1;
output[2] = r2;
output[3] = r3;
output[4] = r4;
}
/* Load a little-endian 64-bit number */
static limb
load_limb(const u8 *in) {
return
((limb)in[0]) |
(((limb)in[1]) << 8) |
(((limb)in[2]) << 16) |
(((limb)in[3]) << 24) |
(((limb)in[4]) << 32) |
(((limb)in[5]) << 40) |
(((limb)in[6]) << 48) |
(((limb)in[7]) << 56);
}
static void
store_limb(u8 *out, limb in) {
out[0] = in & 0xff;
out[1] = (in >> 8) & 0xff;
out[2] = (in >> 16) & 0xff;
out[3] = (in >> 24) & 0xff;
out[4] = (in >> 32) & 0xff;
out[5] = (in >> 40) & 0xff;
out[6] = (in >> 48) & 0xff;
out[7] = (in >> 56) & 0xff;
}
/* Take a little-endian, 32-byte number and expand it into polynomial form */
static void
fexpand(limb *output, const u8 *in) {
output[0] = load_limb(in) & 0x7ffffffffffff;
output[1] = (load_limb(in+6) >> 3) & 0x7ffffffffffff;
output[2] = (load_limb(in+12) >> 6) & 0x7ffffffffffff;
output[3] = (load_limb(in+19) >> 1) & 0x7ffffffffffff;
output[4] = (load_limb(in+24) >> 12) & 0x7ffffffffffff;
}
/* Take a fully reduced polynomial form number and contract it into a
* little-endian, 32-byte array
*/
static void
fcontract(u8 *output, const felem input) {
uint128_t t[5];
t[0] = input[0];
t[1] = input[1];
t[2] = input[2];
t[3] = input[3];
t[4] = input[4];
t[1] += t[0] >> 51; t[0] &= 0x7ffffffffffff;
t[2] += t[1] >> 51; t[1] &= 0x7ffffffffffff;
t[3] += t[2] >> 51; t[2] &= 0x7ffffffffffff;
t[4] += t[3] >> 51; t[3] &= 0x7ffffffffffff;
t[0] += 19 * (t[4] >> 51); t[4] &= 0x7ffffffffffff;
t[1] += t[0] >> 51; t[0] &= 0x7ffffffffffff;
t[2] += t[1] >> 51; t[1] &= 0x7ffffffffffff;
t[3] += t[2] >> 51; t[2] &= 0x7ffffffffffff;
t[4] += t[3] >> 51; t[3] &= 0x7ffffffffffff;
t[0] += 19 * (t[4] >> 51); t[4] &= 0x7ffffffffffff;
/* now t is between 0 and 2^255-1, properly carried. */
/* case 1: between 0 and 2^255-20. case 2: between 2^255-19 and 2^255-1. */
t[0] += 19;
t[1] += t[0] >> 51; t[0] &= 0x7ffffffffffff;
t[2] += t[1] >> 51; t[1] &= 0x7ffffffffffff;
t[3] += t[2] >> 51; t[2] &= 0x7ffffffffffff;
t[4] += t[3] >> 51; t[3] &= 0x7ffffffffffff;
t[0] += 19 * (t[4] >> 51); t[4] &= 0x7ffffffffffff;
/* now between 19 and 2^255-1 in both cases, and offset by 19. */
t[0] += 0x8000000000000 - 19;
t[1] += 0x8000000000000 - 1;
t[2] += 0x8000000000000 - 1;
t[3] += 0x8000000000000 - 1;
t[4] += 0x8000000000000 - 1;
/* now between 2^255 and 2^256-20, and offset by 2^255. */
t[1] += t[0] >> 51; t[0] &= 0x7ffffffffffff;
t[2] += t[1] >> 51; t[1] &= 0x7ffffffffffff;
t[3] += t[2] >> 51; t[2] &= 0x7ffffffffffff;
t[4] += t[3] >> 51; t[3] &= 0x7ffffffffffff;
t[4] &= 0x7ffffffffffff;
store_limb(output, t[0] | (t[1] << 51));
store_limb(output+8, (t[1] >> 13) | (t[2] << 38));
store_limb(output+16, (t[2] >> 26) | (t[3] << 25));
store_limb(output+24, (t[3] >> 39) | (t[4] << 12));
}
/* Input: Q, Q', Q-Q'
* Output: 2Q, Q+Q'
*
* x2 z3: long form
* x3 z3: long form
* x z: short form, destroyed
* xprime zprime: short form, destroyed
* qmqp: short form, preserved
*/
static void
fmonty(limb *x2, limb *z2, /* output 2Q */
limb *x3, limb *z3, /* output Q + Q' */
limb *x, limb *z, /* input Q */
limb *xprime, limb *zprime, /* input Q' */
const limb *qmqp /* input Q - Q' */) {
limb origx[5], origxprime[5], zzz[5], xx[5], zz[5], xxprime[5],
zzprime[5], zzzprime[5];
memcpy(origx, x, 5 * sizeof(limb));
fsum(x, z);
fdifference_backwards(z, origx); // does x - z
memcpy(origxprime, xprime, sizeof(limb) * 5);
fsum(xprime, zprime);
fdifference_backwards(zprime, origxprime);
fmul(xxprime, xprime, z);
fmul(zzprime, x, zprime);
memcpy(origxprime, xxprime, sizeof(limb) * 5);
fsum(xxprime, zzprime);
fdifference_backwards(zzprime, origxprime);
fsquare_times(x3, xxprime, 1);
fsquare_times(zzzprime, zzprime, 1);
fmul(z3, zzzprime, qmqp);
fsquare_times(xx, x, 1);
fsquare_times(zz, z, 1);
fmul(x2, xx, zz);
fdifference_backwards(zz, xx); // does zz = xx - zz
fscalar_product(zzz, zz, 121665);
fsum(zzz, xx);
fmul(z2, zz, zzz);
}
// -----------------------------------------------------------------------------
// Maybe swap the contents of two limb arrays (@a and @b), each @len elements
// long. Perform the swap iff @swap is non-zero.
//
// This function performs the swap without leaking any side-channel
// information.
// -----------------------------------------------------------------------------
static void
swap_conditional(limb a[5], limb b[5], limb iswap) {
unsigned i;
const limb swap = -iswap;
for (i = 0; i < 5; ++i) {
const limb x = swap & (a[i] ^ b[i]);
a[i] ^= x;
b[i] ^= x;
}
}
/* Calculates nQ where Q is the x-coordinate of a point on the curve
*
* resultx/resultz: the x coordinate of the resulting curve point (short form)
* n: a little endian, 32-byte number
* q: a point of the curve (short form)
*/
static void
cmult(limb *resultx, limb *resultz, const u8 *n, const limb *q) {
limb a[5] = {0}, b[5] = {1}, c[5] = {1}, d[5] = {0};
limb *nqpqx = a, *nqpqz = b, *nqx = c, *nqz = d, *t;
limb e[5] = {0}, f[5] = {1}, g[5] = {0}, h[5] = {1};
limb *nqpqx2 = e, *nqpqz2 = f, *nqx2 = g, *nqz2 = h;
unsigned i, j;
memcpy(nqpqx, q, sizeof(limb) * 5);
for (i = 0; i < 32; ++i) {
u8 byte = n[31 - i];
for (j = 0; j < 8; ++j) {
const limb bit = byte >> 7;
swap_conditional(nqx, nqpqx, bit);
swap_conditional(nqz, nqpqz, bit);
fmonty(nqx2, nqz2,
nqpqx2, nqpqz2,
nqx, nqz,
nqpqx, nqpqz,
q);
swap_conditional(nqx2, nqpqx2, bit);
swap_conditional(nqz2, nqpqz2, bit);
t = nqx;
nqx = nqx2;
nqx2 = t;
t = nqz;
nqz = nqz2;
nqz2 = t;
t = nqpqx;
nqpqx = nqpqx2;
nqpqx2 = t;
t = nqpqz;
nqpqz = nqpqz2;
nqpqz2 = t;
byte <<= 1;
}
}
memcpy(resultx, nqx, sizeof(limb) * 5);
memcpy(resultz, nqz, sizeof(limb) * 5);
}
// -----------------------------------------------------------------------------
// Shamelessly copied from djb's code, tightened a little
// -----------------------------------------------------------------------------
static void
crecip(felem out, const felem z) {
felem a,t0,b,c;
/* 2 */ fsquare_times(a, z, 1); // a = 2
/* 8 */ fsquare_times(t0, a, 2);
/* 9 */ fmul(b, t0, z); // b = 9
/* 11 */ fmul(a, b, a); // a = 11
/* 22 */ fsquare_times(t0, a, 1);
/* 2^5 - 2^0 = 31 */ fmul(b, t0, b);
/* 2^10 - 2^5 */ fsquare_times(t0, b, 5);
/* 2^10 - 2^0 */ fmul(b, t0, b);
/* 2^20 - 2^10 */ fsquare_times(t0, b, 10);
/* 2^20 - 2^0 */ fmul(c, t0, b);
/* 2^40 - 2^20 */ fsquare_times(t0, c, 20);
/* 2^40 - 2^0 */ fmul(t0, t0, c);
/* 2^50 - 2^10 */ fsquare_times(t0, t0, 10);
/* 2^50 - 2^0 */ fmul(b, t0, b);
/* 2^100 - 2^50 */ fsquare_times(t0, b, 50);
/* 2^100 - 2^0 */ fmul(c, t0, b);
/* 2^200 - 2^100 */ fsquare_times(t0, c, 100);
/* 2^200 - 2^0 */ fmul(t0, t0, c);
/* 2^250 - 2^50 */ fsquare_times(t0, t0, 50);
/* 2^250 - 2^0 */ fmul(t0, t0, b);
/* 2^255 - 2^5 */ fsquare_times(t0, t0, 5);
/* 2^255 - 21 */ fmul(out, t0, a);
}
int curve25519_donna(u8 *, const u8 *, const u8 *);
int
curve25519_donna(u8 *mypublic, const u8 *secret, const u8 *basepoint) {
limb bp[5], x[5], z[5], zmone[5];
uint8_t e[32];
int i;
for (i = 0;i < 32;++i) e[i] = secret[i];
e[0] &= 248;
e[31] &= 127;
e[31] |= 64;
fexpand(bp, basepoint);
cmult(x, z, e, bp);
crecip(zmone, z);
fmul(z, x, zmone);
fcontract(mypublic, z);
return 0;
}

860
curve25519-donna.c Normal file
View file

@ -0,0 +1,860 @@
/* Copyright 2008, Google Inc.
* All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions are
* met:
*
* * Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
* * Redistributions in binary form must reproduce the above
* copyright notice, this list of conditions and the following disclaimer
* in the documentation and/or other materials provided with the
* distribution.
* * Neither the name of Google Inc. nor the names of its
* contributors may be used to endorse or promote products derived from
* this software without specific prior written permission.
*
* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
* "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
* LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
* A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
* OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
* SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
* LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
* DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
* THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
* (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
* OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*
* curve25519-donna: Curve25519 elliptic curve, public key function
*
* http://code.google.com/p/curve25519-donna/
*
* Adam Langley <agl@imperialviolet.org>
*
* Derived from public domain C code by Daniel J. Bernstein <djb@cr.yp.to>
*
* More information about curve25519 can be found here
* http://cr.yp.to/ecdh.html
*
* djb's sample implementation of curve25519 is written in a special assembly
* language called qhasm and uses the floating point registers.
*
* This is, almost, a clean room reimplementation from the curve25519 paper. It
* uses many of the tricks described therein. Only the crecip function is taken
* from the sample implementation. */
#include <string.h>
#include <stdint.h>
#ifdef _MSC_VER
#define inline __inline
#endif
typedef uint8_t u8;
typedef int32_t s32;
typedef int64_t limb;
/* Field element representation:
*
* Field elements are written as an array of signed, 64-bit limbs, least
* significant first. The value of the field element is:
* x[0] + 2^26·x[1] + x^51·x[2] + 2^102·x[3] + ...
*
* i.e. the limbs are 26, 25, 26, 25, ... bits wide. */
/* Sum two numbers: output += in */
static void fsum(limb *output, const limb *in) {
unsigned i;
for (i = 0; i < 10; i += 2) {
output[0+i] = output[0+i] + in[0+i];
output[1+i] = output[1+i] + in[1+i];
}
}
/* Find the difference of two numbers: output = in - output
* (note the order of the arguments!). */
static void fdifference(limb *output, const limb *in) {
unsigned i;
for (i = 0; i < 10; ++i) {
output[i] = in[i] - output[i];
}
}
/* Multiply a number by a scalar: output = in * scalar */
static void fscalar_product(limb *output, const limb *in, const limb scalar) {
unsigned i;
for (i = 0; i < 10; ++i) {
output[i] = in[i] * scalar;
}
}
/* Multiply two numbers: output = in2 * in
*
* output must be distinct to both inputs. The inputs are reduced coefficient
* form, the output is not.
*
* output[x] <= 14 * the largest product of the input limbs. */
static void fproduct(limb *output, const limb *in2, const limb *in) {
output[0] = ((limb) ((s32) in2[0])) * ((s32) in[0]);
output[1] = ((limb) ((s32) in2[0])) * ((s32) in[1]) +
((limb) ((s32) in2[1])) * ((s32) in[0]);
output[2] = 2 * ((limb) ((s32) in2[1])) * ((s32) in[1]) +
((limb) ((s32) in2[0])) * ((s32) in[2]) +
((limb) ((s32) in2[2])) * ((s32) in[0]);
output[3] = ((limb) ((s32) in2[1])) * ((s32) in[2]) +
((limb) ((s32) in2[2])) * ((s32) in[1]) +
((limb) ((s32) in2[0])) * ((s32) in[3]) +
((limb) ((s32) in2[3])) * ((s32) in[0]);
output[4] = ((limb) ((s32) in2[2])) * ((s32) in[2]) +
2 * (((limb) ((s32) in2[1])) * ((s32) in[3]) +
((limb) ((s32) in2[3])) * ((s32) in[1])) +
((limb) ((s32) in2[0])) * ((s32) in[4]) +
((limb) ((s32) in2[4])) * ((s32) in[0]);
output[5] = ((limb) ((s32) in2[2])) * ((s32) in[3]) +
((limb) ((s32) in2[3])) * ((s32) in[2]) +
((limb) ((s32) in2[1])) * ((s32) in[4]) +
((limb) ((s32) in2[4])) * ((s32) in[1]) +
((limb) ((s32) in2[0])) * ((s32) in[5]) +
((limb) ((s32) in2[5])) * ((s32) in[0]);
output[6] = 2 * (((limb) ((s32) in2[3])) * ((s32) in[3]) +
((limb) ((s32) in2[1])) * ((s32) in[5]) +
((limb) ((s32) in2[5])) * ((s32) in[1])) +
((limb) ((s32) in2[2])) * ((s32) in[4]) +
((limb) ((s32) in2[4])) * ((s32) in[2]) +
((limb) ((s32) in2[0])) * ((s32) in[6]) +
((limb) ((s32) in2[6])) * ((s32) in[0]);
output[7] = ((limb) ((s32) in2[3])) * ((s32) in[4]) +
((limb) ((s32) in2[4])) * ((s32) in[3]) +
((limb) ((s32) in2[2])) * ((s32) in[5]) +
((limb) ((s32) in2[5])) * ((s32) in[2]) +
((limb) ((s32) in2[1])) * ((s32) in[6]) +
((limb) ((s32) in2[6])) * ((s32) in[1]) +
((limb) ((s32) in2[0])) * ((s32) in[7]) +
((limb) ((s32) in2[7])) * ((s32) in[0]);
output[8] = ((limb) ((s32) in2[4])) * ((s32) in[4]) +
2 * (((limb) ((s32) in2[3])) * ((s32) in[5]) +
((limb) ((s32) in2[5])) * ((s32) in[3]) +
((limb) ((s32) in2[1])) * ((s32) in[7]) +
((limb) ((s32) in2[7])) * ((s32) in[1])) +
((limb) ((s32) in2[2])) * ((s32) in[6]) +
((limb) ((s32) in2[6])) * ((s32) in[2]) +
((limb) ((s32) in2[0])) * ((s32) in[8]) +
((limb) ((s32) in2[8])) * ((s32) in[0]);
output[9] = ((limb) ((s32) in2[4])) * ((s32) in[5]) +
((limb) ((s32) in2[5])) * ((s32) in[4]) +
((limb) ((s32) in2[3])) * ((s32) in[6]) +
((limb) ((s32) in2[6])) * ((s32) in[3]) +
((limb) ((s32) in2[2])) * ((s32) in[7]) +
((limb) ((s32) in2[7])) * ((s32) in[2]) +
((limb) ((s32) in2[1])) * ((s32) in[8]) +
((limb) ((s32) in2[8])) * ((s32) in[1]) +
((limb) ((s32) in2[0])) * ((s32) in[9]) +
((limb) ((s32) in2[9])) * ((s32) in[0]);
output[10] = 2 * (((limb) ((s32) in2[5])) * ((s32) in[5]) +
((limb) ((s32) in2[3])) * ((s32) in[7]) +
((limb) ((s32) in2[7])) * ((s32) in[3]) +
((limb) ((s32) in2[1])) * ((s32) in[9]) +
((limb) ((s32) in2[9])) * ((s32) in[1])) +
((limb) ((s32) in2[4])) * ((s32) in[6]) +
((limb) ((s32) in2[6])) * ((s32) in[4]) +
((limb) ((s32) in2[2])) * ((s32) in[8]) +
((limb) ((s32) in2[8])) * ((s32) in[2]);
output[11] = ((limb) ((s32) in2[5])) * ((s32) in[6]) +
((limb) ((s32) in2[6])) * ((s32) in[5]) +
((limb) ((s32) in2[4])) * ((s32) in[7]) +
((limb) ((s32) in2[7])) * ((s32) in[4]) +
((limb) ((s32) in2[3])) * ((s32) in[8]) +
((limb) ((s32) in2[8])) * ((s32) in[3]) +
((limb) ((s32) in2[2])) * ((s32) in[9]) +
((limb) ((s32) in2[9])) * ((s32) in[2]);
output[12] = ((limb) ((s32) in2[6])) * ((s32) in[6]) +
2 * (((limb) ((s32) in2[5])) * ((s32) in[7]) +
((limb) ((s32) in2[7])) * ((s32) in[5]) +
((limb) ((s32) in2[3])) * ((s32) in[9]) +
((limb) ((s32) in2[9])) * ((s32) in[3])) +
((limb) ((s32) in2[4])) * ((s32) in[8]) +
((limb) ((s32) in2[8])) * ((s32) in[4]);
output[13] = ((limb) ((s32) in2[6])) * ((s32) in[7]) +
((limb) ((s32) in2[7])) * ((s32) in[6]) +
((limb) ((s32) in2[5])) * ((s32) in[8]) +
((limb) ((s32) in2[8])) * ((s32) in[5]) +
((limb) ((s32) in2[4])) * ((s32) in[9]) +
((limb) ((s32) in2[9])) * ((s32) in[4]);
output[14] = 2 * (((limb) ((s32) in2[7])) * ((s32) in[7]) +
((limb) ((s32) in2[5])) * ((s32) in[9]) +
((limb) ((s32) in2[9])) * ((s32) in[5])) +
((limb) ((s32) in2[6])) * ((s32) in[8]) +
((limb) ((s32) in2[8])) * ((s32) in[6]);
output[15] = ((limb) ((s32) in2[7])) * ((s32) in[8]) +
((limb) ((s32) in2[8])) * ((s32) in[7]) +
((limb) ((s32) in2[6])) * ((s32) in[9]) +
((limb) ((s32) in2[9])) * ((s32) in[6]);
output[16] = ((limb) ((s32) in2[8])) * ((s32) in[8]) +
2 * (((limb) ((s32) in2[7])) * ((s32) in[9]) +
((limb) ((s32) in2[9])) * ((s32) in[7]));
output[17] = ((limb) ((s32) in2[8])) * ((s32) in[9]) +
((limb) ((s32) in2[9])) * ((s32) in[8]);
output[18] = 2 * ((limb) ((s32) in2[9])) * ((s32) in[9]);
}
/* Reduce a long form to a short form by taking the input mod 2^255 - 19.
*
* On entry: |output[i]| < 14*2^54
* On exit: |output[0..8]| < 280*2^54 */
static void freduce_degree(limb *output) {
/* Each of these shifts and adds ends up multiplying the value by 19.
*
* For output[0..8], the absolute entry value is < 14*2^54 and we add, at
* most, 19*14*2^54 thus, on exit, |output[0..8]| < 280*2^54. */
output[8] += output[18] << 4;
output[8] += output[18] << 1;
output[8] += output[18];
output[7] += output[17] << 4;
output[7] += output[17] << 1;
output[7] += output[17];
output[6] += output[16] << 4;
output[6] += output[16] << 1;
output[6] += output[16];
output[5] += output[15] << 4;
output[5] += output[15] << 1;
output[5] += output[15];
output[4] += output[14] << 4;
output[4] += output[14] << 1;
output[4] += output[14];
output[3] += output[13] << 4;
output[3] += output[13] << 1;
output[3] += output[13];
output[2] += output[12] << 4;
output[2] += output[12] << 1;
output[2] += output[12];
output[1] += output[11] << 4;
output[1] += output[11] << 1;
output[1] += output[11];
output[0] += output[10] << 4;
output[0] += output[10] << 1;
output[0] += output[10];
}
#if (-1 & 3) != 3
#error "This code only works on a two's complement system"
#endif
/* return v / 2^26, using only shifts and adds.
*
* On entry: v can take any value. */
static inline limb
div_by_2_26(const limb v)
{
/* High word of v; no shift needed. */
const uint32_t highword = (uint32_t) (((uint64_t) v) >> 32);
/* Set to all 1s if v was negative; else set to 0s. */
const int32_t sign = ((int32_t) highword) >> 31;
/* Set to 0x3ffffff if v was negative; else set to 0. */
const int32_t roundoff = ((uint32_t) sign) >> 6;
/* Should return v / (1<<26) */
return (v + roundoff) >> 26;
}
/* return v / (2^25), using only shifts and adds.
*
* On entry: v can take any value. */
static inline limb
div_by_2_25(const limb v)
{
/* High word of v; no shift needed*/
const uint32_t highword = (uint32_t) (((uint64_t) v) >> 32);
/* Set to all 1s if v was negative; else set to 0s. */
const int32_t sign = ((int32_t) highword) >> 31;
/* Set to 0x1ffffff if v was negative; else set to 0. */
const int32_t roundoff = ((uint32_t) sign) >> 7;
/* Should return v / (1<<25) */
return (v + roundoff) >> 25;
}
/* Reduce all coefficients of the short form input so that |x| < 2^26.
*
* On entry: |output[i]| < 280*2^54 */
static void freduce_coefficients(limb *output) {
unsigned i;
output[10] = 0;
for (i = 0; i < 10; i += 2) {
limb over = div_by_2_26(output[i]);
/* The entry condition (that |output[i]| < 280*2^54) means that over is, at
* most, 280*2^28 in the first iteration of this loop. This is added to the
* next limb and we can approximate the resulting bound of that limb by
* 281*2^54. */
output[i] -= over << 26;
output[i+1] += over;
/* For the first iteration, |output[i+1]| < 281*2^54, thus |over| <
* 281*2^29. When this is added to the next limb, the resulting bound can
* be approximated as 281*2^54.
*
* For subsequent iterations of the loop, 281*2^54 remains a conservative
* bound and no overflow occurs. */
over = div_by_2_25(output[i+1]);
output[i+1] -= over << 25;
output[i+2] += over;
}
/* Now |output[10]| < 281*2^29 and all other coefficients are reduced. */
output[0] += output[10] << 4;
output[0] += output[10] << 1;
output[0] += output[10];
output[10] = 0;
/* Now output[1..9] are reduced, and |output[0]| < 2^26 + 19*281*2^29
* So |over| will be no more than 2^16. */
{
limb over = div_by_2_26(output[0]);
output[0] -= over << 26;
output[1] += over;
}
/* Now output[0,2..9] are reduced, and |output[1]| < 2^25 + 2^16 < 2^26. The
* bound on |output[1]| is sufficient to meet our needs. */
}
/* A helpful wrapper around fproduct: output = in * in2.
*
* On entry: |in[i]| < 2^27 and |in2[i]| < 2^27.
*
* output must be distinct to both inputs. The output is reduced degree
* (indeed, one need only provide storage for 10 limbs) and |output[i]| < 2^26. */
static void
fmul(limb *output, const limb *in, const limb *in2) {
limb t[19];
fproduct(t, in, in2);
/* |t[i]| < 14*2^54 */
freduce_degree(t);
freduce_coefficients(t);
/* |t[i]| < 2^26 */
memcpy(output, t, sizeof(limb) * 10);
}
/* Square a number: output = in**2
*
* output must be distinct from the input. The inputs are reduced coefficient
* form, the output is not.
*
* output[x] <= 14 * the largest product of the input limbs. */
static void fsquare_inner(limb *output, const limb *in) {
output[0] = ((limb) ((s32) in[0])) * ((s32) in[0]);
output[1] = 2 * ((limb) ((s32) in[0])) * ((s32) in[1]);
output[2] = 2 * (((limb) ((s32) in[1])) * ((s32) in[1]) +
((limb) ((s32) in[0])) * ((s32) in[2]));
output[3] = 2 * (((limb) ((s32) in[1])) * ((s32) in[2]) +
((limb) ((s32) in[0])) * ((s32) in[3]));
output[4] = ((limb) ((s32) in[2])) * ((s32) in[2]) +
4 * ((limb) ((s32) in[1])) * ((s32) in[3]) +
2 * ((limb) ((s32) in[0])) * ((s32) in[4]);
output[5] = 2 * (((limb) ((s32) in[2])) * ((s32) in[3]) +
((limb) ((s32) in[1])) * ((s32) in[4]) +
((limb) ((s32) in[0])) * ((s32) in[5]));
output[6] = 2 * (((limb) ((s32) in[3])) * ((s32) in[3]) +
((limb) ((s32) in[2])) * ((s32) in[4]) +
((limb) ((s32) in[0])) * ((s32) in[6]) +
2 * ((limb) ((s32) in[1])) * ((s32) in[5]));
output[7] = 2 * (((limb) ((s32) in[3])) * ((s32) in[4]) +
((limb) ((s32) in[2])) * ((s32) in[5]) +
((limb) ((s32) in[1])) * ((s32) in[6]) +
((limb) ((s32) in[0])) * ((s32) in[7]));
output[8] = ((limb) ((s32) in[4])) * ((s32) in[4]) +
2 * (((limb) ((s32) in[2])) * ((s32) in[6]) +
((limb) ((s32) in[0])) * ((s32) in[8]) +
2 * (((limb) ((s32) in[1])) * ((s32) in[7]) +
((limb) ((s32) in[3])) * ((s32) in[5])));
output[9] = 2 * (((limb) ((s32) in[4])) * ((s32) in[5]) +
((limb) ((s32) in[3])) * ((s32) in[6]) +
((limb) ((s32) in[2])) * ((s32) in[7]) +
((limb) ((s32) in[1])) * ((s32) in[8]) +
((limb) ((s32) in[0])) * ((s32) in[9]));
output[10] = 2 * (((limb) ((s32) in[5])) * ((s32) in[5]) +
((limb) ((s32) in[4])) * ((s32) in[6]) +
((limb) ((s32) in[2])) * ((s32) in[8]) +
2 * (((limb) ((s32) in[3])) * ((s32) in[7]) +
((limb) ((s32) in[1])) * ((s32) in[9])));
output[11] = 2 * (((limb) ((s32) in[5])) * ((s32) in[6]) +
((limb) ((s32) in[4])) * ((s32) in[7]) +
((limb) ((s32) in[3])) * ((s32) in[8]) +
((limb) ((s32) in[2])) * ((s32) in[9]));
output[12] = ((limb) ((s32) in[6])) * ((s32) in[6]) +
2 * (((limb) ((s32) in[4])) * ((s32) in[8]) +
2 * (((limb) ((s32) in[5])) * ((s32) in[7]) +
((limb) ((s32) in[3])) * ((s32) in[9])));
output[13] = 2 * (((limb) ((s32) in[6])) * ((s32) in[7]) +
((limb) ((s32) in[5])) * ((s32) in[8]) +
((limb) ((s32) in[4])) * ((s32) in[9]));
output[14] = 2 * (((limb) ((s32) in[7])) * ((s32) in[7]) +
((limb) ((s32) in[6])) * ((s32) in[8]) +
2 * ((limb) ((s32) in[5])) * ((s32) in[9]));
output[15] = 2 * (((limb) ((s32) in[7])) * ((s32) in[8]) +
((limb) ((s32) in[6])) * ((s32) in[9]));
output[16] = ((limb) ((s32) in[8])) * ((s32) in[8]) +
4 * ((limb) ((s32) in[7])) * ((s32) in[9]);
output[17] = 2 * ((limb) ((s32) in[8])) * ((s32) in[9]);
output[18] = 2 * ((limb) ((s32) in[9])) * ((s32) in[9]);
}
/* fsquare sets output = in^2.
*
* On entry: The |in| argument is in reduced coefficients form and |in[i]| <
* 2^27.
*
* On exit: The |output| argument is in reduced coefficients form (indeed, one
* need only provide storage for 10 limbs) and |out[i]| < 2^26. */
static void
fsquare(limb *output, const limb *in) {
limb t[19];
fsquare_inner(t, in);
/* |t[i]| < 14*2^54 because the largest product of two limbs will be <
* 2^(27+27) and fsquare_inner adds together, at most, 14 of those
* products. */
freduce_degree(t);
freduce_coefficients(t);
/* |t[i]| < 2^26 */
memcpy(output, t, sizeof(limb) * 10);
}
/* Take a little-endian, 32-byte number and expand it into polynomial form */
static void
fexpand(limb *output, const u8 *input) {
#define F(n,start,shift,mask) \
output[n] = ((((limb) input[start + 0]) | \
((limb) input[start + 1]) << 8 | \
((limb) input[start + 2]) << 16 | \
((limb) input[start + 3]) << 24) >> shift) & mask;
F(0, 0, 0, 0x3ffffff);
F(1, 3, 2, 0x1ffffff);
F(2, 6, 3, 0x3ffffff);
F(3, 9, 5, 0x1ffffff);
F(4, 12, 6, 0x3ffffff);
F(5, 16, 0, 0x1ffffff);
F(6, 19, 1, 0x3ffffff);
F(7, 22, 3, 0x1ffffff);
F(8, 25, 4, 0x3ffffff);
F(9, 28, 6, 0x1ffffff);
#undef F
}
#if (-32 >> 1) != -16
#error "This code only works when >> does sign-extension on negative numbers"
#endif
/* s32_eq returns 0xffffffff iff a == b and zero otherwise. */
static s32 s32_eq(s32 a, s32 b) {
a = ~(a ^ b);
a &= a << 16;
a &= a << 8;
a &= a << 4;
a &= a << 2;
a &= a << 1;
return a >> 31;
}
/* s32_gte returns 0xffffffff if a >= b and zero otherwise, where a and b are
* both non-negative. */
static s32 s32_gte(s32 a, s32 b) {
a -= b;
/* a >= 0 iff a >= b. */
return ~(a >> 31);
}
/* Take a fully reduced polynomial form number and contract it into a
* little-endian, 32-byte array.
*
* On entry: |input_limbs[i]| < 2^26 */
static void
fcontract(u8 *output, limb *input_limbs) {
int i;
int j;
s32 input[10];
s32 mask;
/* |input_limbs[i]| < 2^26, so it's valid to convert to an s32. */
for (i = 0; i < 10; i++) {
input[i] = input_limbs[i];
}
for (j = 0; j < 2; ++j) {
for (i = 0; i < 9; ++i) {
if ((i & 1) == 1) {
/* This calculation is a time-invariant way to make input[i]
* non-negative by borrowing from the next-larger limb. */
const s32 mask = input[i] >> 31;
const s32 carry = -((input[i] & mask) >> 25);
input[i] = input[i] + (carry << 25);
input[i+1] = input[i+1] - carry;
} else {
const s32 mask = input[i] >> 31;
const s32 carry = -((input[i] & mask) >> 26);
input[i] = input[i] + (carry << 26);
input[i+1] = input[i+1] - carry;
}
}
/* There's no greater limb for input[9] to borrow from, but we can multiply
* by 19 and borrow from input[0], which is valid mod 2^255-19. */
{
const s32 mask = input[9] >> 31;
const s32 carry = -((input[9] & mask) >> 25);
input[9] = input[9] + (carry << 25);
input[0] = input[0] - (carry * 19);
}
/* After the first iteration, input[1..9] are non-negative and fit within
* 25 or 26 bits, depending on position. However, input[0] may be
* negative. */
}
/* The first borrow-propagation pass above ended with every limb
except (possibly) input[0] non-negative.
If input[0] was negative after the first pass, then it was because of a
carry from input[9]. On entry, input[9] < 2^26 so the carry was, at most,
one, since (2**26-1) >> 25 = 1. Thus input[0] >= -19.
In the second pass, each limb is decreased by at most one. Thus the second
borrow-propagation pass could only have wrapped around to decrease
input[0] again if the first pass left input[0] negative *and* input[1]
through input[9] were all zero. In that case, input[1] is now 2^25 - 1,
and this last borrow-propagation step will leave input[1] non-negative. */
{
const s32 mask = input[0] >> 31;
const s32 carry = -((input[0] & mask) >> 26);
input[0] = input[0] + (carry << 26);
input[1] = input[1] - carry;
}
/* All input[i] are now non-negative. However, there might be values between
* 2^25 and 2^26 in a limb which is, nominally, 25 bits wide. */
for (j = 0; j < 2; j++) {
for (i = 0; i < 9; i++) {
if ((i & 1) == 1) {
const s32 carry = input[i] >> 25;
input[i] &= 0x1ffffff;
input[i+1] += carry;
} else {
const s32 carry = input[i] >> 26;
input[i] &= 0x3ffffff;
input[i+1] += carry;
}
}
{
const s32 carry = input[9] >> 25;
input[9] &= 0x1ffffff;
input[0] += 19*carry;
}
}
/* If the first carry-chain pass, just above, ended up with a carry from
* input[9], and that caused input[0] to be out-of-bounds, then input[0] was
* < 2^26 + 2*19, because the carry was, at most, two.
*
* If the second pass carried from input[9] again then input[0] is < 2*19 and
* the input[9] -> input[0] carry didn't push input[0] out of bounds. */
/* It still remains the case that input might be between 2^255-19 and 2^255.
* In this case, input[1..9] must take their maximum value and input[0] must
* be >= (2^255-19) & 0x3ffffff, which is 0x3ffffed. */
mask = s32_gte(input[0], 0x3ffffed);
for (i = 1; i < 10; i++) {
if ((i & 1) == 1) {
mask &= s32_eq(input[i], 0x1ffffff);
} else {
mask &= s32_eq(input[i], 0x3ffffff);
}
}
/* mask is either 0xffffffff (if input >= 2^255-19) and zero otherwise. Thus
* this conditionally subtracts 2^255-19. */
input[0] -= mask & 0x3ffffed;
for (i = 1; i < 10; i++) {
if ((i & 1) == 1) {
input[i] -= mask & 0x1ffffff;
} else {
input[i] -= mask & 0x3ffffff;
}
}
input[1] <<= 2;
input[2] <<= 3;
input[3] <<= 5;
input[4] <<= 6;
input[6] <<= 1;
input[7] <<= 3;
input[8] <<= 4;
input[9] <<= 6;
#define F(i, s) \
output[s+0] |= input[i] & 0xff; \
output[s+1] = (input[i] >> 8) & 0xff; \
output[s+2] = (input[i] >> 16) & 0xff; \
output[s+3] = (input[i] >> 24) & 0xff;
output[0] = 0;
output[16] = 0;
F(0,0);
F(1,3);
F(2,6);
F(3,9);
F(4,12);
F(5,16);
F(6,19);
F(7,22);
F(8,25);
F(9,28);
#undef F
}
/* Input: Q, Q', Q-Q'
* Output: 2Q, Q+Q'
*
* x2 z3: long form
* x3 z3: long form
* x z: short form, destroyed
* xprime zprime: short form, destroyed
* qmqp: short form, preserved
*
* On entry and exit, the absolute value of the limbs of all inputs and outputs
* are < 2^26. */
static void fmonty(limb *x2, limb *z2, /* output 2Q */
limb *x3, limb *z3, /* output Q + Q' */
limb *x, limb *z, /* input Q */
limb *xprime, limb *zprime, /* input Q' */
const limb *qmqp /* input Q - Q' */) {
limb origx[10], origxprime[10], zzz[19], xx[19], zz[19], xxprime[19],
zzprime[19], zzzprime[19], xxxprime[19];
memcpy(origx, x, 10 * sizeof(limb));
fsum(x, z);
/* |x[i]| < 2^27 */
fdifference(z, origx); /* does x - z */
/* |z[i]| < 2^27 */
memcpy(origxprime, xprime, sizeof(limb) * 10);
fsum(xprime, zprime);
/* |xprime[i]| < 2^27 */
fdifference(zprime, origxprime);
/* |zprime[i]| < 2^27 */
fproduct(xxprime, xprime, z);
/* |xxprime[i]| < 14*2^54: the largest product of two limbs will be <
* 2^(27+27) and fproduct adds together, at most, 14 of those products.
* (Approximating that to 2^58 doesn't work out.) */
fproduct(zzprime, x, zprime);
/* |zzprime[i]| < 14*2^54 */
freduce_degree(xxprime);
freduce_coefficients(xxprime);
/* |xxprime[i]| < 2^26 */
freduce_degree(zzprime);
freduce_coefficients(zzprime);
/* |zzprime[i]| < 2^26 */
memcpy(origxprime, xxprime, sizeof(limb) * 10);
fsum(xxprime, zzprime);
/* |xxprime[i]| < 2^27 */
fdifference(zzprime, origxprime);
/* |zzprime[i]| < 2^27 */
fsquare(xxxprime, xxprime);
/* |xxxprime[i]| < 2^26 */
fsquare(zzzprime, zzprime);
/* |zzzprime[i]| < 2^26 */
fproduct(zzprime, zzzprime, qmqp);
/* |zzprime[i]| < 14*2^52 */
freduce_degree(zzprime);
freduce_coefficients(zzprime);
/* |zzprime[i]| < 2^26 */
memcpy(x3, xxxprime, sizeof(limb) * 10);
memcpy(z3, zzprime, sizeof(limb) * 10);
fsquare(xx, x);
/* |xx[i]| < 2^26 */
fsquare(zz, z);
/* |zz[i]| < 2^26 */
fproduct(x2, xx, zz);
/* |x2[i]| < 14*2^52 */
freduce_degree(x2);
freduce_coefficients(x2);
/* |x2[i]| < 2^26 */
fdifference(zz, xx); // does zz = xx - zz
/* |zz[i]| < 2^27 */
memset(zzz + 10, 0, sizeof(limb) * 9);
fscalar_product(zzz, zz, 121665);
/* |zzz[i]| < 2^(27+17) */
/* No need to call freduce_degree here:
fscalar_product doesn't increase the degree of its input. */
freduce_coefficients(zzz);
/* |zzz[i]| < 2^26 */
fsum(zzz, xx);
/* |zzz[i]| < 2^27 */
fproduct(z2, zz, zzz);
/* |z2[i]| < 14*2^(26+27) */
freduce_degree(z2);
freduce_coefficients(z2);
/* |z2|i| < 2^26 */
}
/* Conditionally swap two reduced-form limb arrays if 'iswap' is 1, but leave
* them unchanged if 'iswap' is 0. Runs in data-invariant time to avoid
* side-channel attacks.
*
* NOTE that this function requires that 'iswap' be 1 or 0; other values give
* wrong results. Also, the two limb arrays must be in reduced-coefficient,
* reduced-degree form: the values in a[10..19] or b[10..19] aren't swapped,
* and all all values in a[0..9],b[0..9] must have magnitude less than
* INT32_MAX. */
static void
swap_conditional(limb a[19], limb b[19], limb iswap) {
unsigned i;
const s32 swap = (s32) -iswap;
for (i = 0; i < 10; ++i) {
const s32 x = swap & ( ((s32)a[i]) ^ ((s32)b[i]) );
a[i] = ((s32)a[i]) ^ x;
b[i] = ((s32)b[i]) ^ x;
}
}
/* Calculates nQ where Q is the x-coordinate of a point on the curve
*
* resultx/resultz: the x coordinate of the resulting curve point (short form)
* n: a little endian, 32-byte number
* q: a point of the curve (short form) */
static void
cmult(limb *resultx, limb *resultz, const u8 *n, const limb *q) {
limb a[19] = {0}, b[19] = {1}, c[19] = {1}, d[19] = {0};
limb *nqpqx = a, *nqpqz = b, *nqx = c, *nqz = d, *t;
limb e[19] = {0}, f[19] = {1}, g[19] = {0}, h[19] = {1};
limb *nqpqx2 = e, *nqpqz2 = f, *nqx2 = g, *nqz2 = h;
unsigned i, j;
memcpy(nqpqx, q, sizeof(limb) * 10);
for (i = 0; i < 32; ++i) {
u8 byte = n[31 - i];
for (j = 0; j < 8; ++j) {
const limb bit = byte >> 7;
swap_conditional(nqx, nqpqx, bit);
swap_conditional(nqz, nqpqz, bit);
fmonty(nqx2, nqz2,
nqpqx2, nqpqz2,
nqx, nqz,
nqpqx, nqpqz,
q);
swap_conditional(nqx2, nqpqx2, bit);
swap_conditional(nqz2, nqpqz2, bit);
t = nqx;
nqx = nqx2;
nqx2 = t;
t = nqz;
nqz = nqz2;
nqz2 = t;
t = nqpqx;
nqpqx = nqpqx2;
nqpqx2 = t;
t = nqpqz;
nqpqz = nqpqz2;
nqpqz2 = t;
byte <<= 1;
}
}
memcpy(resultx, nqx, sizeof(limb) * 10);
memcpy(resultz, nqz, sizeof(limb) * 10);
}
// -----------------------------------------------------------------------------
// Shamelessly copied from djb's code
// -----------------------------------------------------------------------------
static void
crecip(limb *out, const limb *z) {
limb z2[10];
limb z9[10];
limb z11[10];
limb z2_5_0[10];
limb z2_10_0[10];
limb z2_20_0[10];
limb z2_50_0[10];
limb z2_100_0[10];
limb t0[10];
limb t1[10];
int i;
/* 2 */ fsquare(z2,z);
/* 4 */ fsquare(t1,z2);
/* 8 */ fsquare(t0,t1);
/* 9 */ fmul(z9,t0,z);
/* 11 */ fmul(z11,z9,z2);
/* 22 */ fsquare(t0,z11);
/* 2^5 - 2^0 = 31 */ fmul(z2_5_0,t0,z9);
/* 2^6 - 2^1 */ fsquare(t0,z2_5_0);
/* 2^7 - 2^2 */ fsquare(t1,t0);
/* 2^8 - 2^3 */ fsquare(t0,t1);
/* 2^9 - 2^4 */ fsquare(t1,t0);
/* 2^10 - 2^5 */ fsquare(t0,t1);
/* 2^10 - 2^0 */ fmul(z2_10_0,t0,z2_5_0);
/* 2^11 - 2^1 */ fsquare(t0,z2_10_0);
/* 2^12 - 2^2 */ fsquare(t1,t0);
/* 2^20 - 2^10 */ for (i = 2;i < 10;i += 2) { fsquare(t0,t1); fsquare(t1,t0); }
/* 2^20 - 2^0 */ fmul(z2_20_0,t1,z2_10_0);
/* 2^21 - 2^1 */ fsquare(t0,z2_20_0);
/* 2^22 - 2^2 */ fsquare(t1,t0);
/* 2^40 - 2^20 */ for (i = 2;i < 20;i += 2) { fsquare(t0,t1); fsquare(t1,t0); }
/* 2^40 - 2^0 */ fmul(t0,t1,z2_20_0);
/* 2^41 - 2^1 */ fsquare(t1,t0);
/* 2^42 - 2^2 */ fsquare(t0,t1);
/* 2^50 - 2^10 */ for (i = 2;i < 10;i += 2) { fsquare(t1,t0); fsquare(t0,t1); }
/* 2^50 - 2^0 */ fmul(z2_50_0,t0,z2_10_0);
/* 2^51 - 2^1 */ fsquare(t0,z2_50_0);
/* 2^52 - 2^2 */ fsquare(t1,t0);
/* 2^100 - 2^50 */ for (i = 2;i < 50;i += 2) { fsquare(t0,t1); fsquare(t1,t0); }
/* 2^100 - 2^0 */ fmul(z2_100_0,t1,z2_50_0);
/* 2^101 - 2^1 */ fsquare(t1,z2_100_0);
/* 2^102 - 2^2 */ fsquare(t0,t1);
/* 2^200 - 2^100 */ for (i = 2;i < 100;i += 2) { fsquare(t1,t0); fsquare(t0,t1); }
/* 2^200 - 2^0 */ fmul(t1,t0,z2_100_0);
/* 2^201 - 2^1 */ fsquare(t0,t1);
/* 2^202 - 2^2 */ fsquare(t1,t0);
/* 2^250 - 2^50 */ for (i = 2;i < 50;i += 2) { fsquare(t0,t1); fsquare(t1,t0); }
/* 2^250 - 2^0 */ fmul(t0,t1,z2_50_0);
/* 2^251 - 2^1 */ fsquare(t1,t0);
/* 2^252 - 2^2 */ fsquare(t0,t1);
/* 2^253 - 2^3 */ fsquare(t1,t0);
/* 2^254 - 2^4 */ fsquare(t0,t1);
/* 2^255 - 2^5 */ fsquare(t1,t0);
/* 2^255 - 21 */ fmul(out,t1,z11);
}
int
curve25519_donna(u8 *mypublic, const u8 *secret, const u8 *basepoint) {
limb bp[10], x[10], z[11], zmone[10];
uint8_t e[32];
int i;
for (i = 0; i < 32; ++i) e[i] = secret[i];
e[0] &= 248;
e[31] &= 127;
e[31] |= 64;
fexpand(bp, basepoint);
cmult(x, z, e, bp);
crecip(zmone, z);
fmul(z, x, zmone);
fcontract(mypublic, z);
return 0;
}

13
curve25519-donna.podspec Normal file
View file

@ -0,0 +1,13 @@
Pod::Spec.new do |s|
s.name = "curve25519-donna"
s.version = "1.2.1"
s.summary = "Implementations of a fast elliptic-curve, Diffie-Hellman primitive"
s.description = <<-DESC
Curve25519 is a state-of-the-art Diffie-Hellman function suitable for a wide variety of applications.
DESC
s.homepage = "http://code.google.com/p/curve25519-donna"
s.license = 'BSD 3-Clause'
s.author = 'Dan Bernstein'
s.source = { :git => "https://github.com/agl/curve25519-donna.git", :tag => "1.2.1" }
s.source_files = 'curve25519-donna.c'
end

4
python-src/curve25519/__init__.py Normal file
View file

@ -0,0 +1,4 @@
from .keys import Private, Public
hush_pyflakes = [Private, Public]; del hush_pyflakes

105
python-src/curve25519/curve25519module.c Normal file
View file

@ -0,0 +1,105 @@
/* tell python that PyArg_ParseTuple(t#) means Py_ssize_t, not int */
#define PY_SSIZE_T_CLEAN
#include <Python.h>
#if (PY_VERSION_HEX < 0x02050000)
typedef int Py_ssize_t;
#endif
/* This is required for compatibility with Python 2. */
#if PY_MAJOR_VERSION >= 3
#include <bytesobject.h>
#define y "y"
#else
#define PyBytes_FromStringAndSize PyString_FromStringAndSize
#define y "t"
#endif
int curve25519_donna(char *mypublic,
const char *secret, const char *basepoint);
static PyObject *
pycurve25519_makeprivate(PyObject *self, PyObject *args)
{
char *in1;
Py_ssize_t in1len;
if (!PyArg_ParseTuple(args, y"#:clamp", &in1, &in1len))
return NULL;
if (in1len != 32) {
PyErr_SetString(PyExc_ValueError, "input must be 32-byte string");
return NULL;
}
in1[0] &= 248;
in1[31] &= 127;
in1[31] |= 64;
return PyBytes_FromStringAndSize((char *)in1, 32);
}
static PyObject *
pycurve25519_makepublic(PyObject *self, PyObject *args)
{
const char *private;
char mypublic[32];
char basepoint[32] = {9};
Py_ssize_t privatelen;
if (!PyArg_ParseTuple(args, y"#:makepublic", &private, &privatelen))
return NULL;
if (privatelen != 32) {
PyErr_SetString(PyExc_ValueError, "input must be 32-byte string");
return NULL;
}
curve25519_donna(mypublic, private, basepoint);
return PyBytes_FromStringAndSize((char *)mypublic, 32);
}
static PyObject *
pycurve25519_makeshared(PyObject *self, PyObject *args)
{
const char *myprivate, *theirpublic;
char shared_key[32];
Py_ssize_t myprivatelen, theirpubliclen;
if (!PyArg_ParseTuple(args, y"#"y"#:generate",
&myprivate, &myprivatelen, &theirpublic, &theirpubliclen))
return NULL;
if (myprivatelen != 32) {
PyErr_SetString(PyExc_ValueError, "input must be 32-byte string");
return NULL;
}
if (theirpubliclen != 32) {
PyErr_SetString(PyExc_ValueError, "input must be 32-byte string");
return NULL;
}
curve25519_donna(shared_key, myprivate, theirpublic);
return PyBytes_FromStringAndSize((char *)shared_key, 32);
}
static PyMethodDef
curve25519_functions[] = {
{"make_private", pycurve25519_makeprivate, METH_VARARGS, "data->private"},
{"make_public", pycurve25519_makepublic, METH_VARARGS, "private->public"},
{"make_shared", pycurve25519_makeshared, METH_VARARGS, "private+public->shared"},
{NULL, NULL, 0, NULL},
};
#if PY_MAJOR_VERSION >= 3
static struct PyModuleDef
curve25519_module = {
PyModuleDef_HEAD_INIT,
"_curve25519",
NULL,
NULL,
curve25519_functions,
};
PyObject *
PyInit__curve25519(void)
{
return PyModule_Create(&curve25519_module);
}
#else
PyMODINIT_FUNC
init_curve25519(void)
{
(void)Py_InitModule("_curve25519", curve25519_functions);
}
#endif

46
python-src/curve25519/keys.py Normal file
View file

@ -0,0 +1,46 @@
from . import _curve25519
from hashlib import sha256
import os
# the curve25519 functions are really simple, and could be used without an
# OOP layer, but it's a bit too easy to accidentally swap the private and
# public keys that way.
def _hash_shared(shared):
return sha256(b"curve25519-shared:"+shared).digest()
class Private:
def __init__(self, secret=None, seed=None):
if secret is None:
if seed is None:
secret = os.urandom(32)
else:
secret = sha256(b"curve25519-private:"+seed).digest()
else:
assert seed is None, "provide secret, seed, or neither, not both"
if not isinstance(secret, bytes) or len(secret) != 32:
raise TypeError("secret= must be 32-byte string")
self.private = _curve25519.make_private(secret)
def serialize(self):
return self.private
def get_public(self):
return Public(_curve25519.make_public(self.private))
def get_shared_key(self, public, hashfunc=None):
if not isinstance(public, Public):
raise ValueError("'public' must be an instance of Public")
if hashfunc is None:
hashfunc = _hash_shared
shared = _curve25519.make_shared(self.private, public.public)
return hashfunc(shared)
class Public:
def __init__(self, public):
assert isinstance(public, bytes)
assert len(public) == 32
self.public = public
def serialize(self):
return self.public

0
python-src/curve25519/test/__init__.py Normal file
View file

99
python-src/curve25519/test/test_curve25519.py Executable file
View file

@ -0,0 +1,99 @@
#! /usr/bin/python
import unittest
from curve25519 import Private, Public
from hashlib import sha1, sha256
from binascii import hexlify
class Basic(unittest.TestCase):
def test_basic(self):
secret1 = b"abcdefghijklmnopqrstuvwxyz123456"
self.assertEqual(len(secret1), 32)
secret2 = b"654321zyxwvutsrqponmlkjihgfedcba"
self.assertEqual(len(secret2), 32)
priv1 = Private(secret=secret1)
pub1 = priv1.get_public()
priv2 = Private(secret=secret2)
pub2 = priv2.get_public()
shared12 = priv1.get_shared_key(pub2)
e = b"b0818125eab42a8ac1af5e8b9b9c15ed2605c2bbe9675de89e5e6e7f442b9598"
self.assertEqual(hexlify(shared12), e)
shared21 = priv2.get_shared_key(pub1)
self.assertEqual(shared12, shared21)
pub2a = Public(pub2.serialize())
shared12a = priv1.get_shared_key(pub2a)
self.assertEqual(hexlify(shared12a), e)
def test_errors(self):
priv1 = Private()
self.assertRaises(ValueError, priv1.get_shared_key, priv1)
def test_seed(self):
# use 32-byte secret
self.assertRaises(TypeError, Private, secret=123)
self.assertRaises(TypeError, Private, secret=b"too short")
secret1 = b"abcdefghijklmnopqrstuvwxyz123456"
assert len(secret1) == 32
priv1 = Private(secret=secret1)
priv1a = Private(secret=secret1)
priv1b = Private(priv1.serialize())
self.assertEqual(priv1.serialize(), priv1a.serialize())
self.assertEqual(priv1.serialize(), priv1b.serialize())
e = b"6062636465666768696a6b6c6d6e6f707172737475767778797a313233343576"
self.assertEqual(hexlify(priv1.serialize()), e)
# the private key is a clamped form of the secret, so they won't
# quite be the same
p = Private(secret=b"\x00"*32)
self.assertEqual(hexlify(p.serialize()), b"00"*31+b"40")
p = Private(secret=b"\xff"*32)
self.assertEqual(hexlify(p.serialize()), b"f8"+b"ff"*30+b"7f")
# use arbitrary-length seed
self.assertRaises(TypeError, Private, seed=123)
priv1 = Private(seed=b"abc")
priv1a = Private(seed=b"abc")
priv1b = Private(priv1.serialize())
self.assertEqual(priv1.serialize(), priv1a.serialize())
self.assertEqual(priv1.serialize(), priv1b.serialize())
self.assertRaises(AssertionError, Private, seed=b"abc", secret=b"no")
priv1 = Private(seed=b"abc")
priv1a = Private(priv1.serialize())
self.assertEqual(priv1.serialize(), priv1a.serialize())
self.assertRaises(AssertionError, Private, seed=b"abc", secret=b"no")
# use built-in os.urandom
priv2 = Private()
priv2a = Private(priv2.private)
self.assertEqual(priv2.serialize(), priv2a.serialize())
# attempt to use both secret= and seed=, not allowed
self.assertRaises(AssertionError, Private, seed=b"abc", secret=b"no")
def test_hashfunc(self):
priv1 = Private(seed=b"abc")
priv2 = Private(seed=b"def")
shared_sha256 = priv1.get_shared_key(priv2.get_public())
e = b"da959ffe77ebeb4757fe5ba310e28ede425ae0d0ff5ec9c884e2d08f311cf5e5"
self.assertEqual(hexlify(shared_sha256), e)
# confirm the hash function remains what we think it is
def myhash(shared_key):
return sha256(b"curve25519-shared:"+shared_key).digest()
shared_myhash = priv1.get_shared_key(priv2.get_public(), myhash)
self.assertEqual(hexlify(shared_myhash), e)
def hexhash(shared_key):
return sha1(shared_key).hexdigest().encode()
shared_hexhash = priv1.get_shared_key(priv2.get_public(), hexhash)
self.assertEqual(shared_hexhash,
b"80eec98222c8edc4324fb9477a3c775ce7c6c93a")
if __name__ == "__main__":
unittest.main()

46
python-src/curve25519/test/test_speed.py Executable file
View file

@ -0,0 +1,46 @@
#! /usr/bin/python
from time import time
from curve25519 import Private
count = 10000
elapsed_get_public = 0.0
elapsed_get_shared = 0.0
def abbreviate_time(data):
# 1.23s, 790ms, 132us
if data is None:
return ""
s = float(data)
if s >= 10:
#return abbreviate.abbreviate_time(data)
return "%d" % s
if s >= 1.0:
return "%.2fs" % s
if s >= 0.01:
return "%dms" % (1000*s)
if s >= 0.001:
return "%.1fms" % (1000*s)
if s >= 0.000001:
return "%.1fus" % (1000000*s)
return "%dns" % (1000000000*s)
def nohash(key): return key
for i in range(count):
p = Private()
start = time()
pub = p.get_public()
elapsed_get_public += time() - start
pub2 = Private().get_public()
start = time()
shared = p.get_shared_key(pub2) #, hashfunc=nohash)
elapsed_get_shared += time() - start
print("get_public: %s" % abbreviate_time(elapsed_get_public / count))
print("get_shared: %s" % abbreviate_time(elapsed_get_shared / count))
# these take about 560us-570us each (with the default compiler settings, -Os)
# on my laptop, same with -O2
# of which the python overhead is about 5us
# and the get_shared_key() hash step adds about 5us

38
setup.py Executable file
View file

@ -0,0 +1,38 @@
#! /usr/bin/python
from subprocess import Popen, PIPE
from distutils.core import setup, Extension
version = Popen(["git", "describe", "--tags"], stdout=PIPE).communicate()[0]\
.strip().decode("utf8")
ext_modules = [Extension("curve25519._curve25519",
["python-src/curve25519/curve25519module.c",
"curve25519-donna.c"],
)]
short_description="Python wrapper for the Curve25519 cryptographic library"
long_description="""\
Curve25519 is a fast elliptic-curve key-agreement protocol, in which two
parties Alice and Bob each generate a (public,private) keypair, exchange
public keys, and can then compute the same shared key. Specifically, Alice
computes F(Aprivate, Bpublic), Bob computes F(Bprivate, Apublic), and both
get the same value (and nobody else can guess that shared value, even if they
know Apublic and Bpublic).
This is a Python wrapper for the portable 'curve25519-donna' implementation
of this algorithm, written by Adam Langley, hosted at
http://code.google.com/p/curve25519-donna/
"""
setup(name="curve25519-donna",
version=version,
description=short_description,
long_description=long_description,
author="Brian Warner",
author_email="warner-pycurve25519-donna@lothar.com",
license="BSD",
packages=["curve25519", "curve25519.test"],
package_dir={"curve25519": "python-src/curve25519"},
ext_modules=ext_modules,
)

50
speed-curve25519.c Normal file
View file

@ -0,0 +1,50 @@
#include <stdio.h>
#include <string.h>
#include <sys/time.h>
#include <time.h>
#include <stdint.h>
typedef uint8_t u8;
extern void curve25519_donna(u8 *output, const u8 *secret, const u8 *bp);
static uint64_t
time_now() {
struct timeval tv;
uint64_t ret;
gettimeofday(&tv, NULL);
ret = tv.tv_sec;
ret *= 1000000;
ret += tv.tv_usec;
return ret;
}
int
main() {
static const unsigned char basepoint[32] = {9};
unsigned char mysecret[32], mypublic[32];
unsigned i;
uint64_t start, end;
memset(mysecret, 42, 32);
mysecret[0] &= 248;
mysecret[31] &= 127;
mysecret[31] |= 64;
// Load the caches
for (i = 0; i < 1000; ++i) {
curve25519_donna(mypublic, mysecret, basepoint);
}
start = time_now();
for (i = 0; i < 30000; ++i) {
curve25519_donna(mypublic, mysecret, basepoint);
}
end = time_now();
printf("%luus\n", (unsigned long) ((end - start) / 30000));
return 0;
}

54
test-curve25519.c Normal file
View file

@ -0,0 +1,54 @@
/*
test-curve25519 version 20050915
D. J. Bernstein
Public domain.
Tiny modifications by agl
*/
#include <stdio.h>
extern void curve25519_donna(unsigned char *output, const unsigned char *a,
const unsigned char *b);
void doit(unsigned char *ek,unsigned char *e,unsigned char *k);
void doit(unsigned char *ek,unsigned char *e,unsigned char *k)
{
int i;
for (i = 0;i < 32;++i) printf("%02x",(unsigned int) e[i]); printf(" ");
for (i = 0;i < 32;++i) printf("%02x",(unsigned int) k[i]); printf(" ");
curve25519_donna(ek,e,k);
for (i = 0;i < 32;++i) printf("%02x",(unsigned int) ek[i]); printf("\n");
}
unsigned char e1k[32];
unsigned char e2k[32];
unsigned char e1e2k[32];
unsigned char e2e1k[32];
unsigned char e1[32] = {3};
unsigned char e2[32] = {5};
unsigned char k[32] = {9};
int
main()
{
int loop;
int i;
for (loop = 0;loop < 10000;++loop) {
doit(e1k,e1,k);
doit(e2e1k,e2,e1k);
doit(e2k,e2,k);
doit(e1e2k,e1,e2k);
for (i = 0;i < 32;++i) if (e1e2k[i] != e2e1k[i]) {
printf("fail\n");
return 1;
}
for (i = 0;i < 32;++i) e1[i] ^= e2k[i];
for (i = 0;i < 32;++i) e2[i] ^= e1k[i];
for (i = 0;i < 32;++i) k[i] ^= e1e2k[i];
}
return 0;
}

39
test-noncanon.c Normal file
View file

@ -0,0 +1,39 @@
/* This file can be used to test whether the code handles non-canonical curve
* points (i.e. points with the 256th bit set) in the same way as the reference
* implementation. */
#include <stdint.h>
#include <stdio.h>
#include <string.h>
extern void curve25519_donna(unsigned char *output, const unsigned char *a,
const unsigned char *b);
int
main()
{
static const uint8_t point1[32] = {
0x25,0x00,0x00,0x00,0x00,0x00,0x00,0x00,
0x00,0x00,0x00,0x00,0x00,0x00,0x00,0x00,
0x00,0x00,0x00,0x00,0x00,0x00,0x00,0x00,
0x00,0x00,0x00,0x00,0x00,0x00,0x00,0x00,
};
static const uint8_t point2[32] = {
0xff,0xff,0xff,0xff,0xff,0xff,0xff,0xff,
0xff,0xff,0xff,0xff,0xff,0xff,0xff,0xff,
0xff,0xff,0xff,0xff,0xff,0xff,0xff,0xff,
0xff,0xff,0xff,0xff,0xff,0xff,0xff,0xff,
};
static const uint8_t scalar[32] = { 1 };
uint8_t out1[32], out2[32];
curve25519_donna(out1, scalar, point1);
curve25519_donna(out2, scalar, point2);
if (0 == memcmp(out1, out2, sizeof(out1))) {
fprintf(stderr, "Top bit not ignored.\n");
return 1;
}
fprintf(stderr, "Top bit correctly ignored.\n");
return 0;
}

72
test-sc-curve25519.c Normal file
View file

@ -0,0 +1,72 @@
#define _GNU_SOURCE
#include <stdio.h>
#include <string.h>
#include <stdint.h>
#include <math.h>
extern void curve25519_donna(uint8_t *, const uint8_t *, const uint8_t *);
extern uint64_t tsc_read();
int
main(int argc, char **argv) {
uint8_t private_key[32], public[32], peer1[32], peer2[32], output[32];
static const uint8_t basepoint[32] = {9};
unsigned i;
uint64_t sum = 0, sum_squares = 0, skipped = 0, mean;
static const unsigned count = 200000;
memset(private_key, 42, sizeof(private_key));
private_key[0] &= 248;
private_key[31] &= 127;
private_key[31] |= 64;
curve25519_donna(public, private_key, basepoint);
memset(peer1, 0, sizeof(peer1));
memset(peer2, 255, sizeof(peer2));
for (i = 0; i < count; ++i) {
const uint64_t start = tsc_read();
curve25519_donna(output, peer1, public);
const uint64_t end = tsc_read();
const uint64_t delta = end - start;
if (delta > 650000) {
// something terrible happened (task switch etc)
skipped++;
continue;
}
sum += delta;
sum_squares += (delta * delta);
}
mean = sum / ((uint64_t) count);
printf("all 0: mean:%lu sd:%f skipped:%lu\n",
mean,
sqrt((double)(sum_squares/((uint64_t) count) - mean*mean)),
skipped);
sum = sum_squares = skipped = 0;
for (i = 0; i < count; ++i) {
const uint64_t start = tsc_read();
curve25519_donna(output, peer2, public);
const uint64_t end = tsc_read();
const uint64_t delta = end - start;
if (delta > 650000) {
// something terrible happened (task switch etc)
skipped++;
continue;
}
sum += delta;
sum_squares += (delta * delta);
}
mean = sum / ((uint64_t) count);
printf("all 1: mean:%lu sd:%f skipped:%lu\n",
mean,
sqrt((double)(sum_squares/((uint64_t) count) - mean*mean)),
skipped);
return 0;
}

8
test-sc-curve25519.s Normal file
View file

@ -0,0 +1,8 @@
.text
.globl tsc_read
tsc_read:
rdtsc
shl $32,%rdx
or %rdx,%rax
ret