1 /*
2 version 20081011
3 Matthew Dempsky
4 Public domain.
5 Derived from public domain code by D. J. Bernstein.
6 */
7 #include "curve25519.h"
8
add(unsigned int out[32],const unsigned int a[32],const unsigned int b[32])9 static void add(unsigned int out[32],const unsigned int a[32],const unsigned int b[32])
10 {
11 unsigned int j;
12 unsigned int u;
13 u = 0;
14 for (j = 0;j < 31;++j) { u += a[j] + b[j]; out[j] = u & 255; u >>= 8; }
15 u += a[31] + b[31]; out[31] = u;
16 }
17
sub(unsigned int out[32],const unsigned int a[32],const unsigned int b[32])18 static void sub(unsigned int out[32],const unsigned int a[32],const unsigned int b[32])
19 {
20 unsigned int j;
21 unsigned int u;
22 u = 218;
23 for (j = 0;j < 31;++j) {
24 u += a[j] + 65280 - b[j];
25 out[j] = u & 255;
26 u >>= 8;
27 }
28 u += a[31] - b[31];
29 out[31] = u;
30 }
31
squeeze(unsigned int a[32])32 static void squeeze(unsigned int a[32])
33 {
34 unsigned int j;
35 unsigned int u;
36 u = 0;
37 for (j = 0;j < 31;++j) { u += a[j]; a[j] = u & 255; u >>= 8; }
38 u += a[31]; a[31] = u & 127;
39 u = 19 * (u >> 7);
40 for (j = 0;j < 31;++j) { u += a[j]; a[j] = u & 255; u >>= 8; }
41 u += a[31]; a[31] = u;
42 }
43
44 static const unsigned int minusp[32] = {
45 19, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 128
46 } ;
47
freeze(unsigned int a[32])48 static void freeze(unsigned int a[32])
49 {
50 unsigned int aorig[32];
51 unsigned int j;
52 unsigned int negative;
53
54 for (j = 0;j < 32;++j) aorig[j] = a[j];
55 add(a,a,minusp);
56 negative = -((a[31] >> 7) & 1);
57 for (j = 0;j < 32;++j) a[j] ^= negative & (aorig[j] ^ a[j]);
58 }
59
mult(unsigned int out[32],const unsigned int a[32],const unsigned int b[32])60 static void mult(unsigned int out[32],const unsigned int a[32],const unsigned int b[32])
61 {
62 unsigned int i;
63 unsigned int j;
64 unsigned int u;
65
66 for (i = 0;i < 32;++i) {
67 u = 0;
68 for (j = 0;j <= i;++j) u += a[j] * b[i - j];
69 for (j = i + 1;j < 32;++j) u += 38 * a[j] * b[i + 32 - j];
70 out[i] = u;
71 }
72 squeeze(out);
73 }
74
mult121665(unsigned int out[32],const unsigned int a[32])75 static void mult121665(unsigned int out[32],const unsigned int a[32])
76 {
77 unsigned int j;
78 unsigned int u;
79
80 u = 0;
81 for (j = 0;j < 31;++j) { u += 121665 * a[j]; out[j] = u & 255; u >>= 8; }
82 u += 121665 * a[31]; out[31] = u & 127;
83 u = 19 * (u >> 7);
84 for (j = 0;j < 31;++j) { u += out[j]; out[j] = u & 255; u >>= 8; }
85 u += out[j]; out[j] = u;
86 }
87
square(unsigned int out[32],const unsigned int a[32])88 static void square(unsigned int out[32],const unsigned int a[32])
89 {
90 unsigned int i;
91 unsigned int j;
92 unsigned int u;
93
94 for (i = 0;i < 32;++i) {
95 u = 0;
96 for (j = 0;j < i - j;++j) u += a[j] * a[i - j];
97 for (j = i + 1;j < i + 32 - j;++j) u += 38 * a[j] * a[i + 32 - j];
98 u *= 2;
99 if ((i & 1) == 0) {
100 u += a[i / 2] * a[i / 2];
101 u += 38 * a[i / 2 + 16] * a[i / 2 + 16];
102 }
103 out[i] = u;
104 }
105 squeeze(out);
106 }
107
select(unsigned int p[64],unsigned int q[64],const unsigned int r[64],const unsigned int s[64],unsigned int b)108 static void select(unsigned int p[64],unsigned int q[64],const unsigned int r[64],const unsigned int s[64],unsigned int b)
109 {
110 unsigned int j;
111 unsigned int t;
112 unsigned int bminus1;
113
114 bminus1 = b - 1;
115 for (j = 0;j < 64;++j) {
116 t = bminus1 & (r[j] ^ s[j]);
117 p[j] = s[j] ^ t;
118 q[j] = r[j] ^ t;
119 }
120 }
121
mainloop(unsigned int work[64],const unsigned char e[32])122 static void mainloop(unsigned int work[64],const unsigned char e[32])
123 {
124 unsigned int xzm1[64];
125 unsigned int xzm[64];
126 unsigned int xzmb[64];
127 unsigned int xzm1b[64];
128 unsigned int xznb[64];
129 unsigned int xzn1b[64];
130 unsigned int a0[64];
131 unsigned int a1[64];
132 unsigned int b0[64];
133 unsigned int b1[64];
134 unsigned int c1[64];
135 unsigned int r[32];
136 unsigned int s[32];
137 unsigned int t[32];
138 unsigned int u[32];
139 unsigned int j;
140 unsigned int b;
141 int pos;
142
143 for (j = 0;j < 32;++j) xzm1[j] = work[j];
144 xzm1[32] = 1;
145 for (j = 33;j < 64;++j) xzm1[j] = 0;
146
147 xzm[0] = 1;
148 for (j = 1;j < 64;++j) xzm[j] = 0;
149
150 for (pos = 254;pos >= 0;--pos) {
151 b = e[pos / 8] >> (pos & 7);
152 b &= 1;
153 select(xzmb,xzm1b,xzm,xzm1,b);
154 add(a0,xzmb,xzmb + 32);
155 sub(a0 + 32,xzmb,xzmb + 32);
156 add(a1,xzm1b,xzm1b + 32);
157 sub(a1 + 32,xzm1b,xzm1b + 32);
158 square(b0,a0);
159 square(b0 + 32,a0 + 32);
160 mult(b1,a1,a0 + 32);
161 mult(b1 + 32,a1 + 32,a0);
162 add(c1,b1,b1 + 32);
163 sub(c1 + 32,b1,b1 + 32);
164 square(r,c1 + 32);
165 sub(s,b0,b0 + 32);
166 mult121665(t,s);
167 add(u,t,b0);
168 mult(xznb,b0,b0 + 32);
169 mult(xznb + 32,s,u);
170 square(xzn1b,c1);
171 mult(xzn1b + 32,r,work);
172 select(xzm,xzm1,xznb,xzn1b,b);
173 }
174
175 for (j = 0;j < 64;++j) work[j] = xzm[j];
176 }
177
recip(unsigned int out[32],const unsigned int z[32])178 static void recip(unsigned int out[32],const unsigned int z[32])
179 {
180 unsigned int z2[32];
181 unsigned int z9[32];
182 unsigned int z11[32];
183 unsigned int z2_5_0[32];
184 unsigned int z2_10_0[32];
185 unsigned int z2_20_0[32];
186 unsigned int z2_50_0[32];
187 unsigned int z2_100_0[32];
188 unsigned int t0[32];
189 unsigned int t1[32];
190 int i;
191
192 /* 2 */ square(z2,z);
193 /* 4 */ square(t1,z2);
194 /* 8 */ square(t0,t1);
195 /* 9 */ mult(z9,t0,z);
196 /* 11 */ mult(z11,z9,z2);
197 /* 22 */ square(t0,z11);
198 /* 2^5 - 2^0 = 31 */ mult(z2_5_0,t0,z9);
199
200 /* 2^6 - 2^1 */ square(t0,z2_5_0);
201 /* 2^7 - 2^2 */ square(t1,t0);
202 /* 2^8 - 2^3 */ square(t0,t1);
203 /* 2^9 - 2^4 */ square(t1,t0);
204 /* 2^10 - 2^5 */ square(t0,t1);
205 /* 2^10 - 2^0 */ mult(z2_10_0,t0,z2_5_0);
206
207 /* 2^11 - 2^1 */ square(t0,z2_10_0);
208 /* 2^12 - 2^2 */ square(t1,t0);
209 /* 2^20 - 2^10 */ for (i = 2;i < 10;i += 2) { square(t0,t1); square(t1,t0); }
210 /* 2^20 - 2^0 */ mult(z2_20_0,t1,z2_10_0);
211
212 /* 2^21 - 2^1 */ square(t0,z2_20_0);
213 /* 2^22 - 2^2 */ square(t1,t0);
214 /* 2^40 - 2^20 */ for (i = 2;i < 20;i += 2) { square(t0,t1); square(t1,t0); }
215 /* 2^40 - 2^0 */ mult(t0,t1,z2_20_0);
216
217 /* 2^41 - 2^1 */ square(t1,t0);
218 /* 2^42 - 2^2 */ square(t0,t1);
219 /* 2^50 - 2^10 */ for (i = 2;i < 10;i += 2) { square(t1,t0); square(t0,t1); }
220 /* 2^50 - 2^0 */ mult(z2_50_0,t0,z2_10_0);
221
222 /* 2^51 - 2^1 */ square(t0,z2_50_0);
223 /* 2^52 - 2^2 */ square(t1,t0);
224 /* 2^100 - 2^50 */ for (i = 2;i < 50;i += 2) { square(t0,t1); square(t1,t0); }
225 /* 2^100 - 2^0 */ mult(z2_100_0,t1,z2_50_0);
226
227 /* 2^101 - 2^1 */ square(t1,z2_100_0);
228 /* 2^102 - 2^2 */ square(t0,t1);
229 /* 2^200 - 2^100 */ for (i = 2;i < 100;i += 2) { square(t1,t0); square(t0,t1); }
230 /* 2^200 - 2^0 */ mult(t1,t0,z2_100_0);
231
232 /* 2^201 - 2^1 */ square(t0,t1);
233 /* 2^202 - 2^2 */ square(t1,t0);
234 /* 2^250 - 2^50 */ for (i = 2;i < 50;i += 2) { square(t0,t1); square(t1,t0); }
235 /* 2^250 - 2^0 */ mult(t0,t1,z2_50_0);
236
237 /* 2^251 - 2^1 */ square(t1,t0);
238 /* 2^252 - 2^2 */ square(t0,t1);
239 /* 2^253 - 2^3 */ square(t1,t0);
240 /* 2^254 - 2^4 */ square(t0,t1);
241 /* 2^255 - 2^5 */ square(t1,t0);
242 /* 2^255 - 21 */ mult(out,t1,z11);
243 }
244
crypto_scalarmult(unsigned char * q,const unsigned char * n,const unsigned char * p)245 static void crypto_scalarmult(unsigned char *q,
246 const unsigned char *n,
247 const unsigned char *p)
248 {
249 unsigned int work[96];
250 unsigned char e[32];
251 unsigned int i;
252 for (i = 0;i < 32;++i) e[i] = n[i];
253 e[0] &= 248;
254 e[31] &= 127;
255 e[31] |= 64;
256 for (i = 0;i < 32;++i) work[i] = p[i];
257 mainloop(work,e);
258 recip(work + 32,work + 32);
259 mult(work + 64,work,work + 32);
260 freeze(work + 64);
261 for (i = 0;i < 32;++i) q[i] = work[64 + i];
262 }
263
cf_curve25519_mul(uint8_t out[32],const uint8_t scalar[32],const uint8_t point[32])264 void cf_curve25519_mul(uint8_t out[32], const uint8_t scalar[32], const uint8_t point[32])
265 {
266 crypto_scalarmult(out, scalar, point);
267 }
268
cf_curve25519_mul_base(uint8_t out[32],const uint8_t scalar[32])269 void cf_curve25519_mul_base(uint8_t out[32], const uint8_t scalar[32])
270 {
271 uint8_t base_point[32] = { 9 };
272 cf_curve25519_mul(out, scalar, base_point);
273 }
274