diff options
Diffstat (limited to 'src/math')
59 files changed, 5677 insertions, 0 deletions
diff --git a/src/math/bigint/big_code.cpp b/src/math/bigint/big_code.cpp new file mode 100644 index 000000000..e7a5e4946 --- /dev/null +++ b/src/math/bigint/big_code.cpp @@ -0,0 +1,152 @@ +/************************************************* +* BigInt Encoding/Decoding Source File * +* (C) 1999-2007 Jack Lloyd * +*************************************************/ + +#include <botan/bigint.h> +#include <botan/divide.h> +#include <botan/charset.h> +#include <botan/hex.h> + +namespace Botan { + +/************************************************* +* Encode a BigInt * +*************************************************/ +void BigInt::encode(byte output[], const BigInt& n, Base base) + { + if(base == Binary) + n.binary_encode(output); + else if(base == Hexadecimal) + { + SecureVector<byte> binary(n.encoded_size(Binary)); + n.binary_encode(binary); + for(u32bit j = 0; j != binary.size(); ++j) + Hex_Encoder::encode(binary[j], output + 2*j); + } + else if(base == Octal) + { + BigInt copy = n; + const u32bit output_size = n.encoded_size(Octal); + for(u32bit j = 0; j != output_size; ++j) + { + output[output_size - 1 - j] = Charset::digit2char(copy % 8); + copy /= 8; + } + } + else if(base == Decimal) + { + BigInt copy = n; + BigInt remainder; + copy.set_sign(Positive); + const u32bit output_size = n.encoded_size(Decimal); + for(u32bit j = 0; j != output_size; ++j) + { + divide(copy, 10, copy, remainder); + output[output_size - 1 - j] = + Charset::digit2char(remainder.word_at(0)); + if(copy.is_zero()) + break; + } + } + else + throw Invalid_Argument("Unknown BigInt encoding method"); + } + +/************************************************* +* Encode a BigInt * +*************************************************/ +SecureVector<byte> BigInt::encode(const BigInt& n, Base base) + { + SecureVector<byte> output(n.encoded_size(base)); + encode(output, n, base); + if(base != Binary) + for(u32bit j = 0; j != output.size(); ++j) + if(output[j] == 0) + output[j] = '0'; + return output; + } + +/************************************************* +* Encode a BigInt, with leading 0s if needed * +*************************************************/ +SecureVector<byte> BigInt::encode_1363(const BigInt& n, u32bit bytes) + { + const u32bit n_bytes = n.bytes(); + if(n_bytes > bytes) + throw Encoding_Error("encode_1363: n is too large to encode properly"); + + const u32bit leading_0s = bytes - n_bytes; + + SecureVector<byte> output(bytes); + encode(output + leading_0s, n, Binary); + return output; + } + +/************************************************* +* Decode a BigInt * +*************************************************/ +BigInt BigInt::decode(const MemoryRegion<byte>& buf, Base base) + { + return BigInt::decode(buf, buf.size(), base); + } + +/************************************************* +* Decode a BigInt * +*************************************************/ +BigInt BigInt::decode(const byte buf[], u32bit length, Base base) + { + BigInt r; + if(base == Binary) + r.binary_decode(buf, length); + else if(base == Hexadecimal) + { + SecureVector<byte> hex; + for(u32bit j = 0; j != length; ++j) + if(Hex_Decoder::is_valid(buf[j])) + hex.append(buf[j]); + + u32bit offset = (hex.size() % 2); + SecureVector<byte> binary(hex.size() / 2 + offset); + + if(offset) + { + byte temp[2] = { '0', hex[0] }; + binary[0] = Hex_Decoder::decode(temp); + } + + for(u32bit j = offset; j != binary.size(); ++j) + binary[j] = Hex_Decoder::decode(hex+2*j-offset); + r.binary_decode(binary, binary.size()); + } + else if(base == Decimal || base == Octal) + { + const u32bit RADIX = ((base == Decimal) ? 10 : 8); + for(u32bit j = 0; j != length; ++j) + { + if(Charset::is_space(buf[j])) + continue; + + if(!Charset::is_digit(buf[j])) + throw Invalid_Argument("BigInt::decode: " + "Invalid character in decimal input"); + + byte x = Charset::char2digit(buf[j]); + if(x >= RADIX) + { + if(RADIX == 10) + throw Invalid_Argument("BigInt: Invalid decimal string"); + else + throw Invalid_Argument("BigInt: Invalid octal string"); + } + + r *= RADIX; + r += x; + } + } + else + throw Invalid_Argument("Unknown BigInt decoding method"); + return r; + } + +} diff --git a/src/math/bigint/big_io.cpp b/src/math/bigint/big_io.cpp new file mode 100644 index 000000000..3c201e8b2 --- /dev/null +++ b/src/math/bigint/big_io.cpp @@ -0,0 +1,53 @@ +/************************************************* +* BigInt Input/Output Source File * +* (C) 1999-2007 Jack Lloyd * +*************************************************/ + +#include <botan/bigint.h> +#include <iostream> + +namespace Botan { + +/************************************************* +* Write the BigInt into a stream * +*************************************************/ +std::ostream& operator<<(std::ostream& stream, const BigInt& n) + { + BigInt::Base base = BigInt::Decimal; + if(stream.flags() & std::ios::hex) + base = BigInt::Hexadecimal; + else if(stream.flags() & std::ios::oct) + base = BigInt::Octal; + + if(n == 0) + stream.write("0", 1); + else + { + if(n < 0) + stream.write("-", 1); + SecureVector<byte> buffer = BigInt::encode(n, base); + u32bit skip = 0; + while(buffer[skip] == '0' && skip < buffer.size()) + ++skip; + stream.write(reinterpret_cast<const char*>(buffer.begin()) + skip, + buffer.size() - skip); + } + if(!stream.good()) + throw Stream_IO_Error("BigInt output operator has failed"); + return stream; + } + +/************************************************* +* Read the BigInt from a stream * +*************************************************/ +std::istream& operator>>(std::istream& stream, BigInt& n) + { + std::string str; + std::getline(stream, str); + if(stream.bad() || (stream.fail() && !stream.eof())) + throw Stream_IO_Error("BigInt input operator has failed"); + n = BigInt(str); + return stream; + } + +} diff --git a/src/math/bigint/big_ops2.cpp b/src/math/bigint/big_ops2.cpp new file mode 100644 index 000000000..ef083f394 --- /dev/null +++ b/src/math/bigint/big_ops2.cpp @@ -0,0 +1,222 @@ +/************************************************* +* BigInt Assignment Operators Source File * +* (C) 1999-2007 Jack Lloyd * +*************************************************/ + +#include <botan/bigint.h> +#include <botan/mp_core.h> +#include <botan/bit_ops.h> +#include <algorithm> + +namespace Botan { + +/************************************************* +* Addition Operator * +*************************************************/ +BigInt& BigInt::operator+=(const BigInt& y) + { + const u32bit x_sw = sig_words(), y_sw = y.sig_words(); + + const u32bit reg_size = std::max(x_sw, y_sw) + 1; + grow_to(reg_size); + + if(sign() == y.sign()) + bigint_add2(get_reg(), reg_size - 1, y.data(), y_sw); + else + { + s32bit relative_size = bigint_cmp(data(), x_sw, y.data(), y_sw); + + if(relative_size < 0) + { + SecureVector<word> z(reg_size - 1); + bigint_sub3(z, y.data(), reg_size - 1, data(), x_sw); + copy_mem(get_reg().begin(), z.begin(), z.size()); + set_sign(y.sign()); + } + else if(relative_size == 0) + { + get_reg().clear(); + set_sign(Positive); + } + else if(relative_size > 0) + bigint_sub2(get_reg(), x_sw, y.data(), y_sw); + } + + return (*this); + } + +/************************************************* +* Subtraction Operator * +*************************************************/ +BigInt& BigInt::operator-=(const BigInt& y) + { + const u32bit x_sw = sig_words(), y_sw = y.sig_words(); + + s32bit relative_size = bigint_cmp(data(), x_sw, y.data(), y_sw); + + const u32bit reg_size = std::max(x_sw, y_sw) + 1; + grow_to(reg_size); + + if(relative_size < 0) + { + if(sign() == y.sign()) + { + SecureVector<word> z(reg_size - 1); + bigint_sub3(z, y.data(), reg_size - 1, data(), x_sw); + copy_mem(get_reg().begin(), z.begin(), z.size()); + } + else + bigint_add2(get_reg(), reg_size - 1, y.data(), y_sw); + + set_sign(y.reverse_sign()); + } + else if(relative_size == 0) + { + if(sign() == y.sign()) + { + get_reg().clear(); + set_sign(Positive); + } + else + bigint_shl1(get_reg(), x_sw, 0, 1); + } + else if(relative_size > 0) + { + if(sign() == y.sign()) + bigint_sub2(get_reg(), x_sw, y.data(), y_sw); + else + bigint_add2(get_reg(), reg_size - 1, y.data(), y_sw); + } + + return (*this); + } + +/************************************************* +* Multiplication Operator * +*************************************************/ +BigInt& BigInt::operator*=(const BigInt& y) + { + const u32bit x_sw = sig_words(), y_sw = y.sig_words(); + set_sign((sign() == y.sign()) ? Positive : Negative); + + if(x_sw == 0 || y_sw == 0) + { + get_reg().clear(); + set_sign(Positive); + } + else if(x_sw == 1 && y_sw) + { + grow_to(y_sw + 2); + bigint_linmul3(get_reg(), y.data(), y_sw, word_at(0)); + } + else if(y_sw == 1 && x_sw) + { + grow_to(x_sw + 2); + bigint_linmul2(get_reg(), x_sw, y.word_at(0)); + } + else + { + grow_to(size() + y.size()); + + SecureVector<word> z(data(), x_sw); + SecureVector<word> workspace(size()); + + bigint_mul(get_reg(), size(), workspace, + z, z.size(), x_sw, + y.data(), y.size(), y_sw); + } + + return (*this); + } + +/************************************************* +* Division Operator * +*************************************************/ +BigInt& BigInt::operator/=(const BigInt& y) + { + if(y.sig_words() == 1 && power_of_2(y.word_at(0))) + (*this) >>= (y.bits() - 1); + else + (*this) = (*this) / y; + return (*this); + } + +/************************************************* +* Modulo Operator * +*************************************************/ +BigInt& BigInt::operator%=(const BigInt& mod) + { + return (*this = (*this) % mod); + } + +/************************************************* +* Modulo Operator * +*************************************************/ +word BigInt::operator%=(word mod) + { + if(mod == 0) + throw BigInt::DivideByZero(); + if(power_of_2(mod)) + { + word result = (word_at(0) & (mod - 1)); + clear(); + grow_to(2); + get_reg()[0] = result; + return result; + } + + word remainder = 0; + + for(u32bit j = sig_words(); j > 0; --j) + remainder = bigint_modop(remainder, word_at(j-1), mod); + clear(); + grow_to(2); + + if(remainder && sign() == BigInt::Negative) + get_reg()[0] = mod - remainder; + else + get_reg()[0] = remainder; + + set_sign(BigInt::Positive); + + return word_at(0); + } + +/************************************************* +* Left Shift Operator * +*************************************************/ +BigInt& BigInt::operator<<=(u32bit shift) + { + if(shift) + { + const u32bit shift_words = shift / MP_WORD_BITS, + shift_bits = shift % MP_WORD_BITS, + words = sig_words(); + + grow_to(words + shift_words + (shift_bits ? 1 : 0)); + bigint_shl1(get_reg(), words, shift_words, shift_bits); + } + + return (*this); + } + +/************************************************* +* Right Shift Operator * +*************************************************/ +BigInt& BigInt::operator>>=(u32bit shift) + { + if(shift) + { + const u32bit shift_words = shift / MP_WORD_BITS, + shift_bits = shift % MP_WORD_BITS; + + bigint_shr1(get_reg(), sig_words(), shift_words, shift_bits); + + if(is_zero()) + set_sign(Positive); + } + + return (*this); + } + +} diff --git a/src/math/bigint/big_ops3.cpp b/src/math/bigint/big_ops3.cpp new file mode 100644 index 000000000..ff24eab1c --- /dev/null +++ b/src/math/bigint/big_ops3.cpp @@ -0,0 +1,188 @@ +/************************************************* +* BigInt Binary Operators Source File * +* (C) 1999-2007 Jack Lloyd * +*************************************************/ + +#include <botan/bigint.h> +#include <botan/divide.h> +#include <botan/mp_core.h> +#include <botan/bit_ops.h> +#include <algorithm> + +namespace Botan { + +/************************************************* +* Addition Operator * +*************************************************/ +BigInt operator+(const BigInt& x, const BigInt& y) + { + const u32bit x_sw = x.sig_words(), y_sw = y.sig_words(); + + BigInt z(x.sign(), std::max(x_sw, y_sw) + 1); + + if((x.sign() == y.sign())) + bigint_add3(z.get_reg(), x.data(), x_sw, y.data(), y_sw); + else + { + s32bit relative_size = bigint_cmp(x.data(), x_sw, y.data(), y_sw); + + if(relative_size < 0) + { + bigint_sub3(z.get_reg(), y.data(), y_sw, x.data(), x_sw); + z.set_sign(y.sign()); + } + else if(relative_size == 0) + z.set_sign(BigInt::Positive); + else if(relative_size > 0) + bigint_sub3(z.get_reg(), x.data(), x_sw, y.data(), y_sw); + } + + return z; + } + +/************************************************* +* Subtraction Operator * +*************************************************/ +BigInt operator-(const BigInt& x, const BigInt& y) + { + const u32bit x_sw = x.sig_words(), y_sw = y.sig_words(); + + s32bit relative_size = bigint_cmp(x.data(), x_sw, y.data(), y_sw); + + BigInt z(BigInt::Positive, std::max(x_sw, y_sw) + 1); + + if(relative_size < 0) + { + if(x.sign() == y.sign()) + bigint_sub3(z.get_reg(), y.data(), y_sw, x.data(), x_sw); + else + bigint_add3(z.get_reg(), x.data(), x_sw, y.data(), y_sw); + z.set_sign(y.reverse_sign()); + } + else if(relative_size == 0) + { + if(x.sign() != y.sign()) + bigint_shl2(z.get_reg(), x.data(), x_sw, 0, 1); + } + else if(relative_size > 0) + { + if(x.sign() == y.sign()) + bigint_sub3(z.get_reg(), x.data(), x_sw, y.data(), y_sw); + else + bigint_add3(z.get_reg(), x.data(), x_sw, y.data(), y_sw); + z.set_sign(x.sign()); + } + return z; + } + +/************************************************* +* Multiplication Operator * +*************************************************/ +BigInt operator*(const BigInt& x, const BigInt& y) + { + const u32bit x_sw = x.sig_words(), y_sw = y.sig_words(); + + BigInt z(BigInt::Positive, x.size() + y.size()); + + if(x_sw == 1 && y_sw) + bigint_linmul3(z.get_reg(), y.data(), y_sw, x.word_at(0)); + else if(y_sw == 1 && x_sw) + bigint_linmul3(z.get_reg(), x.data(), x_sw, y.word_at(0)); + else if(x_sw && y_sw) + { + SecureVector<word> workspace(z.size()); + bigint_mul(z.get_reg(), z.size(), workspace, + x.data(), x.size(), x_sw, + y.data(), y.size(), y_sw); + } + + if(x_sw && y_sw && x.sign() != y.sign()) + z.flip_sign(); + return z; + } + +/************************************************* +* Division Operator * +*************************************************/ +BigInt operator/(const BigInt& x, const BigInt& y) + { + BigInt q, r; + divide(x, y, q, r); + return q; + } + +/************************************************* +* Modulo Operator * +*************************************************/ +BigInt operator%(const BigInt& n, const BigInt& mod) + { + if(mod.is_zero()) + throw BigInt::DivideByZero(); + if(mod.is_negative()) + throw Invalid_Argument("BigInt::operator%: modulus must be > 0"); + if(n.is_positive() && mod.is_positive() && n < mod) + return n; + + BigInt q, r; + divide(n, mod, q, r); + return r; + } + +/************************************************* +* Modulo Operator * +*************************************************/ +word operator%(const BigInt& n, word mod) + { + if(mod == 0) + throw BigInt::DivideByZero(); + if(power_of_2(mod)) + return (n.word_at(0) & (mod - 1)); + + word remainder = 0; + + for(u32bit j = n.sig_words(); j > 0; --j) + remainder = bigint_modop(remainder, n.word_at(j-1), mod); + + if(remainder && n.sign() == BigInt::Negative) + return mod - remainder; + return remainder; + } + +/************************************************* +* Left Shift Operator * +*************************************************/ +BigInt operator<<(const BigInt& x, u32bit shift) + { + if(shift == 0) + return x; + + const u32bit shift_words = shift / MP_WORD_BITS, + shift_bits = shift % MP_WORD_BITS; + + const u32bit x_sw = x.sig_words(); + + BigInt y(x.sign(), x_sw + shift_words + (shift_bits ? 1 : 0)); + bigint_shl2(y.get_reg(), x.data(), x_sw, shift_words, shift_bits); + return y; + } + +/************************************************* +* Right Shift Operator * +*************************************************/ +BigInt operator>>(const BigInt& x, u32bit shift) + { + if(shift == 0) + return x; + if(x.bits() <= shift) + return 0; + + const u32bit shift_words = shift / MP_WORD_BITS, + shift_bits = shift % MP_WORD_BITS, + x_sw = x.sig_words(); + + BigInt y(x.sign(), x_sw - shift_words); + bigint_shr2(y.get_reg(), x.data(), x_sw, shift_words, shift_bits); + return y; + } + +} diff --git a/src/math/bigint/big_rand.cpp b/src/math/bigint/big_rand.cpp new file mode 100644 index 000000000..055873642 --- /dev/null +++ b/src/math/bigint/big_rand.cpp @@ -0,0 +1,59 @@ +/************************************************* +* BigInt Random Generation Source File * +* (C) 1999-2007 Jack Lloyd * +*************************************************/ + +#include <botan/bigint.h> +#include <botan/parsing.h> + +namespace Botan { + +/************************************************* +* Construct a BigInt of a specific form * +*************************************************/ +BigInt::BigInt(NumberType type, u32bit bits) + { + set_sign(Positive); + + if(type == Power2) + set_bit(bits); + else + throw Invalid_Argument("BigInt(NumberType): Unknown type"); + } + +/************************************************* +* Randomize this number * +*************************************************/ +void BigInt::randomize(RandomNumberGenerator& rng, + u32bit bitsize) + { + set_sign(Positive); + + if(bitsize == 0) + clear(); + else + { + SecureVector<byte> array((bitsize + 7) / 8); + rng.randomize(array, array.size()); + if(bitsize % 8) + array[0] &= 0xFF >> (8 - (bitsize % 8)); + array[0] |= 0x80 >> ((bitsize % 8) ? (8 - bitsize % 8) : 0); + binary_decode(array, array.size()); + } + } + +/************************************************* +* Generate a random integer within given range * +*************************************************/ +BigInt BigInt::random_integer(RandomNumberGenerator& rng, + const BigInt& min, const BigInt& max) + { + BigInt range = max - min; + + if(range <= 0) + throw Invalid_Argument("random_integer: invalid min/max values"); + + return (min + (BigInt(rng, range.bits() + 2) % range)); + } + +} diff --git a/src/math/bigint/bigint.cpp b/src/math/bigint/bigint.cpp new file mode 100644 index 000000000..e3c7931e6 --- /dev/null +++ b/src/math/bigint/bigint.cpp @@ -0,0 +1,367 @@ +/************************************************* +* BigInt Base Source File * +* (C) 1999-2008 Jack Lloyd * +*************************************************/ + +#include <botan/bigint.h> +#include <botan/mp_core.h> +#include <botan/loadstor.h> +#include <botan/parsing.h> +#include <botan/util.h> + +namespace Botan { + +/************************************************* +* Construct a BigInt from a regular number * +*************************************************/ +BigInt::BigInt(u64bit n) + { + set_sign(Positive); + + if(n == 0) + return; + + const u32bit limbs_needed = sizeof(u64bit) / sizeof(word); + + reg.create(4*limbs_needed); + for(u32bit j = 0; j != limbs_needed; ++j) + reg[j] = ((n >> (j*MP_WORD_BITS)) & MP_WORD_MASK); + } + +/************************************************* +* Construct a BigInt of the specified size * +*************************************************/ +BigInt::BigInt(Sign s, u32bit size) + { + reg.create(round_up(size, 8)); + signedness = s; + } + +/************************************************* +* Construct a BigInt from a "raw" BigInt * +*************************************************/ +BigInt::BigInt(const BigInt& b) + { + const u32bit b_words = b.sig_words(); + + if(b_words) + { + reg.create(round_up(b_words, 8)); + reg.copy(b.data(), b_words); + set_sign(b.sign()); + } + else + { + reg.create(2); + set_sign(Positive); + } + } + +/************************************************* +* Construct a BigInt from a string * +*************************************************/ +BigInt::BigInt(const std::string& str) + { + Base base = Decimal; + u32bit markers = 0; + bool negative = false; + if(str.length() > 0 && str[0] == '-') { markers += 1; negative = true; } + + if(str.length() > markers + 2 && str[markers ] == '0' && + str[markers + 1] == 'x') + { markers += 2; base = Hexadecimal; } + else if(str.length() > markers + 1 && str[markers] == '0') + { markers += 1; base = Octal; } + + *this = decode(reinterpret_cast<const byte*>(str.data()) + markers, + str.length() - markers, base); + + if(negative) set_sign(Negative); + else set_sign(Positive); + } + +/************************************************* +* Construct a BigInt from an encoded BigInt * +*************************************************/ +BigInt::BigInt(const byte input[], u32bit length, Base base) + { + set_sign(Positive); + *this = decode(input, length, base); + } + +/************************************************* +* Construct a BigInt from an encoded BigInt * +*************************************************/ +BigInt::BigInt(RandomNumberGenerator& rng, u32bit bits) + { + set_sign(Positive); + randomize(rng, bits); + } + +/************************************************* +* Swap this BigInt with another * +*************************************************/ +void BigInt::swap(BigInt& other) + { + reg.swap(other.reg); + std::swap(signedness, other.signedness); + } + +/************************************************* +* Grow the internal storage * +*************************************************/ +void BigInt::grow_reg(u32bit n) + { + reg.grow_to(round_up(size() + n, 8)); + } + +/************************************************* +* Grow the internal storage * +*************************************************/ +void BigInt::grow_to(u32bit n) + { + if(n > size()) + reg.grow_to(round_up(n, 8)); + } + +/************************************************* +* Comparison Function * +*************************************************/ +s32bit BigInt::cmp(const BigInt& n, bool check_signs) const + { + if(check_signs) + { + if(n.is_positive() && this->is_negative()) return -1; + if(n.is_negative() && this->is_positive()) return 1; + if(n.is_negative() && this->is_negative()) + return (-bigint_cmp(data(), sig_words(), n.data(), n.sig_words())); + } + return bigint_cmp(data(), sig_words(), n.data(), n.sig_words()); + } + +/************************************************* +* Convert this number to a u32bit, if possible * +*************************************************/ +u32bit BigInt::to_u32bit() const + { + if(is_negative()) + throw Encoding_Error("BigInt::to_u32bit: Number is negative"); + if(bits() >= 32) + throw Encoding_Error("BigInt::to_u32bit: Number is too big to convert"); + + u32bit out = 0; + for(u32bit j = 0; j != 4; ++j) + out = (out << 8) | byte_at(3-j); + return out; + } + +/************************************************* +* Return byte n of this number * +*************************************************/ +byte BigInt::byte_at(u32bit n) const + { + const u32bit WORD_BYTES = sizeof(word); + u32bit word_num = n / WORD_BYTES, byte_num = n % WORD_BYTES; + if(word_num >= size()) + return 0; + else + return get_byte(WORD_BYTES - byte_num - 1, reg[word_num]); + } + +/************************************************* +* Return bit n of this number * +*************************************************/ +bool BigInt::get_bit(u32bit n) const + { + return ((word_at(n / MP_WORD_BITS) >> (n % MP_WORD_BITS)) & 1); + } + +/************************************************* +* Return bits {offset...offset+length} * +*************************************************/ +u32bit BigInt::get_substring(u32bit offset, u32bit length) const + { + if(length > 32) + throw Invalid_Argument("BigInt::get_substring: Substring size too big"); + + u64bit piece = 0; + for(u32bit j = 0; j != 8; ++j) + piece = (piece << 8) | byte_at((offset / 8) + (7-j)); + + u64bit mask = (1 << length) - 1; + u32bit shift = (offset % 8); + + return static_cast<u32bit>((piece >> shift) & mask); + } + +/************************************************* +* Set bit number n * +*************************************************/ +void BigInt::set_bit(u32bit n) + { + const u32bit which = n / MP_WORD_BITS; + const word mask = static_cast<word>(1) << (n % MP_WORD_BITS); + if(which >= size()) grow_to(which + 1); + reg[which] |= mask; + } + +/************************************************* +* Clear bit number n * +*************************************************/ +void BigInt::clear_bit(u32bit n) + { + const u32bit which = n / MP_WORD_BITS; + const word mask = static_cast<word>(1) << (n % MP_WORD_BITS); + if(which < size()) + reg[which] &= ~mask; + } + +/************************************************* +* Clear all but the lowest n bits * +*************************************************/ +void BigInt::mask_bits(u32bit n) + { + if(n == 0) { clear(); return; } + if(n >= bits()) return; + + const u32bit top_word = n / MP_WORD_BITS; + const word mask = (static_cast<word>(1) << (n % MP_WORD_BITS)) - 1; + + if(top_word < size()) + for(u32bit j = top_word + 1; j != size(); ++j) + reg[j] = 0; + + reg[top_word] &= mask; + } + +/************************************************* +* Count how many bytes are being used * +*************************************************/ +u32bit BigInt::bytes() const + { + return (bits() + 7) / 8; + } + +/************************************************* +* Count how many bits are being used * +*************************************************/ +u32bit BigInt::bits() const + { + if(sig_words() == 0) + return 0; + + u32bit full_words = sig_words() - 1, top_bits = MP_WORD_BITS; + word top_word = word_at(full_words), mask = MP_WORD_TOP_BIT; + + while(top_bits && ((top_word & mask) == 0)) + { mask >>= 1; top_bits--; } + + return (full_words * MP_WORD_BITS + top_bits); + } + +/************************************************* +* Calcluate the size in a certain base * +*************************************************/ +u32bit BigInt::encoded_size(Base base) const + { + static const double LOG_2_BASE_10 = 0.30102999566; + + if(base == Binary) + return bytes(); + else if(base == Hexadecimal) + return 2*bytes(); + else if(base == Octal) + return ((bits() + 2) / 3); + else if(base == Decimal) + return static_cast<u32bit>((bits() * LOG_2_BASE_10) + 1); + else + throw Invalid_Argument("Unknown base for BigInt encoding"); + } + +/************************************************* +* Set the sign * +*************************************************/ +void BigInt::set_sign(Sign s) + { + if(is_zero()) + signedness = Positive; + else + signedness = s; + } + +/************************************************* +* Reverse the value of the sign flag * +*************************************************/ +void BigInt::flip_sign() + { + set_sign(reverse_sign()); + } + +/************************************************* +* Return the opposite value of the current sign * +*************************************************/ +BigInt::Sign BigInt::reverse_sign() const + { + if(sign() == Positive) + return Negative; + return Positive; + } + +/************************************************* +* Return the negation of this number * +*************************************************/ +BigInt BigInt::operator-() const + { + BigInt x = (*this); + x.flip_sign(); + return x; + } + +/************************************************* +* Return the absolute value of this number * +*************************************************/ +BigInt BigInt::abs() const + { + BigInt x = (*this); + x.set_sign(Positive); + return x; + } + +/************************************************* +* Encode this number into bytes * +*************************************************/ +void BigInt::binary_encode(byte output[]) const + { + const u32bit sig_bytes = bytes(); + for(u32bit j = 0; j != sig_bytes; ++j) + output[sig_bytes-j-1] = byte_at(j); + } + +/************************************************* +* Set this number to the value in buf * +*************************************************/ +void BigInt::binary_decode(const byte buf[], u32bit length) + { + const u32bit WORD_BYTES = sizeof(word); + + reg.create(round_up((length / WORD_BYTES) + 1, 8)); + + for(u32bit j = 0; j != length / WORD_BYTES; ++j) + { + u32bit top = length - WORD_BYTES*j; + for(u32bit k = WORD_BYTES; k > 0; --k) + reg[j] = (reg[j] << 8) | buf[top - k]; + } + for(u32bit j = 0; j != length % WORD_BYTES; ++j) + reg[length / WORD_BYTES] = (reg[length / WORD_BYTES] << 8) | buf[j]; + } + +/************************************************* +* Set this number to the value in buf * +*************************************************/ +void BigInt::binary_decode(const MemoryRegion<byte>& buf) + { + binary_decode(buf, buf.size()); + } + +} diff --git a/src/math/bigint/bigint.h b/src/math/bigint/bigint.h new file mode 100644 index 000000000..bb1d7ef7b --- /dev/null +++ b/src/math/bigint/bigint.h @@ -0,0 +1,184 @@ +/************************************************* +* BigInt Header File * +* (C) 1999-2008 Jack Lloyd * +*************************************************/ + +#ifndef BOTAN_BIGINT_H__ +#define BOTAN_BIGINT_H__ + +#include <botan/rng.h> +#include <botan/secmem.h> +#include <botan/mp_types.h> +#include <iosfwd> + +namespace Botan { + +/************************************************* +* BigInt * +*************************************************/ +class BOTAN_DLL BigInt + { + public: + enum Base { Octal = 8, Decimal = 10, Hexadecimal = 16, Binary = 256 }; + enum Sign { Negative = 0, Positive = 1 }; + enum NumberType { Power2 }; + + struct DivideByZero : public Exception + { DivideByZero() : Exception("BigInt divide by zero") {} }; + + BigInt& operator+=(const BigInt&); + BigInt& operator-=(const BigInt&); + + BigInt& operator*=(const BigInt&); + BigInt& operator/=(const BigInt&); + BigInt& operator%=(const BigInt&); + word operator%=(word); + BigInt& operator<<=(u32bit); + BigInt& operator>>=(u32bit); + + BigInt& operator++() { return (*this += 1); } + BigInt& operator--() { return (*this -= 1); } + BigInt operator++(int) { BigInt x = (*this); ++(*this); return x; } + BigInt operator--(int) { BigInt x = (*this); --(*this); return x; } + + BigInt operator-() const; + bool operator !() const { return (!is_nonzero()); } + + s32bit cmp(const BigInt&, bool = true) const; + bool is_even() const { return (get_bit(0) == 0); } + bool is_odd() const { return (get_bit(0) == 1); } + + bool is_zero() const + { + const u32bit sw = sig_words(); + + for(u32bit i = 0; i != sw; ++i) + if(reg[i]) + return false; + return true; + } + + bool is_nonzero() const { return (!is_zero()); } + + void set_bit(u32bit); + void clear_bit(u32bit); + void mask_bits(u32bit); + + bool get_bit(u32bit) const; + u32bit get_substring(u32bit, u32bit) const; + byte byte_at(u32bit) const; + + // same as operator[], remove this + word word_at(u32bit n) const + { return ((n < size()) ? get_reg()[n] : 0); } + + u32bit to_u32bit() const; + + bool is_negative() const { return (sign() == Negative); } + bool is_positive() const { return (sign() == Positive); } + Sign sign() const { return (signedness); } + Sign reverse_sign() const; + void flip_sign(); + void set_sign(Sign); + BigInt abs() const; + + u32bit size() const { return get_reg().size(); } + + u32bit sig_words() const + { + const word* x = reg.begin(); + u32bit sig = reg.size(); + + while(sig && (x[sig-1] == 0)) + sig--; + return sig; + } + + u32bit bytes() const; + u32bit bits() const; + + const word* data() const { return reg.begin(); } + SecureVector<word>& get_reg() { return reg; } + const SecureVector<word>& get_reg() const { return reg; } + + void grow_reg(u32bit); + void grow_to(u32bit); + + word& operator[](u32bit i) { return reg[i]; } + word operator[](u32bit i) const { return reg[i]; } + void clear() { get_reg().clear(); } + + void randomize(RandomNumberGenerator& rng, u32bit n); + + void binary_encode(byte[]) const; + void binary_decode(const byte[], u32bit); + void binary_decode(const MemoryRegion<byte>&); + u32bit encoded_size(Base = Binary) const; + + static BigInt random_integer(RandomNumberGenerator&, + const BigInt&, const BigInt&); + + static SecureVector<byte> encode(const BigInt&, Base = Binary); + static void encode(byte[], const BigInt&, Base = Binary); + static BigInt decode(const byte[], u32bit, Base = Binary); + static BigInt decode(const MemoryRegion<byte>&, Base = Binary); + static SecureVector<byte> encode_1363(const BigInt&, u32bit); + + void swap(BigInt&); + + BigInt() { signedness = Positive; } + BigInt(u64bit); + BigInt(const BigInt&); + BigInt(const std::string&); + BigInt(const byte[], u32bit, Base = Binary); + BigInt(RandomNumberGenerator& rng, u32bit bits); + BigInt(Sign, u32bit); + BigInt(NumberType, u32bit); + private: + SecureVector<word> reg; + Sign signedness; + }; + +/************************************************* +* Arithmetic Operators * +*************************************************/ +BigInt BOTAN_DLL operator+(const BigInt&, const BigInt&); +BigInt BOTAN_DLL operator-(const BigInt&, const BigInt&); +BigInt BOTAN_DLL operator*(const BigInt&, const BigInt&); +BigInt BOTAN_DLL operator/(const BigInt&, const BigInt&); +BigInt BOTAN_DLL operator%(const BigInt&, const BigInt&); +word BOTAN_DLL operator%(const BigInt&, word); +BigInt BOTAN_DLL operator<<(const BigInt&, u32bit); +BigInt BOTAN_DLL operator>>(const BigInt&, u32bit); + +/************************************************* +* Comparison Operators * +*************************************************/ +inline bool operator==(const BigInt& a, const BigInt& b) + { return (a.cmp(b) == 0); } +inline bool operator!=(const BigInt& a, const BigInt& b) + { return (a.cmp(b) != 0); } +inline bool operator<=(const BigInt& a, const BigInt& b) + { return (a.cmp(b) <= 0); } +inline bool operator>=(const BigInt& a, const BigInt& b) + { return (a.cmp(b) >= 0); } +inline bool operator<(const BigInt& a, const BigInt& b) + { return (a.cmp(b) < 0); } +inline bool operator>(const BigInt& a, const BigInt& b) + { return (a.cmp(b) > 0); } + +/************************************************* +* I/O Operators * +*************************************************/ +BOTAN_DLL std::ostream& operator<<(std::ostream&, const BigInt&); +BOTAN_DLL std::istream& operator>>(std::istream&, BigInt&); + +} + +namespace std { + +inline void swap(Botan::BigInt& a, Botan::BigInt& b) { a.swap(b); } + +} + +#endif diff --git a/src/math/bigint/divide.cpp b/src/math/bigint/divide.cpp new file mode 100644 index 000000000..ba088ced4 --- /dev/null +++ b/src/math/bigint/divide.cpp @@ -0,0 +1,104 @@ +/************************************************* +* Division Algorithm Source File * +* (C) 1999-2007 Jack Lloyd * +*************************************************/ + +#include <botan/divide.h> +#include <botan/mp_core.h> + +namespace Botan { + +namespace { + +/************************************************* +* Handle signed operands, if necessary * +*************************************************/ +void sign_fixup(const BigInt& x, const BigInt& y, BigInt& q, BigInt& r) + { + if(x.sign() == BigInt::Negative) + { + q.flip_sign(); + if(r.is_nonzero()) { --q; r = y.abs() - r; } + } + if(y.sign() == BigInt::Negative) + q.flip_sign(); + } + +} + +/************************************************* +* Solve x = q * y + r * +*************************************************/ +void divide(const BigInt& x, const BigInt& y_arg, BigInt& q, BigInt& r) + { + if(y_arg.is_zero()) + throw BigInt::DivideByZero(); + + BigInt y = y_arg; + const u32bit y_words = y.sig_words(); + r = x; + + r.set_sign(BigInt::Positive); + y.set_sign(BigInt::Positive); + + s32bit compare = r.cmp(y); + + if(compare < 0) + q = 0; + else if(compare == 0) + { + q = 1; + r = 0; + } + else + { + u32bit shifts = 0; + word y_top = y[y.sig_words()-1]; + while(y_top < MP_WORD_TOP_BIT) { y_top <<= 1; ++shifts; } + y <<= shifts; + r <<= shifts; + + const u32bit n = r.sig_words() - 1, t = y_words - 1; + + q.get_reg().create(n - t + 1); + if(n <= t) + { + while(r > y) { r -= y; ++q; } + r >>= shifts; + sign_fixup(x, y_arg, q, r); + return; + } + + BigInt temp = y << (MP_WORD_BITS * (n-t)); + + while(r >= temp) { r -= temp; ++q[n-t]; } + + for(u32bit j = n; j != t; --j) + { + const word x_j0 = r.word_at(j); + const word x_j1 = r.word_at(j-1); + const word y_t = y.word_at(t); + + if(x_j0 == y_t) + q[j-t-1] = MP_WORD_MAX; + else + q[j-t-1] = bigint_divop(x_j0, x_j1, y_t); + + while(bigint_divcore(q[j-t-1], y_t, y.word_at(t-1), + x_j0, x_j1, r.word_at(j-2))) + --q[j-t-1]; + + r -= (q[j-t-1] * y) << (MP_WORD_BITS * (j-t-1)); + if(r.is_negative()) + { + r += y << (MP_WORD_BITS * (j-t-1)); + --q[j-t-1]; + } + } + r >>= shifts; + } + + sign_fixup(x, y_arg, q, r); + } + +} diff --git a/src/math/bigint/divide.h b/src/math/bigint/divide.h new file mode 100644 index 000000000..d0de3da75 --- /dev/null +++ b/src/math/bigint/divide.h @@ -0,0 +1,17 @@ +/************************************************* +* Division Header File * +* (C) 1999-2007 Jack Lloyd * +*************************************************/ + +#ifndef BOTAN_DIVISON_ALGORITHM__ +#define BOTAN_DIVISON_ALGORITHM__ + +#include <botan/bigint.h> + +namespace Botan { + +void BOTAN_DLL divide(const BigInt&, const BigInt&, BigInt&, BigInt&); + +} + +#endif diff --git a/src/math/bigint/info.txt b/src/math/bigint/info.txt new file mode 100644 index 000000000..d6c5df763 --- /dev/null +++ b/src/math/bigint/info.txt @@ -0,0 +1,31 @@ +realname "BigInt" + +load_on auto + +define BIGINT + +<requires> +hex +mp_amd64|mp_asm64|mp_ia32|mp_ia32_msvc|mp_generic +monty_generic +mulop_generic +</requires> + +<add> +bigint.h +divide.h +mp_core.h +mp_types.h +big_code.cpp +big_io.cpp +big_ops2.cpp +big_ops3.cpp +big_rand.cpp +bigint.cpp +divide.cpp +mp_asm.cpp +mp_comba.cpp +mp_karat.cpp +mp_misc.cpp +mp_shift.cpp +</add> diff --git a/src/math/bigint/monty_amd64/info.txt b/src/math/bigint/monty_amd64/info.txt new file mode 100644 index 000000000..32308bf41 --- /dev/null +++ b/src/math/bigint/monty_amd64/info.txt @@ -0,0 +1,31 @@ +realname "Montgomery Reduction (x86-64)" + +mp_bits 64 + +load_on request + +<add> +mp_monty.S +</add> + +<requires> +asm_amd64 +</requires> + +<arch> +amd64 +</arch> + +<cc> +gcc +icc +</cc> + +# ELF systems +<os> +linux +freebsd +netbsd +openbsd +solaris +</os> diff --git a/src/math/bigint/monty_amd64/mp_monty.S b/src/math/bigint/monty_amd64/mp_monty.S new file mode 100644 index 000000000..3dd4040bc --- /dev/null +++ b/src/math/bigint/monty_amd64/mp_monty.S @@ -0,0 +1,397 @@ +/************************************************* +* Montgomery Reduction Source File * +* (C) 2008 Jack Lloyd * +*************************************************/ + +#include <botan/asm_macr.h> + +START_LISTING(mp_monty.S) + +START_FUNCTION(bigint_monty_redc) + pushq %r15 # + pushq %r14 # + pushq %r13 # + pushq %r12 # + pushq %rbp # + pushq %rbx # + + movq %rdi, %r14 # z + movq %rdx, %r12 # x + movl %esi, %ebp # z_size + + xorl %esi, %esi # j.76 + movq %r8, -16(%rsp) # u, u + movl %ecx, %ebx # x_size, x_size + movl %ecx, %r8d # x_size, blocks_of_8 + andl $-8, %r8d #, blocks_of_8 + testl %ecx, %ecx # x_size + je .L3 #, + mov %ecx, %eax # x_size, pretmp.71 + leal 1(%rbx), %r15d #, k.73 + salq $3, %rax #, + xorl %r13d, %r13d # j + movq %rax, -8(%rsp) #, pretmp.21 + .p2align 4,,10 + .p2align 3 +.L11: + mov %r13d, %eax # j, j + movq -16(%rsp), %rdi # u, y + leaq (%r14,%rax,8), %r11 #, z_j + xorl %r9d, %r9d # i + imulq (%r11), %rdi #* z_j, y + xorl %r10d, %r10d # carry + testl %r8d, %r8d # blocks_of_8 + je .L7 #, + .p2align 4,,10 + .p2align 3 +.LOOP_MUL_ADD: + mov %r9d, %ecx # i, i + addl $8, %r9d #, i + salq $3, %rcx #, D.2315 + leaq (%r11,%rcx), %rsi #, tmp130 + leaq (%r12,%rcx), %rcx #, tmp131 + + movq 8*0(%rcx), %rax + mulq %rdi # y + addq %r10, %rax # carry + adcq $0,%rdx + addq 8*0(%rsi), %rax + adcq $0,%rdx + movq %rdx,%r10 # carry + movq %rax, 8*0 (%rsi) + + movq 8*1(%rcx), %rax + mulq %rdi # y + addq %r10, %rax # carry + adcq $0,%rdx + addq 8*1(%rsi), %rax + adcq $0,%rdx + movq %rdx,%r10 # carry + movq %rax, 8*1 (%rsi) + + movq 8*2(%rcx), %rax + mulq %rdi # y + addq %r10, %rax # carry + adcq $0,%rdx + addq 8*2(%rsi), %rax + adcq $0,%rdx + movq %rdx,%r10 # carry + movq %rax, 8*2 (%rsi) + + movq 8*3(%rcx), %rax + mulq %rdi # y + addq %r10, %rax # carry + adcq $0,%rdx + addq 8*3(%rsi), %rax + adcq $0,%rdx + movq %rdx,%r10 # carry + movq %rax, 8*3 (%rsi) + + movq 8*4(%rcx), %rax + mulq %rdi # y + addq %r10, %rax # carry + adcq $0,%rdx + addq 8*4(%rsi), %rax + adcq $0,%rdx + movq %rdx,%r10 # carry + movq %rax, 8*4 (%rsi) + + movq 8*5(%rcx), %rax + mulq %rdi # y + addq %r10, %rax # carry + adcq $0,%rdx + addq 8*5(%rsi), %rax + adcq $0,%rdx + movq %rdx,%r10 # carry + movq %rax, 8*5 (%rsi) + + movq 8*6(%rcx), %rax + mulq %rdi # y + addq %r10, %rax # carry + adcq $0,%rdx + addq 8*6(%rsi), %rax + adcq $0,%rdx + movq %rdx,%r10 # carry + movq %rax, 8*6 (%rsi) + + movq 8*7(%rcx), %rax + mulq %rdi # y + addq %r10, %rax # carry + adcq $0,%rdx + addq 8*7(%rsi), %rax + adcq $0,%rdx + movq %rdx,%r10 # carry + movq %rax, 8*7 (%rsi) + + cmpl %r9d, %r8d # i, blocks_of_8 + jne .LOOP_MUL_ADD #, + cmpl %r8d, %ebx # blocks_of_8, x_size + je .L8 #, +.L7: + movl %r8d, %esi # blocks_of_8, i + .p2align 4,,10 + .p2align 3 +.L5: + mov %esi, %eax # i, i + movq %rdi, %rcx # y, b + leaq (%r11, %rax,8), %r9 #, D.2325 + incl %esi # i + movq (%r12, %rax,8), %rax #* x, tmp133 + + mulq %rcx # b + addq (%r9), %rax #* D.2325, a + adcq $0,%rdx # + addq %r10, %rax # carry, a + adcq $0,%rdx # + + cmpl %esi, %ebx # i, x_size + movq %rdx, %r10 #, carry + movq %rax, (%r9) # a,* D.2325 + jne .L5 #, +.L8: + movq -8(%rsp), %rdx # pretmp.21, + leaq (%r11,%rdx), %rax #, D.2332 + movq (%rax), %rcx #* D.2332, D.2333 + leaq (%r10,%rcx), %rdx #, z_sum + movq %rdx, (%rax) # z_sum,* D.2332 + cmpq %rdx, %rcx # z_sum, D.2333 + jbe .L9 #, + cmpl %ebp, %r15d # z_size, k.73 + je .L9 #, + movl %r15d, %ecx # k.73, k + jmp .L10 # + .p2align 4,,10 + .p2align 3 +.L31: + incl %ecx # k + cmpl %ecx, %ebp # k, z_size + .p2align 4,,4 + .p2align 3 + je .L9 #, +.L10: + mov %ecx, %edx # k, k + leaq (%r11,%rdx,8), %rdx #, D.2342 + movq (%rdx), %rax #* D.2342, tmp136 + incq %rax # D.2344 + movq %rax, (%rdx) # D.2344,* D.2342 + testq %rax, %rax # D.2344 + je .L31 #, +.L9: + incl %r13d # j + decl %ebp # z_size + cmpl %r13d, %ebx # j, x_size + jne .L11 #, + movl %ebx, %esi # x_size, j.76 +.L3: + leal (%rbx,%rbx), %eax #, tmp137 + mov %eax, %eax + leaq (%r14, %rax,8), %rdi #, D.2349 + cmpq $0, (%rdi) #,* D.2349 + jne .L12 #, + testl %ebx, %ebx # x_size + je .L12 #, + leal -1(%rbx), %ecx #, j + leal (%rsi,%rcx), %edx #, tmp141 + mov %ecx, %eax # j, j + movq (%r14,%rdx,8), %rbp #* z, + cmpq %rbp, (%r12, %rax,8) #,* x + jb .L12 #, + ja .L_EXIT #, + leal -2(%rsi,%rbx), %edx #, ivtmp.45 + jmp .L14 # + .p2align 4,,10 + .p2align 3 +.L15: + mov %edx, %eax # ivtmp.45, ivtmp.45 + decl %ecx # j + movq (%r14, %rax,8), %rsi #* z, D.2360 + mov %ecx, %eax # j, j + movq (%r12, %rax,8), %rax #* x, temp.68 + cmpq %rax, %rsi + ja .L12 #, + decl %edx # ivtmp.45 + cmpq %rax, %rsi + jb .L_EXIT #, +.L14: + testl %ecx, %ecx # j + jne .L15 #, +.L12: + xorl %ecx, %ecx # j + xorl %r10d, %r10d # carry + mov %ebx, %esi # x_size, pretmp.19 + testl %r8d, %r8d # blocks_of_8 + je .L17 #, + .p2align 4,,10 + .p2align 3 +.L22: + mov %ecx, %edx # j, D.2375 + addl $8, %ecx #, j + leaq (%rdx,%rsi), %rax #, tmp146 + leaq (%r12,%rdx,8), %rdx #, tmp150 + leaq (%r14, %rax,8), %rax #, tmp148 + + rorq %r10 # carry + + movq 8*0(%rdx), %r10 + sbbq %r10, 8*0(%rax) + + movq 8*1(%rdx), %r10 + sbbq %r10, 8*1(%rax) + + movq 8*2(%rdx), %r10 + sbbq %r10, 8*2(%rax) + + movq 8*3(%rdx), %r10 + sbbq %r10, 8*3(%rax) + + movq 8*4(%rdx), %r10 + sbbq %r10, 8*4(%rax) + + movq 8*5(%rdx), %r10 + sbbq %r10, 8*5(%rax) + + movq 8*6(%rdx), %r10 + sbbq %r10, 8*6(%rax) + + movq 8*7(%rdx), %r10 + sbbq %r10, 8*7(%rax) + + sbbq %r10,%r10 # carry + negq %r10 # carry + + cmpl %ecx, %r8d # j, blocks_of_8 + jne .L22 #, +.L17: + cmpl %r8d, %ebx # blocks_of_8, x_size + je .L19 #, + leal (%r8,%rbx), %r9d #, ivtmp.33 + movl %r8d, %esi # blocks_of_8, j + .p2align 4,,10 + .p2align 3 +.L20: + mov %r9d, %eax # ivtmp.33, ivtmp.33 + mov %esi, %ecx # j, j + leaq (%r14, %rax,8), %rax #, D.2387 + incl %esi # j + movq (%rax), %rdx #* D.2387, tmp153 + incl %r9d # ivtmp.33 + + rorq %r10 # carry + sbbq (%r12,%rcx,8),%rdx #* x, x + sbbq %r10,%r10 # carry + negq %r10 # carry + + cmpl %esi, %ebx # j, x_size + movq %rdx, (%rax) # x,* D.2387 + jne .L20 #, +.L19: + testq %r10, %r10 # carry + je .L_EXIT #, + decq (%rdi) #* D.2349 +.L_EXIT: + popq %rbx # + popq %rbp # + popq %r12 # + popq %r13 # + popq %r14 # + popq %r15 # +END_FUNCTION(bigint_monty_redc) + + +#if 0 + #define Z_ARR ARG_1 // rdi +#define Z_SIZE ARG_2_32 // esi +// X_ARR is ARG_3 == rdx, moved b/c needed for multiply +#define X_SIZE ARG_4_32 // ecx +#define U ARG_5 // r8 + +/* + We need all arguments for a while (we can reuse U eventually) + So only temp registers are + TEMP_1 %r10 + TEMP_2 %r11 + TEMP_3 = ARG_6 = %r9 + void return, so also + R0 %rax (aka TEMP_9) + is free (but needed for multiply) + + Can push: + %rbx (base pointer, callee saved) + %rpb (frame pointer, callee saved) + %r12-%r15 (callee saved) + + Can push base/frame pointers since this is a leaf function + and does not reference any data. +*/ + + push %r12 + push %r13 + push %r14 + push %r15 + +#define LOOP_CTR_I %r12 +#define LOOP_CTR_J %r13 + +#define CARRY TEMP_1 +#define Z_WORD TEMP_2 +#define X_ARR TEMP_3 +#define MUL_LO %rax +#define MUL_HI %rdx + + ASSIGN(X_ARR, ARG_3) + + /* + ZEROIZE(CARRY) + + ASSIGN(LOOP_CTR, X_SIZE) + + JUMP_IF_ZERO(LOOP_CTR, .L_MULADD_DONE) + JUMP_IF_LT(LOOP_CTR, 8, .LOOP_MULADD1) + +#define MULADD_OP(N) \ + ASSIGN(MUL_LO, ARRAY8(X_ARR, N)) ; \ + ASSIGN(Z_WORD, ARRAY8(Z_ARR, N)) ; \ + MUL(Y) ; \ + ADD(Z_WORD, CARRY) ; \ + ASSIGN(CARRY, MUL_HI) ; \ + ADD_LAST_CARRY(CARRY) ; \ + ADD(Z_WORD, MUL_LO) ; \ + ADD_LAST_CARRY(CARRY) ; \ + ASSIGN(ARRAY8(Z_ARR, N), Z_WORD) + +ALIGN +.LOOP_MULADD8: + MULADD_OP(0) + MULADD_OP(1) + MULADD_OP(2) + MULADD_OP(3) + MULADD_OP(4) + MULADD_OP(5) + MULADD_OP(6) + MULADD_OP(7) + + SUB_IMM(LOOP_CTR, 8) + ADD_IMM(Z_ARR, 64) + ADD_IMM(X_ARR, 64) + cmp IMM(8), LOOP_CTR + jge .LOOP_MULADD8 + + JUMP_IF_ZERO(LOOP_CTR, .L_MULADD_DONE) + +ALIGN +.LOOP_MULADD1: + MULADD_OP(0) + + SUB_IMM(LOOP_CTR, 1) + ADD_IMM(Z_ARR, 8) + ADD_IMM(X_ARR, 8) + + cmp IMM(0), LOOP_CTR + jne .LOOP_MULADD1 +*/ + + pop %r15 + pop %r14 + pop %r13 + pop %r12 +#endif diff --git a/src/math/bigint/monty_generic/info.txt b/src/math/bigint/monty_generic/info.txt new file mode 100644 index 000000000..6f5f0e722 --- /dev/null +++ b/src/math/bigint/monty_generic/info.txt @@ -0,0 +1,7 @@ +realname "Montgomery Reduction" + +load_on dep + +<add> +mp_monty.cpp +</add> diff --git a/src/math/bigint/monty_generic/mp_monty.cpp b/src/math/bigint/monty_generic/mp_monty.cpp new file mode 100644 index 000000000..c162bfd4f --- /dev/null +++ b/src/math/bigint/monty_generic/mp_monty.cpp @@ -0,0 +1,76 @@ +/************************************************* +* Montgomery Reduction Source File * +* (C) 1999-2008 Jack Lloyd * +* 2006 Luca Piccarreta * +*************************************************/ + +#include <botan/mp_core.h> +#include <botan/mp_asm.h> +#include <botan/mp_asmi.h> + +namespace Botan { + +extern "C" { + +/************************************************* +* Montgomery Reduction Algorithm * +*************************************************/ +void bigint_monty_redc(word z[], u32bit z_size, + const word x[], u32bit x_size, word u) + { + const u32bit blocks_of_8 = x_size - (x_size % 8); + + for(u32bit i = 0; i != x_size; ++i) + { + word* z_i = z + i; + + const word y = z_i[0] * u; + + word carry = 0; + + for(u32bit j = 0; j != blocks_of_8; j += 8) + carry = word8_madd3(z_i + j, x + j, y, carry); + + for(u32bit j = blocks_of_8; j != x_size; ++j) + z_i[j] = word_madd3(x[j], y, z_i[j], &carry); + + word z_sum = z_i[x_size] + carry; + carry = (z_sum < z_i[x_size]); + z_i[x_size] = z_sum; + + for(u32bit j = x_size + 1; carry && j != z_size - i; ++j) + { + ++z_i[j]; + carry = !z_i[j]; + } + } + + // Check if z[x_size...x_size+1] >= x[0...x_size] using bigint_cmp (inlined) + if(!z[x_size + x_size]) + { + for(u32bit i = x_size; i > 0; --i) + { + if(z[x_size + i - 1] > x[i-1]) + break; + + if(z[x_size + i - 1] < x[i-1]) + return; + } + } + + // If the compare above is true, subtract using bigint_sub2 (inlined) + word carry = 0; + + for(u32bit i = 0; i != blocks_of_8; i += 8) + carry = word8_sub2(z + x_size + i, x + i, carry); + + for(u32bit i = blocks_of_8; i != x_size; ++i) + z[x_size + i] = word_sub(z[x_size + i], x[i], &carry); + + if(carry) + --z[x_size+x_size]; + } + +} + +} diff --git a/src/math/bigint/mp_amd64/info.txt b/src/math/bigint/mp_amd64/info.txt new file mode 100644 index 000000000..84a5bcf53 --- /dev/null +++ b/src/math/bigint/mp_amd64/info.txt @@ -0,0 +1,19 @@ +realname "MPI Core (x86-64)" + +mp_bits 64 + +load_on dep + +<add> +mp_asm.h +mp_asmi.h +</add> + +<arch> +amd64 +</arch> + +<cc> +gcc +icc +</cc> diff --git a/src/math/bigint/mp_amd64/mp_asm.h b/src/math/bigint/mp_amd64/mp_asm.h new file mode 100644 index 000000000..eca7bae6c --- /dev/null +++ b/src/math/bigint/mp_amd64/mp_asm.h @@ -0,0 +1,67 @@ +/************************************************* +* Lowest Level MPI Algorithms Header File * +* (C) 1999-2008 Jack Lloyd * +* 2006 Luca Piccarreta * +*************************************************/ + +#ifndef BOTAN_MP_ASM_H__ +#define BOTAN_MP_ASM_H__ + +#include <botan/mp_types.h> + +#if (BOTAN_MP_WORD_BITS != 64) + #error The mp_amd64 module requires that BOTAN_MP_WORD_BITS == 64 +#endif + +namespace Botan { + +extern "C" { + +/************************************************* +* Helper Macros for amd64 Assembly * +*************************************************/ +#define ASM(x) x "\n\t" + +/************************************************* +* Word Multiply * +*************************************************/ +inline word word_madd2(word a, word b, word* c) + { + asm( + ASM("mulq %[b]") + ASM("addq %[c],%[a]") + ASM("adcq $0,%[carry]") + + : [a]"=a"(a), [b]"=rm"(b), [carry]"=&d"(*c) + : "0"(a), "1"(b), [c]"g"(*c) : "cc"); + + return a; + } + +/************************************************* +* Word Multiply/Add * +*************************************************/ +inline word word_madd3(word a, word b, word c, word* d) + { + asm( + ASM("mulq %[b]") + + ASM("addq %[c],%[a]") + ASM("adcq $0,%[carry]") + + ASM("addq %[d],%[a]") + ASM("adcq $0,%[carry]") + + : [a]"=a"(a), [b]"=rm"(b), [carry]"=&d"(*d) + : "0"(a), "1"(b), [c]"g"(c), [d]"g"(*d) : "cc"); + + return a; + } + +#undef ASM + +} + +} + +#endif diff --git a/src/math/bigint/mp_amd64/mp_asmi.h b/src/math/bigint/mp_amd64/mp_asmi.h new file mode 100644 index 000000000..16632a38d --- /dev/null +++ b/src/math/bigint/mp_amd64/mp_asmi.h @@ -0,0 +1,233 @@ +/************************************************* +* Lowest Level MPI Algorithms Header File * +* (C) 1999-2007 Jack Lloyd * +* 2006 Luca Piccarreta * +*************************************************/ + +#ifndef BOTAN_MP_ASM_INTERNAL_H__ +#define BOTAN_MP_ASM_INTERNAL_H__ + +#include <botan/mp_asm.h> + +namespace Botan { + +extern "C" { + +/************************************************* +* Helper Macros for amd64 Assembly * +*************************************************/ +#ifndef ASM + #define ASM(x) x "\n\t" +#endif + +#define ADDSUB2_OP(OPERATION, INDEX) \ + ASM("movq 8*" #INDEX "(%[y]), %[carry]") \ + ASM(OPERATION " %[carry], 8*" #INDEX "(%[x])") \ + +#define ADDSUB3_OP(OPERATION, INDEX) \ + ASM("movq 8*" #INDEX "(%[x]), %[carry]") \ + ASM(OPERATION " 8*" #INDEX "(%[y]), %[carry]") \ + ASM("movq %[carry], 8*" #INDEX "(%[z])") \ + +#define LINMUL_OP(WRITE_TO, INDEX) \ + ASM("movq 8*" #INDEX "(%[x]),%%rax") \ + ASM("mulq %[y]") \ + ASM("addq %[carry],%%rax") \ + ASM("adcq $0,%%rdx") \ + ASM("movq %%rdx,%[carry]") \ + ASM("movq %%rax, 8*" #INDEX "(%[" WRITE_TO "])") + +#define MULADD_OP(IGNORED, INDEX) \ + ASM("movq 8*" #INDEX "(%[x]),%%rax") \ + ASM("mulq %[y]") \ + ASM("addq %[carry],%%rax") \ + ASM("adcq $0,%%rdx") \ + ASM("addq 8*" #INDEX "(%[z]),%%rax") \ + ASM("adcq $0,%%rdx") \ + ASM("movq %%rdx,%[carry]") \ + ASM("movq %%rax, 8*" #INDEX " (%[z])") + +#define DO_8_TIMES(MACRO, ARG) \ + MACRO(ARG, 0) \ + MACRO(ARG, 1) \ + MACRO(ARG, 2) \ + MACRO(ARG, 3) \ + MACRO(ARG, 4) \ + MACRO(ARG, 5) \ + MACRO(ARG, 6) \ + MACRO(ARG, 7) + +#define ADD_OR_SUBTRACT(CORE_CODE) \ + ASM("rorq %[carry]") \ + CORE_CODE \ + ASM("sbbq %[carry],%[carry]") \ + ASM("negq %[carry]") + +/************************************************* +* Word Addition * +*************************************************/ +inline word word_add(word x, word y, word* carry) + { + asm( + ADD_OR_SUBTRACT(ASM("adcq %[y],%[x]")) + : [x]"=r"(x), [carry]"=r"(*carry) + : "0"(x), [y]"rm"(y), "1"(*carry) + : "cc"); + return x; + } + +/************************************************* +* Eight Word Block Addition, Two Argument * +*************************************************/ +inline word word8_add2(word x[8], const word y[8], word carry) + { + asm( + ADD_OR_SUBTRACT(DO_8_TIMES(ADDSUB2_OP, "adcq")) + : [carry]"=r"(carry) + : [x]"r"(x), [y]"r"(y), "0"(carry) + : "cc", "memory"); + return carry; + } + +/************************************************* +* Eight Word Block Addition, Three Argument * +*************************************************/ +inline word word8_add3(word z[8], const word x[8], const word y[8], word carry) + { + asm( + ADD_OR_SUBTRACT(DO_8_TIMES(ADDSUB3_OP, "adcq")) + : [carry]"=r"(carry) + : [x]"r"(x), [y]"r"(y), [z]"r"(z), "0"(carry) + : "cc", "memory"); + return carry; + } + +/************************************************* +* Word Subtraction * +*************************************************/ +inline word word_sub(word x, word y, word* carry) + { + asm( + ADD_OR_SUBTRACT(ASM("sbbq %[y],%[x]")) + : [x]"=r"(x), [carry]"=r"(*carry) + : "0"(x), [y]"rm"(y), "1"(*carry) + : "cc"); + return x; + } + +/************************************************* +* Eight Word Block Subtraction, Two Argument * +*************************************************/ +inline word word8_sub2(word x[8], const word y[8], word carry) + { + asm( + ADD_OR_SUBTRACT(DO_8_TIMES(ADDSUB2_OP, "sbbq")) + : [carry]"=r"(carry) + : [x]"r"(x), [y]"r"(y), "0"(carry) + : "cc", "memory"); + return carry; + } + +/************************************************* +* Eight Word Block Subtraction, Three Argument * +*************************************************/ +inline word word8_sub3(word z[8], const word x[8], const word y[8], word carry) + { + asm( + ADD_OR_SUBTRACT(DO_8_TIMES(ADDSUB3_OP, "sbbq")) + : [carry]"=r"(carry) + : [x]"r"(x), [y]"r"(y), [z]"r"(z), "0"(carry) + : "cc", "memory"); + return carry; + } + +/************************************************* +* Eight Word Block Linear Multiplication * +*************************************************/ +inline word word8_linmul2(word x[8], word y, word carry) + { + asm( + DO_8_TIMES(LINMUL_OP, "x") + : [carry]"=r"(carry) + : [x]"r"(x), [y]"rm"(y), "0"(carry) + : "cc", "%rax", "%rdx"); + return carry; + } + +/************************************************* +* Eight Word Block Linear Multiplication * +*************************************************/ +inline word word8_linmul3(word z[8], const word x[8], word y, word carry) + { + asm( + DO_8_TIMES(LINMUL_OP, "z") + : [carry]"=r"(carry) + : [z]"r"(z), [x]"r"(x), [y]"rm"(y), "0"(carry) + : "cc", "%rax", "%rdx"); + return carry; + } + +/************************************************* +* Eight Word Block Multiply/Add * +*************************************************/ +inline word word8_madd3(word z[8], const word x[8], word y, word carry) + { + asm( + DO_8_TIMES(MULADD_OP, "") + : [carry]"=r"(carry) + : [z]"r"(z), [x]"r"(x), [y]"rm"(y), "0"(carry) + : "cc", "%rax", "%rdx"); + return carry; + } + +/************************************************* +* Multiply-Add Accumulator * +*************************************************/ +inline void word3_muladd(word* w2, word* w1, word* w0, word x, word y) + { + asm( + ASM("mulq %[y]") + + ASM("addq %[x],%[w0]") + ASM("adcq %[y],%[w1]") + ASM("adcq $0,%[w2]") + + : [w0]"=r"(*w0), [w1]"=r"(*w1), [w2]"=r"(*w2) + : [x]"a"(x), [y]"d"(y), "0"(*w0), "1"(*w1), "2"(*w2) + : "cc"); + } + +/************************************************* +* Multiply-Add Accumulator * +*************************************************/ +inline void word3_muladd_2(word* w2, word* w1, word* w0, word x, word y) + { + asm( + ASM("mulq %[y]") + + ASM("addq %[x],%[w0]") + ASM("adcq %[y],%[w1]") + ASM("adcq $0,%[w2]") + + ASM("addq %[x],%[w0]") + ASM("adcq %[y],%[w1]") + ASM("adcq $0,%[w2]") + + : [w0]"=r"(*w0), [w1]"=r"(*w1), [w2]"=r"(*w2) + : [x]"a"(x), [y]"d"(y), "0"(*w0), "1"(*w1), "2"(*w2) + : "cc"); + } + + +#undef ASM +#undef DO_8_TIMES +#undef ADD_OR_SUBTRACT +#undef ADDSUB2_OP +#undef ADDSUB3_OP +#undef LINMUL_OP +#undef MULADD_OP + +} + +} +#endif diff --git a/src/math/bigint/mp_asm.cpp b/src/math/bigint/mp_asm.cpp new file mode 100644 index 000000000..e5d1fe0d6 --- /dev/null +++ b/src/math/bigint/mp_asm.cpp @@ -0,0 +1,177 @@ +/************************************************* +* Lowest Level MPI Algorithms Source File * +* (C) 1999-2008 Jack Lloyd * +* 2006 Luca Piccarreta * +*************************************************/ + +#include <botan/mp_asm.h> +#include <botan/mp_asmi.h> +#include <botan/mp_core.h> +#include <botan/mem_ops.h> + +namespace Botan { + +extern "C" { + +/************************************************* +* Two Operand Addition, No Carry * +*************************************************/ +word bigint_add2_nc(word x[], u32bit x_size, const word y[], u32bit y_size) + { + word carry = 0; + + const u32bit blocks = y_size - (y_size % 8); + + for(u32bit j = 0; j != blocks; j += 8) + carry = word8_add2(x + j, y + j, carry); + + for(u32bit j = blocks; j != y_size; ++j) + x[j] = word_add(x[j], y[j], &carry); + + if(!carry) + return 0; + + for(u32bit j = y_size; j != x_size; ++j) + if(++x[j]) + return 0; + + return 1; + } + +/************************************************* +* Three Operand Addition, No Carry * +*************************************************/ +word bigint_add3_nc(word z[], const word x[], u32bit x_size, + const word y[], u32bit y_size) + { + if(x_size < y_size) + { return bigint_add3_nc(z, y, y_size, x, x_size); } + + word carry = 0; + + const u32bit blocks = y_size - (y_size % 8); + + for(u32bit j = 0; j != blocks; j += 8) + carry = word8_add3(z + j, x + j, y + j, carry); + + for(u32bit j = blocks; j != y_size; ++j) + z[j] = word_add(x[j], y[j], &carry); + + for(u32bit j = y_size; j != x_size; ++j) + { + word x_j = x[j] + carry; + if(carry && x_j) + carry = 0; + z[j] = x_j; + } + + return carry; + } + +/************************************************* +* Two Operand Addition * +*************************************************/ +void bigint_add2(word x[], u32bit x_size, const word y[], u32bit y_size) + { + if(bigint_add2_nc(x, x_size, y, y_size)) + ++x[x_size]; + } + +/************************************************* +* Three Operand Addition * +*************************************************/ +void bigint_add3(word z[], const word x[], u32bit x_size, + const word y[], u32bit y_size) + { + if(bigint_add3_nc(z, x, x_size, y, y_size)) + ++z[(x_size > y_size ? x_size : y_size)]; + } + +/************************************************* +* Two Operand Subtraction * +*************************************************/ +void bigint_sub2(word x[], u32bit x_size, const word y[], u32bit y_size) + { + word carry = 0; + + const u32bit blocks = y_size - (y_size % 8); + + for(u32bit j = 0; j != blocks; j += 8) + carry = word8_sub2(x + j, y + j, carry); + + for(u32bit j = blocks; j != y_size; ++j) + x[j] = word_sub(x[j], y[j], &carry); + + if(!carry) return; + + for(u32bit j = y_size; j != x_size; ++j) + { + --x[j]; + if(x[j] != MP_WORD_MAX) return; + } + } + +/************************************************* +* Three Operand Subtraction * +*************************************************/ +void bigint_sub3(word z[], const word x[], u32bit x_size, + const word y[], u32bit y_size) + { + word carry = 0; + + const u32bit blocks = y_size - (y_size % 8); + + for(u32bit j = 0; j != blocks; j += 8) + carry = word8_sub3(z + j, x + j, y + j, carry); + + for(u32bit j = blocks; j != y_size; ++j) + z[j] = word_sub(x[j], y[j], &carry); + + for(u32bit j = y_size; j != x_size; ++j) + { + word x_j = x[j] - carry; + if(carry && x_j != MP_WORD_MAX) + carry = 0; + z[j] = x_j; + } + } + +/************************************************* +* Two Operand Linear Multiply * +*************************************************/ +void bigint_linmul2(word x[], u32bit x_size, word y) + { + const u32bit blocks = x_size - (x_size % 8); + + word carry = 0; + + for(u32bit j = 0; j != blocks; j += 8) + carry = word8_linmul2(x + j, y, carry); + + for(u32bit j = blocks; j != x_size; ++j) + x[j] = word_madd2(x[j], y, &carry); + + x[x_size] = carry; + } + +/************************************************* +* Three Operand Linear Multiply * +*************************************************/ +void bigint_linmul3(word z[], const word x[], u32bit x_size, word y) + { + const u32bit blocks = x_size - (x_size % 8); + + word carry = 0; + + for(u32bit j = 0; j != blocks; j += 8) + carry = word8_linmul3(z + j, x + j, y, carry); + + for(u32bit j = blocks; j != x_size; ++j) + z[j] = word_madd2(x[j], y, &carry); + + z[x_size] = carry; + } + +} + +} diff --git a/src/math/bigint/mp_asm64/info.txt b/src/math/bigint/mp_asm64/info.txt new file mode 100644 index 000000000..85a4391b9 --- /dev/null +++ b/src/math/bigint/mp_asm64/info.txt @@ -0,0 +1,26 @@ +realname "MPI Core (Alpha/IA-64/MIPS64/PowerPC-64/SPARC64)" + +mp_bits 64 + +load_on dep + +<add> +mp_asm.h +mp_generic:mp_asmi.h +</add> + +<arch> +alpha +ia64 +mips64 +ppc64 +sparc64 +</arch> + +# The inline asm only works with gcc, but it looks like (at least on +# UltraSPARC), using 64-bit words and the sythensized multiply is a 5 to 25% +# win, so it's probably worth using elsewhere. +<cc> +gcc +sunwspro +</cc> diff --git a/src/math/bigint/mp_asm64/mp_asm.h b/src/math/bigint/mp_asm64/mp_asm.h new file mode 100644 index 000000000..e455b3616 --- /dev/null +++ b/src/math/bigint/mp_asm64/mp_asm.h @@ -0,0 +1,110 @@ +/************************************************* +* MPI Multiply-Add Core Header File * +* (C) 1999-2007 Jack Lloyd * +*************************************************/ + +#ifndef BOTAN_MP_MADD_H__ +#define BOTAN_MP_MADD_H__ + +#include <botan/mp_types.h> + +namespace Botan { + +#if (BOTAN_MP_WORD_BITS != 64) + #error The mp_asm64 module requires that BOTAN_MP_WORD_BITS == 64 +#endif + +#if defined(BOTAN_TARGET_ARCH_IS_ALPHA) + +#define BOTAN_WORD_MUL(a,b,z1,z0) do { \ + asm("umulh %1,%2,%0" : "=r" (z0) : "r" (a), "r" (b)); \ + z1 = a * b; \ +} while(0); + +#elif defined(BOTAN_TARGET_ARCH_IS_IA64) + +#define BOTAN_WORD_MUL(a,b,z1,z0) do { \ + asm("xmpy.hu %0=%1,%2" : "=f" (z0) : "f" (a), "f" (b)); \ + z1 = a * b; \ +} while(0); + +#elif defined(BOTAN_TARGET_ARCH_IS_PPC64) + +#define BOTAN_WORD_MUL(a,b,z1,z0) do { \ + asm("mulhdu %0,%1,%2" : "=r" (z0) : "r" (a), "r" (b) : "cc"); \ + z1 = a * b; \ +} while(0); + +#elif defined(BOTAN_TARGET_ARCH_IS_MIPS64) + +#define BOTAN_WORD_MUL(a,b,z1,z0) do { \ + asm("dmultu %2,%3" : "=h" (z0), "=l" (z1) : "r" (a), "r" (b)); \ +} while(0); + +#else + +// Do a 64x64->128 multiply using four 64x64->64 multiplies +// plus some adds and shifts. Last resort for CPUs like UltraSPARC, +// with 64-bit registers/ALU, but no 64x64->128 multiply. +inline void bigint_2word_mul(word a, word b, word* z1, word* z0) + { + const u32bit MP_HWORD_BITS = MP_WORD_BITS / 2; + const word MP_HWORD_MASK = ((word)1 << MP_HWORD_BITS) - 1; + + const word a_hi = (a >> MP_HWORD_BITS); + const word a_lo = (a & MP_HWORD_MASK); + const word b_hi = (b >> MP_HWORD_BITS); + const word b_lo = (b & MP_HWORD_MASK); + + word x0 = a_hi * b_hi; + word x1 = a_lo * b_hi; + word x2 = a_hi * b_lo; + word x3 = a_lo * b_lo; + + x2 += x3 >> (MP_HWORD_BITS); + x2 += x1; + if(x2 < x1) + x0 += ((word)1 << MP_HWORD_BITS); + + *z0 = x0 + (x2 >> MP_HWORD_BITS); + *z1 = ((x2 & MP_HWORD_MASK) << MP_HWORD_BITS) + (x3 & MP_HWORD_MASK); + } + +#define BOTAN_WORD_MUL(a,b,z1,z0) bigint_2word_mul(a, b, &z1, &z0) + +#endif + +/************************************************* +* Word Multiply/Add * +*************************************************/ +inline word word_madd2(word a, word b, word* c) + { + word z0 = 0, z1 = 0; + + BOTAN_WORD_MUL(a, b, z1, z0); + + z1 += *c; if(z1 < *c) z0++; + + *c = z0; + return z1; + } + +/************************************************* +* Word Multiply/Add * +*************************************************/ +inline word word_madd3(word a, word b, word c, word* d) + { + word z0 = 0, z1 = 0; + + BOTAN_WORD_MUL(a, b, z1, z0); + + z1 += c; if(z1 < c) z0++; + z1 += *d; if(z1 < *d) z0++; + + *d = z0; + return z1; + } + +} + +#endif diff --git a/src/math/bigint/mp_comba.cpp b/src/math/bigint/mp_comba.cpp new file mode 100644 index 000000000..c7a9c964c --- /dev/null +++ b/src/math/bigint/mp_comba.cpp @@ -0,0 +1,918 @@ +/************************************************* +* Comba Multiplication and Squaring Source File * +* (C) 1999-2007 Jack Lloyd * +*************************************************/ + +#include <botan/mp_core.h> +#include <botan/mp_asmi.h> + +namespace Botan { + +extern "C" { + +/************************************************* +* Comba 4x4 Squaring * +*************************************************/ +void bigint_comba_sqr4(word z[8], const word x[4]) + { + word w2 = 0, w1 = 0, w0 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], x[0]); + z[0] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[1]); + z[1] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[2]); + word3_muladd(&w2, &w1, &w0, x[1], x[1]); + z[2] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[3]); + word3_muladd_2(&w2, &w1, &w0, x[1], x[2]); + z[3] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[1], x[3]); + word3_muladd(&w2, &w1, &w0, x[2], x[2]); + z[4] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[2], x[3]); + z[5] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[3], x[3]); + z[6] = w0; + z[7] = w1; + } + +/************************************************* +* Comba 4x4 Multiplication * +*************************************************/ +void bigint_comba_mul4(word z[8], const word x[4], const word y[4]) + { + word w2 = 0, w1 = 0, w0 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[0]); + z[0] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[1]); + word3_muladd(&w2, &w1, &w0, x[1], y[0]); + z[1] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[2]); + word3_muladd(&w2, &w1, &w0, x[1], y[1]); + word3_muladd(&w2, &w1, &w0, x[2], y[0]); + z[2] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[3]); + word3_muladd(&w2, &w1, &w0, x[1], y[2]); + word3_muladd(&w2, &w1, &w0, x[2], y[1]); + word3_muladd(&w2, &w1, &w0, x[3], y[0]); + z[3] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[1], y[3]); + word3_muladd(&w2, &w1, &w0, x[2], y[2]); + word3_muladd(&w2, &w1, &w0, x[3], y[1]); + z[4] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[2], y[3]); + word3_muladd(&w2, &w1, &w0, x[3], y[2]); + z[5] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[3], y[3]); + z[6] = w0; + z[7] = w1; + } + +/************************************************* +* Comba 6x6 Squaring * +*************************************************/ +void bigint_comba_sqr6(word z[12], const word x[6]) + { + word w2 = 0, w1 = 0, w0 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], x[0]); + z[0] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[1]); + z[1] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[2]); + word3_muladd(&w2, &w1, &w0, x[1], x[1]); + z[2] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[3]); + word3_muladd_2(&w2, &w1, &w0, x[1], x[2]); + z[3] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[4]); + word3_muladd_2(&w2, &w1, &w0, x[1], x[3]); + word3_muladd(&w2, &w1, &w0, x[2], x[2]); + z[4] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[5]); + word3_muladd_2(&w2, &w1, &w0, x[1], x[4]); + word3_muladd_2(&w2, &w1, &w0, x[2], x[3]); + z[5] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[1], x[5]); + word3_muladd_2(&w2, &w1, &w0, x[2], x[4]); + word3_muladd(&w2, &w1, &w0, x[3], x[3]); + z[6] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[2], x[5]); + word3_muladd_2(&w2, &w1, &w0, x[3], x[4]); + z[7] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[3], x[5]); + word3_muladd(&w2, &w1, &w0, x[4], x[4]); + z[8] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[4], x[5]); + z[9] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[5], x[5]); + z[10] = w0; + z[11] = w1; + } + +/************************************************* +* Comba 6x6 Multiplication * +*************************************************/ +void bigint_comba_mul6(word z[12], const word x[6], const word y[6]) + { + word w2 = 0, w1 = 0, w0 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[0]); + z[0] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[1]); + word3_muladd(&w2, &w1, &w0, x[1], y[0]); + z[1] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[2]); + word3_muladd(&w2, &w1, &w0, x[1], y[1]); + word3_muladd(&w2, &w1, &w0, x[2], y[0]); + z[2] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[3]); + word3_muladd(&w2, &w1, &w0, x[1], y[2]); + word3_muladd(&w2, &w1, &w0, x[2], y[1]); + word3_muladd(&w2, &w1, &w0, x[3], y[0]); + z[3] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[4]); + word3_muladd(&w2, &w1, &w0, x[1], y[3]); + word3_muladd(&w2, &w1, &w0, x[2], y[2]); + word3_muladd(&w2, &w1, &w0, x[3], y[1]); + word3_muladd(&w2, &w1, &w0, x[4], y[0]); + z[4] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[5]); + word3_muladd(&w2, &w1, &w0, x[1], y[4]); + word3_muladd(&w2, &w1, &w0, x[2], y[3]); + word3_muladd(&w2, &w1, &w0, x[3], y[2]); + word3_muladd(&w2, &w1, &w0, x[4], y[1]); + word3_muladd(&w2, &w1, &w0, x[5], y[0]); + z[5] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[1], y[5]); + word3_muladd(&w2, &w1, &w0, x[2], y[4]); + word3_muladd(&w2, &w1, &w0, x[3], y[3]); + word3_muladd(&w2, &w1, &w0, x[4], y[2]); + word3_muladd(&w2, &w1, &w0, x[5], y[1]); + z[6] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[2], y[5]); + word3_muladd(&w2, &w1, &w0, x[3], y[4]); + word3_muladd(&w2, &w1, &w0, x[4], y[3]); + word3_muladd(&w2, &w1, &w0, x[5], y[2]); + z[7] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[3], y[5]); + word3_muladd(&w2, &w1, &w0, x[4], y[4]); + word3_muladd(&w2, &w1, &w0, x[5], y[3]); + z[8] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[4], y[5]); + word3_muladd(&w2, &w1, &w0, x[5], y[4]); + z[9] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[5], y[5]); + z[10] = w0; + z[11] = w1; + } + +/************************************************* +* Comba 8x8 Squaring * +*************************************************/ +void bigint_comba_sqr8(word z[16], const word x[8]) + { + word w2 = 0, w1 = 0, w0 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], x[0]); + z[0] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[1]); + z[1] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[2]); + word3_muladd(&w2, &w1, &w0, x[1], x[1]); + z[2] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[3]); + word3_muladd_2(&w2, &w1, &w0, x[1], x[2]); + z[3] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[4]); + word3_muladd_2(&w2, &w1, &w0, x[1], x[3]); + word3_muladd(&w2, &w1, &w0, x[2], x[2]); + z[4] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[5]); + word3_muladd_2(&w2, &w1, &w0, x[1], x[4]); + word3_muladd_2(&w2, &w1, &w0, x[2], x[3]); + z[5] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[6]); + word3_muladd_2(&w2, &w1, &w0, x[1], x[5]); + word3_muladd_2(&w2, &w1, &w0, x[2], x[4]); + word3_muladd(&w2, &w1, &w0, x[3], x[3]); + z[6] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[0], x[7]); + word3_muladd_2(&w2, &w1, &w0, x[1], x[6]); + word3_muladd_2(&w2, &w1, &w0, x[2], x[5]); + word3_muladd_2(&w2, &w1, &w0, x[3], x[4]); + z[7] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[1], x[7]); + word3_muladd_2(&w2, &w1, &w0, x[2], x[6]); + word3_muladd_2(&w2, &w1, &w0, x[3], x[5]); + word3_muladd(&w2, &w1, &w0, x[4], x[4]); + z[8] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[2], x[7]); + word3_muladd_2(&w2, &w1, &w0, x[3], x[6]); + word3_muladd_2(&w2, &w1, &w0, x[4], x[5]); + z[9] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[3], x[7]); + word3_muladd_2(&w2, &w1, &w0, x[4], x[6]); + word3_muladd(&w2, &w1, &w0, x[5], x[5]); + z[10] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[4], x[7]); + word3_muladd_2(&w2, &w1, &w0, x[5], x[6]); + z[11] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[5], x[7]); + word3_muladd(&w2, &w1, &w0, x[6], x[6]); + z[12] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[6], x[7]); + z[13] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[7], x[7]); + z[14] = w0; + z[15] = w1; + } + +/************************************************* +* Comba 8x8 Multiplication * +*************************************************/ +void bigint_comba_mul8(word z[16], const word x[8], const word y[8]) + { + word w2 = 0, w1 = 0, w0 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[0]); + z[0] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[1]); + word3_muladd(&w2, &w1, &w0, x[1], y[0]); + z[1] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[2]); + word3_muladd(&w2, &w1, &w0, x[1], y[1]); + word3_muladd(&w2, &w1, &w0, x[2], y[0]); + z[2] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[3]); + word3_muladd(&w2, &w1, &w0, x[1], y[2]); + word3_muladd(&w2, &w1, &w0, x[2], y[1]); + word3_muladd(&w2, &w1, &w0, x[3], y[0]); + z[3] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[4]); + word3_muladd(&w2, &w1, &w0, x[1], y[3]); + word3_muladd(&w2, &w1, &w0, x[2], y[2]); + word3_muladd(&w2, &w1, &w0, x[3], y[1]); + word3_muladd(&w2, &w1, &w0, x[4], y[0]); + z[4] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[5]); + word3_muladd(&w2, &w1, &w0, x[1], y[4]); + word3_muladd(&w2, &w1, &w0, x[2], y[3]); + word3_muladd(&w2, &w1, &w0, x[3], y[2]); + word3_muladd(&w2, &w1, &w0, x[4], y[1]); + word3_muladd(&w2, &w1, &w0, x[5], y[0]); + z[5] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[6]); + word3_muladd(&w2, &w1, &w0, x[1], y[5]); + word3_muladd(&w2, &w1, &w0, x[2], y[4]); + word3_muladd(&w2, &w1, &w0, x[3], y[3]); + word3_muladd(&w2, &w1, &w0, x[4], y[2]); + word3_muladd(&w2, &w1, &w0, x[5], y[1]); + word3_muladd(&w2, &w1, &w0, x[6], y[0]); + z[6] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[0], y[7]); + word3_muladd(&w2, &w1, &w0, x[1], y[6]); + word3_muladd(&w2, &w1, &w0, x[2], y[5]); + word3_muladd(&w2, &w1, &w0, x[3], y[4]); + word3_muladd(&w2, &w1, &w0, x[4], y[3]); + word3_muladd(&w2, &w1, &w0, x[5], y[2]); + word3_muladd(&w2, &w1, &w0, x[6], y[1]); + word3_muladd(&w2, &w1, &w0, x[7], y[0]); + z[7] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[1], y[7]); + word3_muladd(&w2, &w1, &w0, x[2], y[6]); + word3_muladd(&w2, &w1, &w0, x[3], y[5]); + word3_muladd(&w2, &w1, &w0, x[4], y[4]); + word3_muladd(&w2, &w1, &w0, x[5], y[3]); + word3_muladd(&w2, &w1, &w0, x[6], y[2]); + word3_muladd(&w2, &w1, &w0, x[7], y[1]); + z[8] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[2], y[7]); + word3_muladd(&w2, &w1, &w0, x[3], y[6]); + word3_muladd(&w2, &w1, &w0, x[4], y[5]); + word3_muladd(&w2, &w1, &w0, x[5], y[4]); + word3_muladd(&w2, &w1, &w0, x[6], y[3]); + word3_muladd(&w2, &w1, &w0, x[7], y[2]); + z[9] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[3], y[7]); + word3_muladd(&w2, &w1, &w0, x[4], y[6]); + word3_muladd(&w2, &w1, &w0, x[5], y[5]); + word3_muladd(&w2, &w1, &w0, x[6], y[4]); + word3_muladd(&w2, &w1, &w0, x[7], y[3]); + z[10] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[4], y[7]); + word3_muladd(&w2, &w1, &w0, x[5], y[6]); + word3_muladd(&w2, &w1, &w0, x[6], y[5]); + word3_muladd(&w2, &w1, &w0, x[7], y[4]); + z[11] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[5], y[7]); + word3_muladd(&w2, &w1, &w0, x[6], y[6]); + word3_muladd(&w2, &w1, &w0, x[7], y[5]); + z[12] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[6], y[7]); + word3_muladd(&w2, &w1, &w0, x[7], y[6]); + z[13] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[7], y[7]); + z[14] = w0; + z[15] = w1; + } + +/************************************************* +* Comba 16x16 Squaring * +*************************************************/ +void bigint_comba_sqr16(word z[32], const word x[16]) + { + word w2 = 0, w1 = 0, w0 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], x[ 0]); + z[ 0] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[ 1]); + z[ 1] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[ 2]); + word3_muladd(&w2, &w1, &w0, x[ 1], x[ 1]); + z[ 2] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[ 3]); + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[ 2]); + z[ 3] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[ 4]); + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[ 3]); + word3_muladd(&w2, &w1, &w0, x[ 2], x[ 2]); + z[ 4] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[ 5]); + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[ 4]); + word3_muladd_2(&w2, &w1, &w0, x[ 2], x[ 3]); + z[ 5] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[ 6]); + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[ 5]); + word3_muladd_2(&w2, &w1, &w0, x[ 2], x[ 4]); + word3_muladd(&w2, &w1, &w0, x[ 3], x[ 3]); + z[ 6] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[ 7]); + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[ 6]); + word3_muladd_2(&w2, &w1, &w0, x[ 2], x[ 5]); + word3_muladd_2(&w2, &w1, &w0, x[ 3], x[ 4]); + z[ 7] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[ 8]); + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[ 7]); + word3_muladd_2(&w2, &w1, &w0, x[ 2], x[ 6]); + word3_muladd_2(&w2, &w1, &w0, x[ 3], x[ 5]); + word3_muladd(&w2, &w1, &w0, x[ 4], x[ 4]); + z[ 8] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[ 9]); + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[ 8]); + word3_muladd_2(&w2, &w1, &w0, x[ 2], x[ 7]); + word3_muladd_2(&w2, &w1, &w0, x[ 3], x[ 6]); + word3_muladd_2(&w2, &w1, &w0, x[ 4], x[ 5]); + z[ 9] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[10]); + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[ 9]); + word3_muladd_2(&w2, &w1, &w0, x[ 2], x[ 8]); + word3_muladd_2(&w2, &w1, &w0, x[ 3], x[ 7]); + word3_muladd_2(&w2, &w1, &w0, x[ 4], x[ 6]); + word3_muladd(&w2, &w1, &w0, x[ 5], x[ 5]); + z[10] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[11]); + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[10]); + word3_muladd_2(&w2, &w1, &w0, x[ 2], x[ 9]); + word3_muladd_2(&w2, &w1, &w0, x[ 3], x[ 8]); + word3_muladd_2(&w2, &w1, &w0, x[ 4], x[ 7]); + word3_muladd_2(&w2, &w1, &w0, x[ 5], x[ 6]); + z[11] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[12]); + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[11]); + word3_muladd_2(&w2, &w1, &w0, x[ 2], x[10]); + word3_muladd_2(&w2, &w1, &w0, x[ 3], x[ 9]); + word3_muladd_2(&w2, &w1, &w0, x[ 4], x[ 8]); + word3_muladd_2(&w2, &w1, &w0, x[ 5], x[ 7]); + word3_muladd(&w2, &w1, &w0, x[ 6], x[ 6]); + z[12] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[13]); + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[12]); + word3_muladd_2(&w2, &w1, &w0, x[ 2], x[11]); + word3_muladd_2(&w2, &w1, &w0, x[ 3], x[10]); + word3_muladd_2(&w2, &w1, &w0, x[ 4], x[ 9]); + word3_muladd_2(&w2, &w1, &w0, x[ 5], x[ 8]); + word3_muladd_2(&w2, &w1, &w0, x[ 6], x[ 7]); + z[13] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[14]); + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[13]); + word3_muladd_2(&w2, &w1, &w0, x[ 2], x[12]); + word3_muladd_2(&w2, &w1, &w0, x[ 3], x[11]); + word3_muladd_2(&w2, &w1, &w0, x[ 4], x[10]); + word3_muladd_2(&w2, &w1, &w0, x[ 5], x[ 9]); + word3_muladd_2(&w2, &w1, &w0, x[ 6], x[ 8]); + word3_muladd(&w2, &w1, &w0, x[ 7], x[ 7]); + z[14] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 0], x[15]); + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[14]); + word3_muladd_2(&w2, &w1, &w0, x[ 2], x[13]); + word3_muladd_2(&w2, &w1, &w0, x[ 3], x[12]); + word3_muladd_2(&w2, &w1, &w0, x[ 4], x[11]); + word3_muladd_2(&w2, &w1, &w0, x[ 5], x[10]); + word3_muladd_2(&w2, &w1, &w0, x[ 6], x[ 9]); + word3_muladd_2(&w2, &w1, &w0, x[ 7], x[ 8]); + z[15] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 1], x[15]); + word3_muladd_2(&w2, &w1, &w0, x[ 2], x[14]); + word3_muladd_2(&w2, &w1, &w0, x[ 3], x[13]); + word3_muladd_2(&w2, &w1, &w0, x[ 4], x[12]); + word3_muladd_2(&w2, &w1, &w0, x[ 5], x[11]); + word3_muladd_2(&w2, &w1, &w0, x[ 6], x[10]); + word3_muladd_2(&w2, &w1, &w0, x[ 7], x[ 9]); + word3_muladd(&w2, &w1, &w0, x[ 8], x[ 8]); + z[16] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 2], x[15]); + word3_muladd_2(&w2, &w1, &w0, x[ 3], x[14]); + word3_muladd_2(&w2, &w1, &w0, x[ 4], x[13]); + word3_muladd_2(&w2, &w1, &w0, x[ 5], x[12]); + word3_muladd_2(&w2, &w1, &w0, x[ 6], x[11]); + word3_muladd_2(&w2, &w1, &w0, x[ 7], x[10]); + word3_muladd_2(&w2, &w1, &w0, x[ 8], x[ 9]); + z[17] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 3], x[15]); + word3_muladd_2(&w2, &w1, &w0, x[ 4], x[14]); + word3_muladd_2(&w2, &w1, &w0, x[ 5], x[13]); + word3_muladd_2(&w2, &w1, &w0, x[ 6], x[12]); + word3_muladd_2(&w2, &w1, &w0, x[ 7], x[11]); + word3_muladd_2(&w2, &w1, &w0, x[ 8], x[10]); + word3_muladd(&w2, &w1, &w0, x[ 9], x[ 9]); + z[18] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 4], x[15]); + word3_muladd_2(&w2, &w1, &w0, x[ 5], x[14]); + word3_muladd_2(&w2, &w1, &w0, x[ 6], x[13]); + word3_muladd_2(&w2, &w1, &w0, x[ 7], x[12]); + word3_muladd_2(&w2, &w1, &w0, x[ 8], x[11]); + word3_muladd_2(&w2, &w1, &w0, x[ 9], x[10]); + z[19] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 5], x[15]); + word3_muladd_2(&w2, &w1, &w0, x[ 6], x[14]); + word3_muladd_2(&w2, &w1, &w0, x[ 7], x[13]); + word3_muladd_2(&w2, &w1, &w0, x[ 8], x[12]); + word3_muladd_2(&w2, &w1, &w0, x[ 9], x[11]); + word3_muladd(&w2, &w1, &w0, x[10], x[10]); + z[20] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 6], x[15]); + word3_muladd_2(&w2, &w1, &w0, x[ 7], x[14]); + word3_muladd_2(&w2, &w1, &w0, x[ 8], x[13]); + word3_muladd_2(&w2, &w1, &w0, x[ 9], x[12]); + word3_muladd_2(&w2, &w1, &w0, x[10], x[11]); + z[21] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 7], x[15]); + word3_muladd_2(&w2, &w1, &w0, x[ 8], x[14]); + word3_muladd_2(&w2, &w1, &w0, x[ 9], x[13]); + word3_muladd_2(&w2, &w1, &w0, x[10], x[12]); + word3_muladd(&w2, &w1, &w0, x[11], x[11]); + z[22] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 8], x[15]); + word3_muladd_2(&w2, &w1, &w0, x[ 9], x[14]); + word3_muladd_2(&w2, &w1, &w0, x[10], x[13]); + word3_muladd_2(&w2, &w1, &w0, x[11], x[12]); + z[23] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[ 9], x[15]); + word3_muladd_2(&w2, &w1, &w0, x[10], x[14]); + word3_muladd_2(&w2, &w1, &w0, x[11], x[13]); + word3_muladd(&w2, &w1, &w0, x[12], x[12]); + z[24] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[10], x[15]); + word3_muladd_2(&w2, &w1, &w0, x[11], x[14]); + word3_muladd_2(&w2, &w1, &w0, x[12], x[13]); + z[25] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[11], x[15]); + word3_muladd_2(&w2, &w1, &w0, x[12], x[14]); + word3_muladd(&w2, &w1, &w0, x[13], x[13]); + z[26] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[12], x[15]); + word3_muladd_2(&w2, &w1, &w0, x[13], x[14]); + z[27] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[13], x[15]); + word3_muladd(&w2, &w1, &w0, x[14], x[14]); + z[28] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd_2(&w2, &w1, &w0, x[14], x[15]); + z[29] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[15], x[15]); + z[30] = w0; + z[31] = w1; + } + +/************************************************* +* Comba 16x16 Multiplication * +*************************************************/ +void bigint_comba_mul16(word z[32], const word x[16], const word y[16]) + { + word w2 = 0, w1 = 0, w0 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[ 0]); + z[0] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[ 0]); + z[1] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[ 0]); + z[2] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[ 0]); + z[3] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[ 0]); + z[4] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[ 0]); + z[5] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[ 0]); + z[6] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[ 0]); + z[7] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[ 0]); + z[8] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[ 0]); + z[9] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[10]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[10], y[ 0]); + z[10] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[11]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[10]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[10], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[11], y[ 0]); + z[11] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[12]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[11]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[10]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[10], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[11], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[12], y[ 0]); + z[12] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[13]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[12]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[11]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[10]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[10], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[11], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[12], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[13], y[ 0]); + z[13] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[14]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[13]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[12]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[11]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[10]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[10], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[11], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[12], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[13], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[14], y[ 0]); + z[14] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 0], y[15]); + word3_muladd(&w2, &w1, &w0, x[ 1], y[14]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[13]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[12]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[11]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[10]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[10], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[11], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[12], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[13], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[14], y[ 1]); + word3_muladd(&w2, &w1, &w0, x[15], y[ 0]); + z[15] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 1], y[15]); + word3_muladd(&w2, &w1, &w0, x[ 2], y[14]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[13]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[12]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[11]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[10]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[10], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[11], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[12], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[13], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[14], y[ 2]); + word3_muladd(&w2, &w1, &w0, x[15], y[ 1]); + z[16] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 2], y[15]); + word3_muladd(&w2, &w1, &w0, x[ 3], y[14]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[13]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[12]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[11]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[10]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[10], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[11], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[12], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[13], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[14], y[ 3]); + word3_muladd(&w2, &w1, &w0, x[15], y[ 2]); + z[17] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 3], y[15]); + word3_muladd(&w2, &w1, &w0, x[ 4], y[14]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[13]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[12]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[11]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[10]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[10], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[11], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[12], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[13], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[14], y[ 4]); + word3_muladd(&w2, &w1, &w0, x[15], y[ 3]); + z[18] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 4], y[15]); + word3_muladd(&w2, &w1, &w0, x[ 5], y[14]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[13]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[12]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[11]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[10]); + word3_muladd(&w2, &w1, &w0, x[10], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[11], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[12], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[13], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[14], y[ 5]); + word3_muladd(&w2, &w1, &w0, x[15], y[ 4]); + z[19] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 5], y[15]); + word3_muladd(&w2, &w1, &w0, x[ 6], y[14]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[13]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[12]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[11]); + word3_muladd(&w2, &w1, &w0, x[10], y[10]); + word3_muladd(&w2, &w1, &w0, x[11], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[12], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[13], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[14], y[ 6]); + word3_muladd(&w2, &w1, &w0, x[15], y[ 5]); + z[20] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 6], y[15]); + word3_muladd(&w2, &w1, &w0, x[ 7], y[14]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[13]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[12]); + word3_muladd(&w2, &w1, &w0, x[10], y[11]); + word3_muladd(&w2, &w1, &w0, x[11], y[10]); + word3_muladd(&w2, &w1, &w0, x[12], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[13], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[14], y[ 7]); + word3_muladd(&w2, &w1, &w0, x[15], y[ 6]); + z[21] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 7], y[15]); + word3_muladd(&w2, &w1, &w0, x[ 8], y[14]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[13]); + word3_muladd(&w2, &w1, &w0, x[10], y[12]); + word3_muladd(&w2, &w1, &w0, x[11], y[11]); + word3_muladd(&w2, &w1, &w0, x[12], y[10]); + word3_muladd(&w2, &w1, &w0, x[13], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[14], y[ 8]); + word3_muladd(&w2, &w1, &w0, x[15], y[ 7]); + z[22] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 8], y[15]); + word3_muladd(&w2, &w1, &w0, x[ 9], y[14]); + word3_muladd(&w2, &w1, &w0, x[10], y[13]); + word3_muladd(&w2, &w1, &w0, x[11], y[12]); + word3_muladd(&w2, &w1, &w0, x[12], y[11]); + word3_muladd(&w2, &w1, &w0, x[13], y[10]); + word3_muladd(&w2, &w1, &w0, x[14], y[ 9]); + word3_muladd(&w2, &w1, &w0, x[15], y[ 8]); + z[23] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[ 9], y[15]); + word3_muladd(&w2, &w1, &w0, x[10], y[14]); + word3_muladd(&w2, &w1, &w0, x[11], y[13]); + word3_muladd(&w2, &w1, &w0, x[12], y[12]); + word3_muladd(&w2, &w1, &w0, x[13], y[11]); + word3_muladd(&w2, &w1, &w0, x[14], y[10]); + word3_muladd(&w2, &w1, &w0, x[15], y[ 9]); + z[24] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[10], y[15]); + word3_muladd(&w2, &w1, &w0, x[11], y[14]); + word3_muladd(&w2, &w1, &w0, x[12], y[13]); + word3_muladd(&w2, &w1, &w0, x[13], y[12]); + word3_muladd(&w2, &w1, &w0, x[14], y[11]); + word3_muladd(&w2, &w1, &w0, x[15], y[10]); + z[25] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[11], y[15]); + word3_muladd(&w2, &w1, &w0, x[12], y[14]); + word3_muladd(&w2, &w1, &w0, x[13], y[13]); + word3_muladd(&w2, &w1, &w0, x[14], y[12]); + word3_muladd(&w2, &w1, &w0, x[15], y[11]); + z[26] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[12], y[15]); + word3_muladd(&w2, &w1, &w0, x[13], y[14]); + word3_muladd(&w2, &w1, &w0, x[14], y[13]); + word3_muladd(&w2, &w1, &w0, x[15], y[12]); + z[27] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[13], y[15]); + word3_muladd(&w2, &w1, &w0, x[14], y[14]); + word3_muladd(&w2, &w1, &w0, x[15], y[13]); + z[28] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[14], y[15]); + word3_muladd(&w2, &w1, &w0, x[15], y[14]); + z[29] = w0; w0 = w1; w1 = w2; w2 = 0; + + word3_muladd(&w2, &w1, &w0, x[15], y[15]); + z[30] = w0; + z[31] = w1; + } + +} + +} diff --git a/src/math/bigint/mp_core.h b/src/math/bigint/mp_core.h new file mode 100644 index 000000000..92949cd83 --- /dev/null +++ b/src/math/bigint/mp_core.h @@ -0,0 +1,96 @@ +/************************************************* +* MPI Algorithms Header File * +* (C) 1999-2007 Jack Lloyd * +*************************************************/ + +#ifndef BOTAN_MP_CORE_H__ +#define BOTAN_MP_CORE_H__ + +#include <botan/mp_types.h> + +namespace Botan { + +/************************************************* +* The size of the word type, in bits * +*************************************************/ +const u32bit MP_WORD_BITS = BOTAN_MP_WORD_BITS; + +extern "C" { + +/************************************************* +* Addition/Subtraction Operations * +*************************************************/ +void bigint_add2(word[], u32bit, const word[], u32bit); +void bigint_add3(word[], const word[], u32bit, const word[], u32bit); + +word bigint_add2_nc(word[], u32bit, const word[], u32bit); +word bigint_add3_nc(word[], const word[], u32bit, const word[], u32bit); + +void bigint_sub2(word[], u32bit, const word[], u32bit); +void bigint_sub3(word[], const word[], u32bit, const word[], u32bit); + +/************************************************* +* Shift Operations * +*************************************************/ +void bigint_shl1(word[], u32bit, u32bit, u32bit); +void bigint_shl2(word[], const word[], u32bit, u32bit, u32bit); +void bigint_shr1(word[], u32bit, u32bit, u32bit); +void bigint_shr2(word[], const word[], u32bit, u32bit, u32bit); + +/************************************************* +* Simple O(N^2) Multiplication and Squaring * +*************************************************/ +void bigint_simple_mul(word z[], const word x[], u32bit x_size, + const word y[], u32bit y_size); +void bigint_simple_sqr(word z[], const word x[], u32bit x_size); + +/************************************************* +* Linear Multiply * +*************************************************/ +void bigint_linmul2(word[], u32bit, word); +void bigint_linmul3(word[], const word[], u32bit, word); +void bigint_linmul_add(word[], u32bit, const word[], u32bit, word); + +/************************************************* +* Montgomery Reduction * +*************************************************/ +void bigint_monty_redc(word[], u32bit, const word[], u32bit, word); + +/************************************************* +* Misc Utility Operations * +*************************************************/ +u32bit bigint_divcore(word, word, word, word, word, word); +s32bit bigint_cmp(const word[], u32bit, const word[], u32bit); +word bigint_divop(word, word, word); +word bigint_modop(word, word, word); +void bigint_wordmul(word, word, word*, word*); + +/************************************************* +* Comba Multiplication / Squaring * +*************************************************/ +void bigint_comba_mul4(word[8], const word[4], const word[4]); +void bigint_comba_mul6(word[12], const word[6], const word[6]); +void bigint_comba_mul8(word[16], const word[8], const word[8]); +void bigint_comba_mul16(word[32], const word[16], const word[16]); + +void bigint_comba_sqr4(word[8], const word[4]); +void bigint_comba_sqr6(word[12], const word[6]); +void bigint_comba_sqr8(word[16], const word[8]); +void bigint_comba_sqr8(word[32], const word[16]); +void bigint_comba_sqr16(word[64], const word[32]); + +} + +/************************************************* +* High Level Multiplication/Squaring Interfaces * +*************************************************/ +void bigint_mul(word[], u32bit, word[], + const word[], u32bit, u32bit, + const word[], u32bit, u32bit); + +void bigint_sqr(word[], u32bit, word[], + const word[], u32bit, u32bit); + +} + +#endif diff --git a/src/math/bigint/mp_generic/info.txt b/src/math/bigint/mp_generic/info.txt new file mode 100644 index 000000000..8bf75fec3 --- /dev/null +++ b/src/math/bigint/mp_generic/info.txt @@ -0,0 +1,8 @@ +realname "MPI Core (C++)" + +load_on dep + +<add> +mp_asm.h +mp_asmi.h +</add> diff --git a/src/math/bigint/mp_generic/mp_asm.h b/src/math/bigint/mp_generic/mp_asm.h new file mode 100644 index 000000000..e62a57110 --- /dev/null +++ b/src/math/bigint/mp_generic/mp_asm.h @@ -0,0 +1,52 @@ +/************************************************* +* Lowest Level MPI Algorithms Header File * +* (C) 1999-2008 Jack Lloyd * +* 2006 Luca Piccarreta * +*************************************************/ + +#ifndef BOTAN_MP_ASM_H__ +#define BOTAN_MP_ASM_H__ + +#include <botan/mp_types.h> + +#if (BOTAN_MP_WORD_BITS == 8) + typedef Botan::u16bit dword; +#elif (BOTAN_MP_WORD_BITS == 16) + typedef Botan::u32bit dword; +#elif (BOTAN_MP_WORD_BITS == 32) + typedef Botan::u64bit dword; +#elif (BOTAN_MP_WORD_BITS == 64) + #error BOTAN_MP_WORD_BITS can be 64 only with assembly support +#else + #error BOTAN_MP_WORD_BITS must be 8, 16, 32, or 64 +#endif + +namespace Botan { + +extern "C" { + +/************************************************* +* Word Multiply/Add * +*************************************************/ +inline word word_madd2(word a, word b, word* c) + { + dword z = (dword)a * b + *c; + *c = (word)(z >> BOTAN_MP_WORD_BITS); + return (word)z; + } + +/************************************************* +* Word Multiply/Add * +*************************************************/ +inline word word_madd3(word a, word b, word c, word* d) + { + dword z = (dword)a * b + c + *d; + *d = (word)(z >> BOTAN_MP_WORD_BITS); + return (word)z; + } + +} + +} + +#endif diff --git a/src/math/bigint/mp_generic/mp_asmi.h b/src/math/bigint/mp_generic/mp_asmi.h new file mode 100644 index 000000000..d15295154 --- /dev/null +++ b/src/math/bigint/mp_generic/mp_asmi.h @@ -0,0 +1,189 @@ +/************************************************* +* Lowest Level MPI Algorithms Header File * +* (C) 1999-2008 Jack Lloyd * +* 2006 Luca Piccarreta * +*************************************************/ + +#ifndef BOTAN_MP_ASM_INTERNAL_H__ +#define BOTAN_MP_ASM_INTERNAL_H__ + +#include <botan/mp_asm.h> + +namespace Botan { + +extern "C" { + +/************************************************* +* Word Addition * +*************************************************/ +inline word word_add(word x, word y, word* carry) + { + word z = x + y; + word c1 = (z < x); + z += *carry; + *carry = c1 | (z < *carry); + return z; + } + +/************************************************* +* Eight Word Block Addition, Two Argument * +*************************************************/ +inline word word8_add2(word x[8], const word y[8], word carry) + { + x[0] = word_add(x[0], y[0], &carry); + x[1] = word_add(x[1], y[1], &carry); + x[2] = word_add(x[2], y[2], &carry); + x[3] = word_add(x[3], y[3], &carry); + x[4] = word_add(x[4], y[4], &carry); + x[5] = word_add(x[5], y[5], &carry); + x[6] = word_add(x[6], y[6], &carry); + x[7] = word_add(x[7], y[7], &carry); + return carry; + } + +/************************************************* +* Eight Word Block Addition, Three Argument * +*************************************************/ +inline word word8_add3(word z[8], const word x[8], + const word y[8], word carry) + { + z[0] = word_add(x[0], y[0], &carry); + z[1] = word_add(x[1], y[1], &carry); + z[2] = word_add(x[2], y[2], &carry); + z[3] = word_add(x[3], y[3], &carry); + z[4] = word_add(x[4], y[4], &carry); + z[5] = word_add(x[5], y[5], &carry); + z[6] = word_add(x[6], y[6], &carry); + z[7] = word_add(x[7], y[7], &carry); + return carry; + } + +/************************************************* +* Word Subtraction * +*************************************************/ +inline word word_sub(word x, word y, word* carry) + { + word t0 = x - y; + word c1 = (t0 > x); + word z = t0 - *carry; + *carry = c1 | (z > t0); + return z; + } + +/************************************************* +* Eight Word Block Subtraction, Two Argument * +*************************************************/ +inline word word8_sub2(word x[4], const word y[4], word carry) + { + x[0] = word_sub(x[0], y[0], &carry); + x[1] = word_sub(x[1], y[1], &carry); + x[2] = word_sub(x[2], y[2], &carry); + x[3] = word_sub(x[3], y[3], &carry); + x[4] = word_sub(x[4], y[4], &carry); + x[5] = word_sub(x[5], y[5], &carry); + x[6] = word_sub(x[6], y[6], &carry); + x[7] = word_sub(x[7], y[7], &carry); + return carry; + } + +/************************************************* +* Eight Word Block Subtraction, Three Argument * +*************************************************/ +inline word word8_sub3(word z[8], const word x[8], + const word y[8], word carry) + { + z[0] = word_sub(x[0], y[0], &carry); + z[1] = word_sub(x[1], y[1], &carry); + z[2] = word_sub(x[2], y[2], &carry); + z[3] = word_sub(x[3], y[3], &carry); + z[4] = word_sub(x[4], y[4], &carry); + z[5] = word_sub(x[5], y[5], &carry); + z[6] = word_sub(x[6], y[6], &carry); + z[7] = word_sub(x[7], y[7], &carry); + return carry; + } + +/************************************************* +* Eight Word Block Linear Multiplication * +*************************************************/ +inline word word8_linmul2(word x[4], word y, word carry) + { + x[0] = word_madd2(x[0], y, &carry); + x[1] = word_madd2(x[1], y, &carry); + x[2] = word_madd2(x[2], y, &carry); + x[3] = word_madd2(x[3], y, &carry); + x[4] = word_madd2(x[4], y, &carry); + x[5] = word_madd2(x[5], y, &carry); + x[6] = word_madd2(x[6], y, &carry); + x[7] = word_madd2(x[7], y, &carry); + return carry; + } + +/************************************************* +* Eight Word Block Linear Multiplication * +*************************************************/ +inline word word8_linmul3(word z[8], const word x[8], word y, word carry) + { + z[0] = word_madd2(x[0], y, &carry); + z[1] = word_madd2(x[1], y, &carry); + z[2] = word_madd2(x[2], y, &carry); + z[3] = word_madd2(x[3], y, &carry); + z[4] = word_madd2(x[4], y, &carry); + z[5] = word_madd2(x[5], y, &carry); + z[6] = word_madd2(x[6], y, &carry); + z[7] = word_madd2(x[7], y, &carry); + return carry; + } + +/************************************************* +* Eight Word Block Multiply/Add * +*************************************************/ +inline word word8_madd3(word z[8], const word x[8], word y, word carry) + { + z[0] = word_madd3(x[0], y, z[0], &carry); + z[1] = word_madd3(x[1], y, z[1], &carry); + z[2] = word_madd3(x[2], y, z[2], &carry); + z[3] = word_madd3(x[3], y, z[3], &carry); + z[4] = word_madd3(x[4], y, z[4], &carry); + z[5] = word_madd3(x[5], y, z[5], &carry); + z[6] = word_madd3(x[6], y, z[6], &carry); + z[7] = word_madd3(x[7], y, z[7], &carry); + return carry; + } + +/************************************************* +* Multiply-Add Accumulator * +*************************************************/ +inline void word3_muladd(word* w2, word* w1, word* w0, word a, word b) + { + word carry = *w0; + *w0 = word_madd2(a, b, &carry); + *w1 += carry; + *w2 += (*w1 < carry) ? 1 : 0; + } + +/************************************************* +* Multiply-Add Accumulator * +*************************************************/ +inline void word3_muladd_2(word* w2, word* w1, word* w0, word a, word b) + { + word carry = 0; + a = word_madd2(a, b, &carry); + b = carry; + + word top = (b >> (BOTAN_MP_WORD_BITS-1)); + b <<= 1; + b |= (a >> (BOTAN_MP_WORD_BITS-1)); + a <<= 1; + + carry = 0; + *w0 = word_add(*w0, a, &carry); + *w1 = word_add(*w1, b, &carry); + *w2 = word_add(*w2, top, &carry); + } + +} + +} + +#endif diff --git a/src/math/bigint/mp_ia32/info.txt b/src/math/bigint/mp_ia32/info.txt new file mode 100644 index 000000000..51f98fda8 --- /dev/null +++ b/src/math/bigint/mp_ia32/info.txt @@ -0,0 +1,19 @@ +realname "MPI Core (IA-32)" + +mp_bits 32 + +load_on asm_ok + +<add> +mp_asm.h +mp_asmi.h +</add> + +<arch> +ia32 +</arch> + +<cc> +gcc +icc +</cc> diff --git a/src/math/bigint/mp_ia32/mp_asm.h b/src/math/bigint/mp_ia32/mp_asm.h new file mode 100644 index 000000000..b45140321 --- /dev/null +++ b/src/math/bigint/mp_ia32/mp_asm.h @@ -0,0 +1,65 @@ +/************************************************* +* Lowest Level MPI Algorithms Header File * +* (C) 1999-2008 Jack Lloyd * +* 2006 Luca Piccarreta * +*************************************************/ + +#ifndef BOTAN_MP_ASM_H__ +#define BOTAN_MP_ASM_H__ + +#include <botan/mp_types.h> + +#if (BOTAN_MP_WORD_BITS != 32) + #error The mp_ia32 module requires that BOTAN_MP_WORD_BITS == 32 +#endif + +namespace Botan { + +extern "C" { + +/************************************************* +* Helper Macros for x86 Assembly * +*************************************************/ +#define ASM(x) x "\n\t" + +/************************************************* +* Word Multiply * +*************************************************/ +inline word word_madd2(word a, word b, word* c) + { + asm( + ASM("mull %[b]") + ASM("addl %[c],%[a]") + ASM("adcl $0,%[carry]") + + : [a]"=a"(a), [b]"=rm"(b), [carry]"=&d"(*c) + : "0"(a), "1"(b), [c]"g"(*c) : "cc"); + + return a; + } + +/************************************************* +* Word Multiply/Add * +*************************************************/ +inline word word_madd3(word a, word b, word c, word* d) + { + asm( + ASM("mull %[b]") + + ASM("addl %[c],%[a]") + ASM("adcl $0,%[carry]") + + ASM("addl %[d],%[a]") + ASM("adcl $0,%[carry]") + + : [a]"=a"(a), [b]"=rm"(b), [carry]"=&d"(*d) + : "0"(a), "1"(b), [c]"g"(c), [d]"g"(*d) : "cc"); + + return a; + } + +} + +} + +#endif diff --git a/src/math/bigint/mp_ia32/mp_asmi.h b/src/math/bigint/mp_ia32/mp_asmi.h new file mode 100644 index 000000000..9de0c11e3 --- /dev/null +++ b/src/math/bigint/mp_ia32/mp_asmi.h @@ -0,0 +1,225 @@ +/************************************************* +* Lowest Level MPI Algorithms Header File * +* (C) 1999-2007 Jack Lloyd * +* 2006 Luca Piccarreta * +*************************************************/ + +#ifndef BOTAN_MP_ASM_INTERNAL_H__ +#define BOTAN_MP_ASM_INTERNAL_H__ + +#include <botan/mp_asm.h> + +namespace Botan { + +extern "C" { + +/************************************************* +* Helper Macros for x86 Assembly * +*************************************************/ +#ifndef ASM + #define ASM(x) x "\n\t" +#endif + +#define ADDSUB2_OP(OPERATION, INDEX) \ + ASM("movl 4*" #INDEX "(%[y]), %[carry]") \ + ASM(OPERATION " %[carry], 4*" #INDEX "(%[x])") \ + +#define ADDSUB3_OP(OPERATION, INDEX) \ + ASM("movl 4*" #INDEX "(%[x]), %[carry]") \ + ASM(OPERATION " 4*" #INDEX "(%[y]), %[carry]") \ + ASM("movl %[carry], 4*" #INDEX "(%[z])") \ + +#define LINMUL_OP(WRITE_TO, INDEX) \ + ASM("movl 4*" #INDEX "(%[x]),%%eax") \ + ASM("mull %[y]") \ + ASM("addl %[carry],%%eax") \ + ASM("adcl $0,%%edx") \ + ASM("movl %%edx,%[carry]") \ + ASM("movl %%eax, 4*" #INDEX "(%[" WRITE_TO "])") + +#define MULADD_OP(IGNORED, INDEX) \ + ASM("movl 4*" #INDEX "(%[x]),%%eax") \ + ASM("mull %[y]") \ + ASM("addl %[carry],%%eax") \ + ASM("adcl $0,%%edx") \ + ASM("addl 4*" #INDEX "(%[z]),%%eax") \ + ASM("adcl $0,%%edx") \ + ASM("movl %%edx,%[carry]") \ + ASM("movl %%eax, 4*" #INDEX " (%[z])") + +#define DO_8_TIMES(MACRO, ARG) \ + MACRO(ARG, 0) \ + MACRO(ARG, 1) \ + MACRO(ARG, 2) \ + MACRO(ARG, 3) \ + MACRO(ARG, 4) \ + MACRO(ARG, 5) \ + MACRO(ARG, 6) \ + MACRO(ARG, 7) + +#define ADD_OR_SUBTRACT(CORE_CODE) \ + ASM("rorl %[carry]") \ + CORE_CODE \ + ASM("sbbl %[carry],%[carry]") \ + ASM("negl %[carry]") + +/************************************************* +* Word Addition * +*************************************************/ +inline word word_add(word x, word y, word* carry) + { + asm( + ADD_OR_SUBTRACT(ASM("adcl %[y],%[x]")) + : [x]"=r"(x), [carry]"=r"(*carry) + : "0"(x), [y]"rm"(y), "1"(*carry) + : "cc"); + return x; + } + +/************************************************* +* Eight Word Block Addition, Two Argument * +*************************************************/ +inline word word8_add2(word x[8], const word y[8], word carry) + { + asm( + ADD_OR_SUBTRACT(DO_8_TIMES(ADDSUB2_OP, "adcl")) + : [carry]"=r"(carry) + : [x]"r"(x), [y]"r"(y), "0"(carry) + : "cc", "memory"); + return carry; + } + +/************************************************* +* Eight Word Block Addition, Three Argument * +*************************************************/ +inline word word8_add3(word z[8], const word x[8], const word y[8], word carry) + { + asm( + ADD_OR_SUBTRACT(DO_8_TIMES(ADDSUB3_OP, "adcl")) + : [carry]"=r"(carry) + : [x]"r"(x), [y]"r"(y), [z]"r"(z), "0"(carry) + : "cc", "memory"); + return carry; + } + +/************************************************* +* Word Subtraction * +*************************************************/ +inline word word_sub(word x, word y, word* carry) + { + asm( + ADD_OR_SUBTRACT(ASM("sbbl %[y],%[x]")) + : [x]"=r"(x), [carry]"=r"(*carry) + : "0"(x), [y]"rm"(y), "1"(*carry) + : "cc"); + return x; + } + +/************************************************* +* Eight Word Block Subtraction, Two Argument * +*************************************************/ +inline word word8_sub2(word x[8], const word y[8], word carry) + { + asm( + ADD_OR_SUBTRACT(DO_8_TIMES(ADDSUB2_OP, "sbbl")) + : [carry]"=r"(carry) + : [x]"r"(x), [y]"r"(y), "0"(carry) + : "cc", "memory"); + return carry; + } + +/************************************************* +* Eight Word Block Subtraction, Three Argument * +*************************************************/ +inline word word8_sub3(word z[8], const word x[8], const word y[8], word carry) + { + asm( + ADD_OR_SUBTRACT(DO_8_TIMES(ADDSUB3_OP, "sbbl")) + : [carry]"=r"(carry) + : [x]"r"(x), [y]"r"(y), [z]"r"(z), "0"(carry) + : "cc", "memory"); + return carry; + } + +/************************************************* +* Eight Word Block Linear Multiplication * +*************************************************/ +inline word word8_linmul2(word x[8], word y, word carry) + { + asm( + DO_8_TIMES(LINMUL_OP, "x") + : [carry]"=r"(carry) + : [x]"r"(x), [y]"rm"(y), "0"(carry) + : "cc", "%eax", "%edx"); + return carry; + } + +/************************************************* +* Eight Word Block Linear Multiplication * +*************************************************/ +inline word word8_linmul3(word z[8], const word x[8], word y, word carry) + { + asm( + DO_8_TIMES(LINMUL_OP, "z") + : [carry]"=r"(carry) + : [z]"r"(z), [x]"r"(x), [y]"rm"(y), "0"(carry) + : "cc", "%eax", "%edx"); + return carry; + } + +/************************************************* +* Eight Word Block Multiply/Add * +*************************************************/ +inline word word8_madd3(word z[8], const word x[8], word y, word carry) + { + asm( + DO_8_TIMES(MULADD_OP, "") + : [carry]"=r"(carry) + : [z]"r"(z), [x]"r"(x), [y]"rm"(y), "0"(carry) + : "cc", "%eax", "%edx"); + return carry; + } + +/************************************************* +* Multiply-Add Accumulator * +*************************************************/ +inline void word3_muladd(word* w2, word* w1, word* w0, word x, word y) + { + asm( + ASM("mull %[y]") + + ASM("addl %[x],%[w0]") + ASM("adcl %[y],%[w1]") + ASM("adcl $0,%[w2]") + + : [w0]"=r"(*w0), [w1]"=r"(*w1), [w2]"=r"(*w2) + : [x]"a"(x), [y]"d"(y), "0"(*w0), "1"(*w1), "2"(*w2) + : "cc"); + } + +/************************************************* +* Multiply-Add Accumulator * +*************************************************/ +inline void word3_muladd_2(word* w2, word* w1, word* w0, word x, word y) + { + asm( + ASM("mull %[y]") + + ASM("addl %[x],%[w0]") + ASM("adcl %[y],%[w1]") + ASM("adcl $0,%[w2]") + + ASM("addl %[x],%[w0]") + ASM("adcl %[y],%[w1]") + ASM("adcl $0,%[w2]") + + : [w0]"=r"(*w0), [w1]"=r"(*w1), [w2]"=r"(*w2) + : [x]"a"(x), [y]"d"(y), "0"(*w0), "1"(*w1), "2"(*w2) + : "cc"); + } + +} + +} + +#endif diff --git a/src/math/bigint/mp_ia32_msvc/info.txt b/src/math/bigint/mp_ia32_msvc/info.txt new file mode 100644 index 000000000..9c7ac9b43 --- /dev/null +++ b/src/math/bigint/mp_ia32_msvc/info.txt @@ -0,0 +1,18 @@ +realname "x86 MPI Assembler Core (MSVC)" + +mp_bits 32 + +load_on dep + +<add> +mp_generic:mp_asm.h +mp_asmi.h +</add> + +<arch> +ia32 +</arch> + +<cc> +msvc +</cc> diff --git a/src/math/bigint/mp_ia32_msvc/mp_asmi.h b/src/math/bigint/mp_ia32_msvc/mp_asmi.h new file mode 100644 index 000000000..5eaa46eb4 --- /dev/null +++ b/src/math/bigint/mp_ia32_msvc/mp_asmi.h @@ -0,0 +1,547 @@ +/************************************************* +* Lowest Level MPI Algorithms Header File * +* (C) 1999-2006 Jack Lloyd * +* 2006 Luca Piccarreta * +*************************************************/ + +#ifndef BOTAN_MP_ASM_INTERNAL_H__ +#define BOTAN_MP_ASM_INTERNAL_H__ + +#include "mp_asm.h" + +namespace Botan { + +extern "C" { + +/************************************************* +* Word Addition * +*************************************************/ +inline word word_add(word x, word y, word* carry) + { + word z = x + y; + word c1 = (z < x); + z += *carry; + *carry = c1 | (z < *carry); + return z; + } + +/************************************************* +* Four Word Block Addition, Two Argument * +*************************************************/ +inline word word4_addcarry(word x[4], word carry) + { + __asm { + mov edx,[x] + xor eax,eax + sub eax,[carry] //force CF=1 iff *carry==1 + adc [edx],0 + mov eax,[esi+4] + adc [edx+4],0 + mov eax,[esi+8] + adc [edx+8],0 + mov eax,[esi+12] + adc [edx+12],0 + sbb eax,eax + neg eax + } + } + +/************************************************* +* Four Word Block Addition, Two Argument * +*************************************************/ +inline word word8_add2(word x[8], const word y[8], word carry) + { + __asm { + mov edx,[x] + mov esi,[y] + xor eax,eax + sub eax,[carry] //force CF=1 iff *carry==1 + mov eax,[esi] + adc [edx],eax + mov eax,[esi+4] + adc [edx+4],eax + mov eax,[esi+8] + adc [edx+8],eax + mov eax,[esi+12] + adc [edx+12],eax + mov eax,[esi+16] + adc [edx+16],eax + mov eax,[esi+20] + adc [edx+20],eax + mov eax,[esi+24] + adc [edx+24],eax + mov eax,[esi+28] + adc [edx+28],eax + sbb eax,eax + neg eax + } + } + +/************************************************* +* Four Word Block Addition, Three Argument * +*************************************************/ +inline word word8_add3(word z[8], const word x[8], const word y[8], word carry) + { + __asm { + mov edi,[x] + mov esi,[y] + mov ebx,[z] + xor eax,eax + sub eax,[carry] //force CF=1 iff *carry==1 + mov eax,[edi] + adc eax,[esi] + mov [ebx],eax + + mov eax,[edi+4] + adc eax,[esi+4] + mov [ebx+4],eax + + mov eax,[edi+8] + adc eax,[esi+8] + mov [ebx+8],eax + + mov eax,[edi+12] + adc eax,[esi+12] + mov [ebx+12],eax + + mov eax,[edi+16] + adc eax,[esi+16] + mov [ebx+16],eax + + mov eax,[edi+20] + adc eax,[esi+20] + mov [ebx+20],eax + + mov eax,[edi+24] + adc eax,[esi+24] + mov [ebx+24],eax + + mov eax,[edi+28] + adc eax,[esi+28] + mov [ebx+28],eax + + sbb eax,eax + neg eax + } + } + +/************************************************* +* Word Subtraction * +*************************************************/ +inline word word_sub(word x, word y, word* carry) + { + word t0 = x - y; + word c1 = (t0 > x); + word z = t0 - *carry; + *carry = c1 | (z > t0); + return z; + } + +/************************************************* +* Four Word Block Subtraction, Two Argument * +*************************************************/ +inline word word8_sub2(word x[8], const word y[8], word carry) + { + _asm { + mov edi,[x] + mov esi,[y] + xor eax,eax + sub eax,[carry] //force CF=1 iff *carry==1 + mov eax,[edi] + sbb eax,[esi] + mov [edi],eax + mov eax,[edi+4] + sbb eax,[esi+4] + mov [edi+4],eax + mov eax,[edi+8] + sbb eax,[esi+8] + mov [edi+8],eax + mov eax,[edi+12] + sbb eax,[esi+12] + mov [edi+12],eax + mov eax,[edi+16] + sbb eax,[esi+16] + mov [edi+16],eax + mov eax,[edi+20] + sbb eax,[esi+20] + mov [edi+20],eax + mov eax,[edi+24] + sbb eax,[esi+24] + mov [edi+24],eax + mov eax,[edi+28] + sbb eax,[esi+28] + mov [edi+28],eax + sbb eax,eax + neg eax + } + } + +/************************************************* +* Four Word Block Subtraction, Three Argument * +*************************************************/ +__forceinline word word8_sub3(word z[8], const word x[8], + const word y[8], word carry) + { + __asm { + mov edi,[x] + mov esi,[y] + xor eax,eax + sub eax,[carry] //force CF=1 iff *carry==1 + mov ebx,[z] + mov eax,[edi] + sbb eax,[esi] + mov [ebx],eax + mov eax,[edi+4] + sbb eax,[esi+4] + mov [ebx+4],eax + mov eax,[edi+8] + sbb eax,[esi+8] + mov [ebx+8],eax + mov eax,[edi+12] + sbb eax,[esi+12] + mov [ebx+12],eax + mov eax,[edi+16] + sbb eax,[esi+16] + mov [ebx+16],eax + mov eax,[edi+20] + sbb eax,[esi+20] + mov [ebx+20],eax + mov eax,[edi+24] + sbb eax,[esi+24] + mov [ebx+24],eax + mov eax,[edi+28] + sbb eax,[esi+28] + mov [ebx+28],eax + sbb eax,eax + neg eax + } + } + +/************************************************* +* Four Word Block Linear Multiplication * +*************************************************/ +inline word word8_linmul2(word x[8], word y, word carry) + { + __asm + { + mov esi,[x] + mov eax,[esi] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,[carry] //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [esi],eax //load a + + mov eax,[esi+4] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [esi+4],eax //load a + + mov eax,[esi+8] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [esi+8],eax //load a + + mov eax,[esi+12] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [esi+12],eax //load a + + mov eax,[esi+16] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [esi+16],eax //load a + + mov eax,[esi+20] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [esi+20],eax //load a + + mov eax,[esi+24] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [esi+24],eax //load a + + mov eax,[esi+28] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov [esi+28],eax //load a + + mov eax,edx //store carry + } + } + +/************************************************* +* Eight Word Block Linear Multiplication * +*************************************************/ +__forceinline word word8_muladd(word z[8], const word x[8], + word y, word carry) + { + __asm + { + mov esi,[x] + mov ebx,[y] + mov edi,[z] + mov eax,[esi] //load a + mul ebx //edx(hi):eax(lo)=a*b + add eax,[carry] //sum lo carry + adc edx,0 //sum hi carry + add eax,[edi] //sum lo z + adc edx,0 //sum hi z + mov ecx,edx //carry for next block = hi z + mov [edi],eax //save lo z + + mov eax,[esi+4] + mul ebx + add eax,ecx + adc edx,0 + add eax,[edi+4] + adc edx,0 + mov ecx,edx + mov [edi+4],eax + + mov eax,[esi+8] + mul ebx + add eax,ecx + adc edx,0 + add eax,[edi+8] + adc edx,0 + mov ecx,edx + mov [edi+8],eax + + mov eax,[esi+12] + mul ebx + add eax,ecx + adc edx,0 + add eax,[edi+12] + adc edx,0 + mov ecx,edx + mov [edi+12],eax + + mov eax,[esi+16] + mul ebx + add eax,ecx + adc edx,0 + add eax,[edi+16] + adc edx,0 + mov ecx,edx + mov [edi+16],eax + + mov eax,[esi+20] + mul ebx + add eax,ecx + adc edx,0 + add eax,[edi+20] + adc edx,0 + mov ecx,edx + mov [edi+20],eax + + mov eax,[esi+24] + mul ebx + add eax,ecx + adc edx,0 + add eax,[edi+24] + adc edx,0 + mov ecx,edx + mov [edi+24],eax + + mov eax,[esi+28] + mul ebx + add eax,ecx + adc edx,0 + add eax,[edi+28] + adc edx,0 + mov [edi+28],eax + mov eax,edx + } + } + +__forceinline word word8_linmul3(word z[4], const word x[4], word y, word carry) + { + __asm + { +#if 0 + //it's slower!!! + mov edx,[z] + mov eax,[x] + movd mm7,[y] + + movd mm0,[eax] + movd mm1,[eax+4] + movd mm2,[eax+8] + pmuludq mm0,mm7 + pmuludq mm1,mm7 + pmuludq mm2,mm7 + + movd mm6,[carry] + paddq mm0,mm6 + movd [edx],mm0 + + psrlq mm0,32 + paddq mm1,mm0 + movd [edx+4],mm1 + + movd mm3,[eax+12] + psrlq mm1,32 + paddq mm2,mm1 + movd [edx+8],mm2 + + pmuludq mm3,mm7 + movd mm4,[eax+16] + psrlq mm2,32 + paddq mm3,mm2 + movd [edx+12],mm3 + + pmuludq mm4,mm7 + movd mm5,[eax+20] + psrlq mm3,32 + paddq mm4,mm3 + movd [edx+16],mm4 + + pmuludq mm5,mm7 + movd mm0,[eax+24] + psrlq mm4,32 + paddq mm5,mm4 + movd [edx+20],mm5 + + pmuludq mm0,mm7 + movd mm1,[eax+28] + psrlq mm5,32 + paddq mm0,mm5 + movd [edx+24],mm0 + + pmuludq mm1,mm7 + psrlq mm0,32 + paddq mm1,mm0 + movd [edx+28],mm1 + + psrlq mm1,32 + movd eax,mm1 + emms +#else + mov edi,[z] + mov esi,[x] + mov eax,[esi] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,[carry] //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [edi],eax //load a + + mov eax,[esi+4] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [edi+4],eax //load a + + mov eax,[esi+8] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [edi+8],eax //load a + + mov eax,[esi+12] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [edi+12],eax //load a + + mov eax,[esi+16] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [edi+16],eax //load a + + mov eax,[esi+20] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [edi+20],eax //load a + + mov eax,[esi+24] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov ecx,edx //store carry + mov [edi+24],eax //load a + + mov eax,[esi+28] //load a + mul [y] //edx(hi):eax(lo)=a*b + add eax,ecx //sum lo carry + adc edx,0 //sum hi carry + mov [edi+28],eax //load a + mov eax,edx //store carry +#endif + } + } + +/************************************************* +* Eight Word Block Multiply/Add * +*************************************************/ +inline word word8_madd3(word z[8], const word x[8], word y, word carry) + { + z[0] = word_madd3(x[0], y, z[0], &carry); + z[1] = word_madd3(x[1], y, z[1], &carry); + z[2] = word_madd3(x[2], y, z[2], &carry); + z[3] = word_madd3(x[3], y, z[3], &carry); + z[4] = word_madd3(x[4], y, z[4], &carry); + z[5] = word_madd3(x[5], y, z[5], &carry); + z[6] = word_madd3(x[6], y, z[6], &carry); + z[7] = word_madd3(x[7], y, z[7], &carry); + return carry; + } + +/************************************************* +* Multiply-Add Accumulator * +*************************************************/ +inline void word3_muladd(word* w2, word* w1, word* w0, word a, word b) + { + dword z = (dword)a * b + (*w0); + *w0 = (word)z; //lo + + word t1 = (word)(z >> BOTAN_MP_WORD_BITS); //hi + *w1 += t1; //w1+=lo + *w2 += (*w1 < t1) ? 1 : 0; //w2+=carry + } + +/************************************************* +* Multiply-Add Accumulator * +*************************************************/ +inline void word3_muladd_2(word* w2, word* w1, word* w0, word a, word b) + { + dword z = (dword)a * b; + word t0 = (word)z; + word t1 = (word)(z >> BOTAN_MP_WORD_BITS); + + *w0 += t0; + *w1 += t1 + ((*w0 < t0) ? 1 : 0); + *w2 += (*w1 < t1) ? 1 : 0; + + *w0 += t0; + *w1 += t1 + ((*w0 < t0) ? 1 : 0); + *w2 += (*w1 < t1) ? 1 : 0; + } + +} + +} + +#endif diff --git a/src/math/bigint/mp_karat.cpp b/src/math/bigint/mp_karat.cpp new file mode 100644 index 000000000..15b0551fd --- /dev/null +++ b/src/math/bigint/mp_karat.cpp @@ -0,0 +1,334 @@ +/************************************************* +* Karatsuba Multiplication/Squaring Source File * +* (C) 1999-2008 Jack Lloyd * +*************************************************/ + +#include <botan/mp_core.h> +#include <botan/mem_ops.h> +#include <botan/mp_asmi.h> + +namespace Botan { + +namespace { + +/************************************************* +* Karatsuba Multiplication Operation * +*************************************************/ +void karatsuba_mul(word z[], const word x[], const word y[], u32bit N, + word workspace[]) + { + if(N == 6) + bigint_comba_mul6(z, x, y); + else if(N == 8) + bigint_comba_mul8(z, x, y); + else if(N == 16) + bigint_comba_mul16(z, x, y); + else if(N < BOTAN_KARAT_MUL_THRESHOLD || N % 2) + bigint_simple_mul(z, x, N, y, N); + else + { + const u32bit N2 = N / 2; + + const word* x0 = x; + const word* x1 = x + N2; + const word* y0 = y; + const word* y1 = y + N2; + word* z0 = z; + word* z1 = z + N; + + const s32bit cmp0 = bigint_cmp(x0, N2, x1, N2); + const s32bit cmp1 = bigint_cmp(y1, N2, y0, N2); + + clear_mem(workspace, 2*N); + + if(cmp0 && cmp1) + { + if(cmp0 > 0) + bigint_sub3(z0, x0, N2, x1, N2); + else + bigint_sub3(z0, x1, N2, x0, N2); + + if(cmp1 > 0) + bigint_sub3(z1, y1, N2, y0, N2); + else + bigint_sub3(z1, y0, N2, y1, N2); + + karatsuba_mul(workspace, z0, z1, N2, workspace+N); + } + + karatsuba_mul(z0, x0, y0, N2, workspace+N); + karatsuba_mul(z1, x1, y1, N2, workspace+N); + + const u32bit blocks_of_8 = N - (N % 8); + + word carry = 0; + + for(u32bit j = 0; j != blocks_of_8; j += 8) + carry = word8_add3(workspace + N + j, z0 + j, z1 + j, carry); + + for(u32bit j = blocks_of_8; j != N; ++j) + workspace[N + j] = word_add(z0[j], z1[j], &carry); + + word carry2 = 0; + + for(u32bit j = 0; j != blocks_of_8; j += 8) + carry2 = word8_add2(z + N2 + j, workspace + N + j, carry2); + + for(u32bit j = blocks_of_8; j != N; ++j) + z[N2 + j] = word_add(z[N2 + j], workspace[N + j], &carry2); + + z[N + N2] = word_add(z[N + N2], carry2, &carry); + + if(carry) + for(u32bit j = 1; j != N2; ++j) + if(++z[N + N2 + j]) + break; + + if((cmp0 == cmp1) || (cmp0 == 0) || (cmp1 == 0)) + bigint_add2(z + N2, 2*N-N2, workspace, N); + else + bigint_sub2(z + N2, 2*N-N2, workspace, N); + } + } + +/************************************************* +* Karatsuba Squaring Operation * +*************************************************/ +void karatsuba_sqr(word z[], const word x[], u32bit N, word workspace[]) + { + if(N == 6) + bigint_comba_sqr6(z, x); + else if(N == 8) + bigint_comba_sqr8(z, x); + else if(N == 16) + bigint_comba_sqr16(z, x); + else if(N < BOTAN_KARAT_SQR_THRESHOLD || N % 2) + bigint_simple_sqr(z, x, N); + else + { + const u32bit N2 = N / 2; + + const word* x0 = x; + const word* x1 = x + N2; + word* z0 = z; + word* z1 = z + N; + + const s32bit cmp = bigint_cmp(x0, N2, x1, N2); + + clear_mem(workspace, 2*N); + + if(cmp) + { + if(cmp > 0) + bigint_sub3(z0, x0, N2, x1, N2); + else + bigint_sub3(z0, x1, N2, x0, N2); + + karatsuba_sqr(workspace, z0, N2, workspace+N); + } + + karatsuba_sqr(z0, x0, N2, workspace+N); + karatsuba_sqr(z1, x1, N2, workspace+N); + + const u32bit blocks_of_8 = N - (N % 8); + + word carry = 0; + + for(u32bit j = 0; j != blocks_of_8; j += 8) + carry = word8_add3(workspace + N + j, z0 + j, z1 + j, carry); + + for(u32bit j = blocks_of_8; j != N; ++j) + workspace[N + j] = word_add(z0[j], z1[j], &carry); + + word carry2 = 0; + + for(u32bit j = 0; j != blocks_of_8; j += 8) + carry2 = word8_add2(z + N2 + j, workspace + N + j, carry2); + + for(u32bit j = blocks_of_8; j != N; ++j) + z[N2 + j] = word_add(z[N2 + j], workspace[N + j], &carry2); + + z[N + N2] = word_add(z[N + N2], carry2, &carry); + + if(carry) + for(u32bit j = 1; j != N2; ++j) + if(++z[N + N2 + j]) + break; + + if(cmp == 0) + bigint_add2(z + N2, 2*N-N2, workspace, N); + else + bigint_sub2(z + N2, 2*N-N2, workspace, N); + } + } + +/************************************************* +* Pick a good size for the Karatsuba multiply * +*************************************************/ +u32bit karatsuba_size(u32bit z_size, + u32bit x_size, u32bit x_sw, + u32bit y_size, u32bit y_sw) + { + if(x_sw > x_size || x_sw > y_size || y_sw > x_size || y_sw > y_size) + return 0; + + if(((x_size == x_sw) && (x_size % 2)) || + ((y_size == y_sw) && (y_size % 2))) + return 0; + + const u32bit start = (x_sw > y_sw) ? x_sw : y_sw; + const u32bit end = (x_size < y_size) ? x_size : y_size; + + if(start == end) + { + if(start % 2) + return 0; + return start; + } + + for(u32bit j = start; j <= end; ++j) + { + if(j % 2) + continue; + + if(2*j > z_size) + return 0; + + if(x_sw <= j && j <= x_size && y_sw <= j && j <= y_size) + { + if(j % 4 == 2 && + (j+2) <= x_size && (j+2) <= y_size && 2*(j+2) <= z_size) + return j+2; + return j; + } + } + + return 0; + } + +/************************************************* +* Pick a good size for the Karatsuba squaring * +*************************************************/ +u32bit karatsuba_size(u32bit z_size, u32bit x_size, u32bit x_sw) + { + if(x_sw == x_size) + { + if(x_sw % 2) + return 0; + return x_sw; + } + + for(u32bit j = x_sw; j <= x_size; ++j) + { + if(j % 2) + continue; + + if(2*j > z_size) + return 0; + + if(j % 4 == 2 && (j+2) <= x_size && 2*(j+2) <= z_size) + return j+2; + return j; + } + + return 0; + } + +} + +/************************************************* +* Multiplication Algorithm Dispatcher * +*************************************************/ +void bigint_mul(word z[], u32bit z_size, word workspace[], + const word x[], u32bit x_size, u32bit x_sw, + const word y[], u32bit y_size, u32bit y_sw) + { + if(x_sw == 1) + { + bigint_linmul3(z, y, y_sw, x[0]); + } + else if(y_sw == 1) + { + bigint_linmul3(z, x, x_sw, y[0]); + } + else if(x_sw <= 4 && x_size >= 4 && + y_sw <= 4 && y_size >= 4 && z_size >= 8) + { + bigint_comba_mul4(z, x, y); + } + else if(x_sw <= 6 && x_size >= 6 && + y_sw <= 6 && y_size >= 6 && z_size >= 12) + { + bigint_comba_mul6(z, x, y); + } + else if(x_sw <= 8 && x_size >= 8 && + y_sw <= 8 && y_size >= 8 && z_size >= 16) + { + bigint_comba_mul8(z, x, y); + } + else if(x_sw <= 16 && x_size >= 16 && + y_sw <= 16 && y_size >= 16 && z_size >= 32) + { + bigint_comba_mul16(z, x, y); + } + else if(x_sw < BOTAN_KARAT_MUL_THRESHOLD || y_sw < BOTAN_KARAT_MUL_THRESHOLD) + bigint_simple_mul(z, x, x_sw, y, y_sw); + else + { + const u32bit N = karatsuba_size(z_size, x_size, x_sw, y_size, y_sw); + + if(N) + { + clear_mem(workspace, 2*N); + karatsuba_mul(z, x, y, N, workspace); + } + else + bigint_simple_mul(z, x, x_sw, y, y_sw); + } + } + +/************************************************* +* Squaring Algorithm Dispatcher * +*************************************************/ +void bigint_sqr(word z[], u32bit z_size, word workspace[], + const word x[], u32bit x_size, u32bit x_sw) + { + if(x_sw == 1) + { + bigint_linmul3(z, x, x_sw, x[0]); + } + else if(x_sw <= 4 && x_size >= 4 && z_size >= 8) + { + bigint_comba_sqr4(z, x); + } + else if(x_sw <= 6 && x_size >= 6 && z_size >= 12) + { + bigint_comba_sqr6(z, x); + } + else if(x_sw <= 8 && x_size >= 8 && z_size >= 16) + { + bigint_comba_sqr8(z, x); + } + else if(x_sw <= 16 && x_size >= 16 && z_size >= 32) + { + bigint_comba_sqr16(z, x); + } + else if(x_size < BOTAN_KARAT_SQR_THRESHOLD) + { + bigint_simple_sqr(z, x, x_sw); + } + else + { + const u32bit N = karatsuba_size(z_size, x_size, x_sw); + + if(N) + { + clear_mem(workspace, 2*N); + karatsuba_sqr(z, x, N, workspace); + } + else + bigint_simple_sqr(z, x, x_sw); + } + } + +} diff --git a/src/math/bigint/mp_misc.cpp b/src/math/bigint/mp_misc.cpp new file mode 100644 index 000000000..db9c8cda0 --- /dev/null +++ b/src/math/bigint/mp_misc.cpp @@ -0,0 +1,92 @@ +/************************************************* +* MP Misc Functions Source File * +* (C) 1999-2008 Jack Lloyd * +*************************************************/ + +#include <botan/mp_core.h> +#include <botan/mp_asm.h> + +namespace Botan { + +extern "C" { + +/************************************************* +* Core Division Operation * +*************************************************/ +u32bit bigint_divcore(word q, word y1, word y2, + word x1, word x2, word x3) + { + word y0 = 0; + y2 = word_madd2(q, y2, &y0); + y1 = word_madd2(q, y1, &y0); + + if(y0 > x1) return 1; + if(y0 < x1) return 0; + if(y1 > x2) return 1; + if(y1 < x2) return 0; + if(y2 > x3) return 1; + if(y2 < x3) return 0; + return 0; + } + +/************************************************* +* Compare two MP integers * +*************************************************/ +s32bit bigint_cmp(const word x[], u32bit x_size, + const word y[], u32bit y_size) + { + if(x_size < y_size) { return (-bigint_cmp(y, y_size, x, x_size)); } + + while(x_size > y_size) + { + if(x[x_size-1]) + return 1; + x_size--; + } + for(u32bit j = x_size; j > 0; --j) + { + if(x[j-1] > y[j-1]) return 1; + if(x[j-1] < y[j-1]) return -1; + } + return 0; + } + +/************************************************* +* Do a 2-word/1-word Division * +*************************************************/ +word bigint_divop(word n1, word n0, word d) + { + word high = n1 % d, quotient = 0; + + for(u32bit j = 0; j != MP_WORD_BITS; ++j) + { + word high_top_bit = (high & MP_WORD_TOP_BIT); + + high <<= 1; + high |= (n0 >> (MP_WORD_BITS-1-j)) & 1; + quotient <<= 1; + + if(high_top_bit || high >= d) + { + high -= d; + quotient |= 1; + } + } + + return quotient; + } + +/************************************************* +* Do a 2-word/1-word Modulo * +*************************************************/ +word bigint_modop(word n1, word n0, word d) + { + word z = bigint_divop(n1, n0, d); + word dummy = 0; + z = word_madd2(z, d, &dummy); + return (n0-z); + } + +} + +} diff --git a/src/math/bigint/mp_shift.cpp b/src/math/bigint/mp_shift.cpp new file mode 100644 index 000000000..033774e46 --- /dev/null +++ b/src/math/bigint/mp_shift.cpp @@ -0,0 +1,136 @@ +/************************************************* +* MP Shift Algorithms Source File * +* (C) 1999-2007 Jack Lloyd * +*************************************************/ + +#include <botan/mp_core.h> +#include <botan/mem_ops.h> + +namespace Botan { + +extern "C" { + +/************************************************* +* Single Operand Left Shift * +*************************************************/ +void bigint_shl1(word x[], u32bit x_size, u32bit word_shift, u32bit bit_shift) + { + if(word_shift) + { + for(u32bit j = 1; j != x_size + 1; ++j) + x[(x_size - j) + word_shift] = x[x_size - j]; + clear_mem(x, word_shift); + } + + if(bit_shift) + { + word carry = 0; + for(u32bit j = word_shift; j != x_size + word_shift + 1; ++j) + { + word temp = x[j]; + x[j] = (temp << bit_shift) | carry; + carry = (temp >> (MP_WORD_BITS - bit_shift)); + } + } + } + +/************************************************* +* Single Operand Right Shift * +*************************************************/ +void bigint_shr1(word x[], u32bit x_size, u32bit word_shift, u32bit bit_shift) + { + if(x_size < word_shift) + { + clear_mem(x, x_size); + return; + } + + if(word_shift) + { + copy_mem(x, x + word_shift, x_size - word_shift); + clear_mem(x + x_size - word_shift, word_shift); + } + + if(bit_shift) + { + word carry = 0; + + u32bit top = x_size - word_shift; + + while(top >= 4) + { + word w = x[top-1]; + x[top-1] = (w >> bit_shift) | carry; + carry = (w << (MP_WORD_BITS - bit_shift)); + + w = x[top-2]; + x[top-2] = (w >> bit_shift) | carry; + carry = (w << (MP_WORD_BITS - bit_shift)); + + w = x[top-3]; + x[top-3] = (w >> bit_shift) | carry; + carry = (w << (MP_WORD_BITS - bit_shift)); + + w = x[top-4]; + x[top-4] = (w >> bit_shift) | carry; + carry = (w << (MP_WORD_BITS - bit_shift)); + + top -= 4; + } + + while(top) + { + word w = x[top-1]; + x[top-1] = (w >> bit_shift) | carry; + carry = (w << (MP_WORD_BITS - bit_shift)); + + top--; + } + } + } + +/************************************************* +* Two Operand Left Shift * +*************************************************/ +void bigint_shl2(word y[], const word x[], u32bit x_size, + u32bit word_shift, u32bit bit_shift) + { + for(u32bit j = 0; j != x_size; ++j) + y[j + word_shift] = x[j]; + if(bit_shift) + { + word carry = 0; + for(u32bit j = word_shift; j != x_size + word_shift + 1; ++j) + { + word w = y[j]; + y[j] = (w << bit_shift) | carry; + carry = (w >> (MP_WORD_BITS - bit_shift)); + } + } + } + +/************************************************* +* Two Operand Right Shift * +*************************************************/ +void bigint_shr2(word y[], const word x[], u32bit x_size, + u32bit word_shift, u32bit bit_shift) + { + if(x_size < word_shift) return; + + for(u32bit j = 0; j != x_size - word_shift; ++j) + y[j] = x[j + word_shift]; + if(bit_shift) + { + word carry = 0; + for(u32bit j = x_size - word_shift; j > 0; --j) + { + word w = y[j-1]; + y[j-1] = (w >> bit_shift) | carry; + carry = (w << (MP_WORD_BITS - bit_shift)); + } + } + } + +} + +} diff --git a/src/math/bigint/mp_types.h b/src/math/bigint/mp_types.h new file mode 100644 index 000000000..81b6d7395 --- /dev/null +++ b/src/math/bigint/mp_types.h @@ -0,0 +1,31 @@ +/************************************************* +* Low Level MPI Types Header File * +* (C) 1999-2007 Jack Lloyd * +*************************************************/ + +#ifndef BOTAN_MPI_TYPES_H__ +#define BOTAN_MPI_TYPES_H__ + +#include <botan/types.h> + +namespace Botan { + +#if (BOTAN_MP_WORD_BITS == 8) + typedef byte word; +#elif (BOTAN_MP_WORD_BITS == 16) + typedef u16bit word; +#elif (BOTAN_MP_WORD_BITS == 32) + typedef u32bit word; +#elif (BOTAN_MP_WORD_BITS == 64) + typedef u64bit word; +#else + #error BOTAN_MP_WORD_BITS must be 8, 16, 32, or 64 +#endif + +const word MP_WORD_MASK = ~static_cast<word>(0); +const word MP_WORD_TOP_BIT = static_cast<word>(1) << (8*sizeof(word) - 1); +const word MP_WORD_MAX = MP_WORD_MASK; + +} + +#endif diff --git a/src/math/bigint/mulop_amd64/info.txt b/src/math/bigint/mulop_amd64/info.txt new file mode 100644 index 000000000..670780d9c --- /dev/null +++ b/src/math/bigint/mulop_amd64/info.txt @@ -0,0 +1,31 @@ +realname "BigInt Multiply-Add (x86-64)" + +mp_bits 64 + +load_on dep + +<add> +mp_mulop_amd64.S +</add> + +<requires> +asm_amd64 +</requires> + +<arch> +amd64 +</arch> + +<cc> +gcc +icc +</cc> + +# ELF systems +<os> +linux +freebsd +netbsd +openbsd +solaris +</os> diff --git a/src/math/bigint/mulop_amd64/mp_mulop.cpp b/src/math/bigint/mulop_amd64/mp_mulop.cpp new file mode 100644 index 000000000..d1aa51489 --- /dev/null +++ b/src/math/bigint/mulop_amd64/mp_mulop.cpp @@ -0,0 +1,94 @@ +/************************************************* +* Simple O(N^2) Multiplication and Squaring * +* (C) 1999-2008 Jack Lloyd * +*************************************************/ + +#include <botan/mp_asm.h> +#include <botan/mp_asmi.h> +#include <botan/mp_core.h> +#include <botan/mem_ops.h> + +namespace Botan { + +extern "C" { + +/************************************************* +* Simple O(N^2) Multiplication * +*************************************************/ +void bigint_simple_mul(word z[], const word x[], u32bit x_size, + const word y[], u32bit y_size) + { + const u32bit blocks = x_size - (x_size % 8); + + clear_mem(z, x_size + y_size); + + for(u32bit i = 0; i != y_size; ++i) + { + word carry = 0; + + for(u32bit j = 0; j != blocks; j += 8) + carry = word8_madd3(z + i + j, x + j, y[i], carry); + + for(u32bit j = blocks; j != x_size; ++j) + z[i+j] = word_madd3(x[j], y[i], z[i+j], &carry); + + z[x_size+i] = carry; + } + } + +inline word word_sqr(word x, + +/************************************************* +* Simple O(N^2) Squaring + +This is exactly the same algorithm as bigint_simple_mul, +however because C/C++ compilers suck at alias analysis it +is good to have the version where the compiler knows +that x == y +*************************************************/ +void bigint_simple_sqr(word z[], const word x[], u32bit x_size) + { + clear_mem(z, 2*x_size); + + for(u32bit i = 0; i != x_size; ++i) + { + const word x_i = x[i]; + + word carry = z[2*i]; + z[2*i] = word_madd2(x_i, x_i, z[2*i], &carry); + + for(u32bit j = i; j != x_size; ++j) + { + // z[i+j] = z[i+j] + 2 * x[j] * x_i + carry; + + /* + load z[i+j] into register + load x[j] into %hi + mulq %[x_i] -> x[i] * x[j] -> %lo:%hi + shlq %lo, $1 + + // put carry bit (cf) from %lo into %temp + xorl %temp + adcq $0, %temp + + // high bit of lo now in cf + shl %hi, $1 + // add in lowest bid from %lo + orl %temp, %hi + + addq %[c], %[lo] + adcq $0, %[hi] + addq %[z_ij], %[lo] + adcq $0, %[hi] + + */ + + } + + z[x_size+i] = carry; + } + } + +} + +} diff --git a/src/math/bigint/mulop_amd64/mp_mulop_amd64.S b/src/math/bigint/mulop_amd64/mp_mulop_amd64.S new file mode 100644 index 000000000..e5bba23fb --- /dev/null +++ b/src/math/bigint/mulop_amd64/mp_mulop_amd64.S @@ -0,0 +1,128 @@ +/************************************************* +* Simple O(N^2) Multiplication and Squaring * +* (C) 1999-2008 Jack Lloyd * +*************************************************/ + +#include <botan/asm_macr.h> + +START_LISTING(mp_mulop.S) + +#if 0 +void bigint_simple_sqr(word z[], const word x[], u32bit x_size) + { + const u32bit blocks = x_size - (x_size % 8); + + clear_mem(z, 2*x_size); + + for(u32bit i = 0; i != x_size; ++i) + { + word carry = 0; + + /* + for(u32bit j = 0; j != blocks; j += 8) + carry = word8_madd3(z + i + j, x + j, x[i], carry); + + for(u32bit j = blocks; j != x_size; ++j) + z[i+j] = word_madd3(x[j], x[i], z[i+j], &carry); + */ + + + for(u32bit j = 0; j != x_size; ++j) + z[i+j] = word_madd3(x[j], x[i], z[i+j], &carry); + + for(u32bit j = 0; j != x_size; ++j) + { + dword z = (dword)a * b + c + *d; + *d = (word)(z >> BOTAN_MP_WORD_BITS); + return (word)z; + } + + + + z[i+j] = word_madd3(x[j], x[i], z[i+j], &carry); + + } + + + + z[x_size+i] = carry; + } + } + +#endif + +START_FUNCTION(bigint_simple_sqr) + +#define Z_ARR ARG_1 +#define X_ARR ARG_2 +//#define X_SIZE ARG_3_32 + +#define CARRY TEMP_1 +#define Z_WORD TEMP_2 +#define LOOP_I TEMP_3 +#define LOOP_J TEMP_4 +#define X_SIZE TEMP_5 +#define MUL_LO %rax +// arg 3, xsize +#define MUL_HI %rdx + +// need arg3 == rdx for multiply + ASSIGN(X_SIZE, ARG3_32) + + ZEROIZE(CARRY) + + ZEROIZE(LOOP_I) + +.LOOP_ZEROIZE_Z: + + cmp LOOP_I, X_SIZE + + + + + JUMP_IF_ZERO(LOOP_CTR, .L_MULADD_DONE) + JUMP_IF_LT(LOOP_CTR, 8, .LOOP_MULADD1) + +#define MULADD_OP(N) \ + ASSIGN(MUL_LO, ARRAY8(X_ARR, N)) ; \ + ASSIGN(Z_WORD, ARRAY8(Z_ARR, N)) ; \ + MUL(Y) ; \ + ADD(Z_WORD, CARRY) ; \ + ASSIGN(CARRY, MUL_HI) ; \ + ADD_LAST_CARRY(CARRY) ; \ + ADD(Z_WORD, MUL_LO) ; \ + ADD_LAST_CARRY(CARRY) ; \ + ASSIGN(ARRAY8(Z_ARR, N), Z_WORD) + +.LOOP_MULADD8: + MULADD_OP(0) + MULADD_OP(1) + MULADD_OP(2) + MULADD_OP(3) + MULADD_OP(4) + MULADD_OP(5) + MULADD_OP(6) + MULADD_OP(7) + + SUB_IMM(LOOP_CTR, 8) + ADD_IMM(Z_ARR, 64) + ADD_IMM(X_ARR, 64) + cmp IMM(8), LOOP_CTR + jge .LOOP_MULADD8 + + JUMP_IF_ZERO(LOOP_CTR, .L_MULADD_DONE) + +ALIGN +.LOOP_MULADD1: + MULADD_OP(0) + + SUB_IMM(LOOP_CTR, 1) + ADD_IMM(Z_ARR, 8) + ADD_IMM(X_ARR, 8) + + cmp IMM(0), LOOP_CTR + jne .LOOP_MULADD1 + +.L_MULADD_DONE: + RETURN_VALUE_IS(CARRY) +END_FUNCTION(bigint_simple_square) diff --git a/src/math/bigint/mulop_generic/info.txt b/src/math/bigint/mulop_generic/info.txt new file mode 100644 index 000000000..28ebe41eb --- /dev/null +++ b/src/math/bigint/mulop_generic/info.txt @@ -0,0 +1,7 @@ +realname "BigInt Multiply-Add" + +load_on dep + +<add> +mp_mulop.cpp +</add> diff --git a/src/math/bigint/mulop_generic/mp_mulop.cpp b/src/math/bigint/mulop_generic/mp_mulop.cpp new file mode 100644 index 000000000..3ab28d306 --- /dev/null +++ b/src/math/bigint/mulop_generic/mp_mulop.cpp @@ -0,0 +1,72 @@ +/************************************************* +* Simple O(N^2) Multiplication and Squaring * +* (C) 1999-2008 Jack Lloyd * +*************************************************/ + +#include <botan/mp_asm.h> +#include <botan/mp_asmi.h> +#include <botan/mp_core.h> +#include <botan/mem_ops.h> + +namespace Botan { + +extern "C" { + +/************************************************* +* Simple O(N^2) Multiplication * +*************************************************/ +void bigint_simple_mul(word z[], const word x[], u32bit x_size, + const word y[], u32bit y_size) + { + const u32bit x_size_8 = x_size - (x_size % 8); + + clear_mem(z, x_size + y_size); + + for(u32bit i = 0; i != y_size; ++i) + { + const word y_i = y[i]; + + word carry = 0; + + for(u32bit j = 0; j != x_size_8; j += 8) + carry = word8_madd3(z + i + j, x + j, y_i, carry); + + for(u32bit j = x_size_8; j != x_size; ++j) + z[i+j] = word_madd3(x[j], y_i, z[i+j], &carry); + + z[x_size+i] = carry; + } + } + +/************************************************* +* Simple O(N^2) Squaring + +This is exactly the same algorithm as bigint_simple_mul, +however because C/C++ compilers suck at alias analysis it +is good to have the version where the compiler knows +that x == y +*************************************************/ +void bigint_simple_sqr(word z[], const word x[], u32bit x_size) + { + const u32bit x_size_8 = x_size - (x_size % 8); + + clear_mem(z, 2*x_size); + + for(u32bit i = 0; i != x_size; ++i) + { + const word x_i = x[i]; + word carry = 0; + + for(u32bit j = 0; j != x_size_8; j += 8) + carry = word8_madd3(z + i + j, x + j, x_i, carry); + + for(u32bit j = x_size_8; j != x_size; ++j) + z[i+j] = word_madd3(x[j], x_i, z[i+j], &carry); + + z[x_size+i] = carry; + } + } + +} + +} diff --git a/src/math/bigint/mulop_ia32/info.txt b/src/math/bigint/mulop_ia32/info.txt new file mode 100644 index 000000000..1c89e95c1 --- /dev/null +++ b/src/math/bigint/mulop_ia32/info.txt @@ -0,0 +1,33 @@ +realname "BigInt Multiply-Add (IA-32)" + +mp_bits 32 + +# Out of date, still implements bigint_mul_add_words + +load_on request + +<add> +mp_mulop.S +</add> + +<requires> +asm_ia32 +</requires> + +<arch> +ia32 +</arch> + +<cc> +gcc +icc +</cc> + +# ELF systems +<os> +linux +freebsd +netbsd +openbsd +solaris +</os> diff --git a/src/math/bigint/mulop_ia32/mp_mulop.S b/src/math/bigint/mulop_ia32/mp_mulop.S new file mode 100644 index 000000000..a5f0d3b27 --- /dev/null +++ b/src/math/bigint/mulop_ia32/mp_mulop.S @@ -0,0 +1,62 @@ +/************************************************* +* Multiply/Add Algorithm Source File * +* (C) 1999-2007 Jack Lloyd * +*************************************************/ + +#include <botan/asm_macr.h> + +START_LISTING(mp_muladd.S) + +START_FUNCTION(bigint_mul_add_words) + SPILL_REGS() +#define PUSHED 4 + +#define LOOP_CTR ESI + ASSIGN(LOOP_CTR, ARG(3)) /* x_size */ + ZEROIZE(EDI) + + ASSIGN(ECX, ARG(1)) /* z[] */ + ASSIGN(EBX, ARG(2)) /* x[] */ + ASSIGN(EBP, ARG(4)) /* y */ + +#define MULADD_OP(N) \ + ASSIGN(EAX, ARRAY4(EBX, N)) ; \ + MUL(EBP) ; \ + ADD_W_CARRY(EAX, EDX, EDI) ; \ + ASSIGN(EDI, EDX) ; \ + ADD_W_CARRY(ARRAY4(ECX, N), EDI, EAX) ; + + JUMP_IF_ZERO(LOOP_CTR, .MUL_ADD_DONE) + JUMP_IF_LT(LOOP_CTR, 8, .MULADD1_LOOP) + +START_LOOP(.MULADD8) + MULADD_OP(0) + MULADD_OP(1) + MULADD_OP(2) + MULADD_OP(3) + MULADD_OP(4) + MULADD_OP(5) + MULADD_OP(6) + MULADD_OP(7) + + SUB_IMM(LOOP_CTR, 8) + ADD_IMM(EBX, 32) + ADD_IMM(ECX, 32) +LOOP_UNTIL_LT(LOOP_CTR, 8, .MULADD8) + + JUMP_IF_ZERO(LOOP_CTR, .MUL_ADD_DONE) + +START_LOOP(.MULADD1) + MULADD_OP(0) + + SUB_IMM(LOOP_CTR, 1) + ADD_IMM(EBX, 4) + ADD_IMM(ECX, 4) +LOOP_UNTIL_EQ(LOOP_CTR, 0, .MULADD1) + +.MUL_ADD_DONE: + + ASSIGN(EAX, EDI) +#undef PUSHED + RESTORE_REGS() +END_FUNCTION(bigint_mul_add_words) diff --git a/src/math/blinding.cpp b/src/math/numbertheory/blinding.cpp index 740904d10..740904d10 100644 --- a/src/math/blinding.cpp +++ b/src/math/numbertheory/blinding.cpp diff --git a/src/math/blinding.h b/src/math/numbertheory/blinding.h index 958686fb1..958686fb1 100644 --- a/src/math/blinding.h +++ b/src/math/numbertheory/blinding.h diff --git a/src/math/def_powm.h b/src/math/numbertheory/def_powm.h index c91ff002c..c91ff002c 100644 --- a/src/math/def_powm.h +++ b/src/math/numbertheory/def_powm.h diff --git a/src/math/dsa_gen.cpp b/src/math/numbertheory/dsa_gen.cpp index baaba66ee..baaba66ee 100644 --- a/src/math/dsa_gen.cpp +++ b/src/math/numbertheory/dsa_gen.cpp diff --git a/src/math/info.txt b/src/math/numbertheory/info.txt index 6ce3583c8..6ce3583c8 100644 --- a/src/math/info.txt +++ b/src/math/numbertheory/info.txt diff --git a/src/math/jacobi.cpp b/src/math/numbertheory/jacobi.cpp index 57c78508a..57c78508a 100644 --- a/src/math/jacobi.cpp +++ b/src/math/numbertheory/jacobi.cpp diff --git a/src/math/make_prm.cpp b/src/math/numbertheory/make_prm.cpp index 30ac9a623..30ac9a623 100644 --- a/src/math/make_prm.cpp +++ b/src/math/numbertheory/make_prm.cpp diff --git a/src/math/mp_numth.cpp b/src/math/numbertheory/mp_numth.cpp index b024d2e2d..b024d2e2d 100644 --- a/src/math/mp_numth.cpp +++ b/src/math/numbertheory/mp_numth.cpp diff --git a/src/math/numthry.cpp b/src/math/numbertheory/numthry.cpp index ffd523e82..ffd523e82 100644 --- a/src/math/numthry.cpp +++ b/src/math/numbertheory/numthry.cpp diff --git a/src/math/numthry.h b/src/math/numbertheory/numthry.h index a271a428a..a271a428a 100644 --- a/src/math/numthry.h +++ b/src/math/numbertheory/numthry.h diff --git a/src/math/pow_mod.cpp b/src/math/numbertheory/pow_mod.cpp index 17ca7b796..17ca7b796 100644 --- a/src/math/pow_mod.cpp +++ b/src/math/numbertheory/pow_mod.cpp diff --git a/src/math/pow_mod.h b/src/math/numbertheory/pow_mod.h index 37e0871da..37e0871da 100644 --- a/src/math/pow_mod.h +++ b/src/math/numbertheory/pow_mod.h diff --git a/src/math/powm_fw.cpp b/src/math/numbertheory/powm_fw.cpp index c29b9f311..c29b9f311 100644 --- a/src/math/powm_fw.cpp +++ b/src/math/numbertheory/powm_fw.cpp diff --git a/src/math/powm_mnt.cpp b/src/math/numbertheory/powm_mnt.cpp index 6091d467a..6091d467a 100644 --- a/src/math/powm_mnt.cpp +++ b/src/math/numbertheory/powm_mnt.cpp diff --git a/src/math/primes.cpp b/src/math/numbertheory/primes.cpp index d005167e5..d005167e5 100644 --- a/src/math/primes.cpp +++ b/src/math/numbertheory/primes.cpp diff --git a/src/math/reducer.cpp b/src/math/numbertheory/reducer.cpp index 47c5c20fc..47c5c20fc 100644 --- a/src/math/reducer.cpp +++ b/src/math/numbertheory/reducer.cpp diff --git a/src/math/reducer.h b/src/math/numbertheory/reducer.h index 48008e73b..48008e73b 100644 --- a/src/math/reducer.h +++ b/src/math/numbertheory/reducer.h diff --git a/src/math/ressol.cpp b/src/math/numbertheory/ressol.cpp index 0cd2b988a..0cd2b988a 100644 --- a/src/math/ressol.cpp +++ b/src/math/numbertheory/ressol.cpp |