|
| 1 | +# include <stdint.h> |
| 2 | + |
| 3 | +typedef unsigned __int128 uint128_t; |
| 4 | + |
| 5 | +typedef uint64_t limb; |
| 6 | +typedef uint128_t widelimb; |
| 7 | + |
| 8 | +typedef limb felem[4]; |
| 9 | +typedef widelimb widefelem[7]; |
| 10 | + |
| 11 | +felem p = {0x1FFFFFFFFFFFFFF, |
| 12 | + 0xFFFFFFFFFFFFFF, |
| 13 | + 0xFFFFE000000000, |
| 14 | + 0x00000000000002}; |
| 15 | + |
| 16 | + |
| 17 | +/*- |
| 18 | + * Reduce seven 128-bit coefficients to four 64-bit coefficients. |
| 19 | + * Requires in[i] < 2^126, |
| 20 | + * ensures out[0] < 2^56, out[1] < 2^56, out[2] < 2^56, out[3] <= 2^56 + 2^16 */ |
| 21 | +void reduce(felem out, const widefelem in) |
| 22 | +{ |
| 23 | + |
| 24 | + __CPROVER_assume(in[0]<(widelimb)((widelimb)1<<126)); |
| 25 | + __CPROVER_assume(in[1]<((widelimb)1<<126)); |
| 26 | + __CPROVER_assume(in[2]<((widelimb)1<<126)); |
| 27 | + __CPROVER_assume(in[3]<((widelimb)1<<126)); |
| 28 | + |
| 29 | + static const widelimb two127p15 = (((widelimb) 1) << 127) + |
| 30 | + (((widelimb) 1) << 15); |
| 31 | + static const widelimb two127m71 = (((widelimb) 1) << 127) - |
| 32 | + (((widelimb) 1) << 71); |
| 33 | + static const widelimb two127m71m55 = (((widelimb) 1) << 127) - |
| 34 | + (((widelimb) 1) << 71) - (((widelimb) 1) << 55); |
| 35 | + widelimb output[5]; |
| 36 | + |
| 37 | + /* Add 0 mod 2^224-2^96+1 to ensure all differences are positive */ |
| 38 | + output[0] = in[0] + two127p15; |
| 39 | + output[1] = in[1] + two127m71m55; |
| 40 | + output[2] = in[2] + two127m71; |
| 41 | + output[3] = in[3]; |
| 42 | + output[4] = in[4]; |
| 43 | + |
| 44 | + /* Eliminate in[4], in[5], in[6] */ |
| 45 | + output[4] += in[6] >> 16; |
| 46 | + output[3] += (in[6] & 0xffff) << 40; |
| 47 | + output[2] -= in[6]; |
| 48 | + |
| 49 | + output[3] += in[5] >> 16; |
| 50 | + output[2] += (in[5] & 0xffff) << 40; |
| 51 | + output[1] -= in[5]; |
| 52 | + |
| 53 | + output[2] += output[4] >> 16; |
| 54 | + output[1] += (output[4] & 0xffff) << 40; |
| 55 | + output[0] -= output[4]; |
| 56 | + |
| 57 | + /* Carry 2 -> 3 -> 4 */ |
| 58 | + output[3] += output[2] >> 56; |
| 59 | + output[2] &= 0x00ffffffffffffff; |
| 60 | + |
| 61 | + output[4] = output[3] >> 56; |
| 62 | + output[3] &= 0x00ffffffffffffff; |
| 63 | + |
| 64 | + /* Now output[2] < 2^56, output[3] < 2^56, output[4] < 2^72 */ |
| 65 | + |
| 66 | + /* Eliminate output[4] */ |
| 67 | + output[2] += output[4] >> 16; |
| 68 | + /* output[2] < 2^56 + 2^56 = 2^57 */ |
| 69 | + output[1] += (output[4] & 0xffff) << 40; |
| 70 | + output[0] -= output[4]; |
| 71 | + |
| 72 | + /* Carry 0 -> 1 -> 2 -> 3 */ |
| 73 | + output[1] += output[0] >> 56; |
| 74 | + out[0] = output[0] & 0x00ffffffffffffff; |
| 75 | + |
| 76 | + output[2] += output[1] >> 56; |
| 77 | + /* output[2] < 2^57 + 2^72 */ |
| 78 | + |
| 79 | + assert(output[2] < (((widelimb)1)<<57)+(((widelimb)1)<<72)); |
| 80 | + |
| 81 | + out[1] = output[1] & 0x00ffffffffffffff; |
| 82 | + output[3] += output[2] >> 56; |
| 83 | + /* output[3] <= 2^56 + 2^16 */ |
| 84 | + assert(output[3] < (((widelimb)1)<<56)+(((widelimb)1)<<16)); |
| 85 | + |
| 86 | + out[2] = output[2] & 0x00ffffffffffffff; |
| 87 | + |
| 88 | + /*- |
| 89 | + * out[0] < 2^56, out[1] < 2^56, out[2] < 2^56, |
| 90 | + * out[3] <= 2^56 + 2^16 (due to final carry), |
| 91 | + * so out < 2*p |
| 92 | + * p[0],p[1],p[2],p[3] |
| 93 | + * p = p[0] + 2**56 * p[1] ... |
| 94 | + * |
| 95 | + * out[3] < p[3] || |
| 96 | + * out[3] == p[3] && ( |
| 97 | + * |
| 98 | + */ |
| 99 | + out[3] = output[3]; |
| 100 | + assert(out[3] < p[0] || |
| 101 | + (out[3] == p[0] && ( |
| 102 | + out[2] < p[1] || (out[2]==p[1] && (out[1]<p[2] || out[1]==p[2] |
| 103 | + && (out[0]<p[3])))))); |
| 104 | +} |
0 commit comments