dlib C++ Library - bigint_kernel_2.cpp

// Copyright (C) 2003 Davis E. King (davis@dlib.net)
// License: Boost Software License See LICENSE.txt for the full license.
#ifndef DLIB_BIGINT_KERNEL_2_CPp_
#define DLIB_BIGINT_KERNEL_2_CPp_
#include "bigint_kernel_2.h"
#include <iostream>
#include <cmath>
namespace dlib
{
// ----------------------------------------------------------------------------------------
// ----------------------------------------------------------------------------------------
 // member/friend function definitions
// ----------------------------------------------------------------------------------------
// ----------------------------------------------------------------------------------------
 bigint_kernel_2::
 bigint_kernel_2 (
 ) :
 slack(25),
 data(new data_record(slack))
 {}
// ----------------------------------------------------------------------------------------
 bigint_kernel_2::
 bigint_kernel_2 (
 uint32 value
 ) :
 slack(25),
 data(new data_record(slack))
 {
 *(data->number) = static_cast<uint16>(value&0xFFFF);
 *(data->number+1) = static_cast<uint16>((value>>16)&0xFFFF);
 if (*(data->number+1) != 0)
 data->digits_used = 2;
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2::
 bigint_kernel_2 (
 const bigint_kernel_2& item
 ) :
 slack(25),
 data(item.data)
 {
 data->references += 1;
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2::
 ~bigint_kernel_2 (
 )
 {
 if (data->references == 1)
 {
 delete data;
 }
 else
 {
 data->references -= 1;
 }
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 bigint_kernel_2::
 operator+ (
 const bigint_kernel_2& rhs
 ) const
 {
 data_record* temp = new data_record (
 std::max(rhs.data->digits_used,data->digits_used) + slack
 );
 long_add(data,rhs.data,temp);
 return bigint_kernel_2(temp,0);
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator+= (
 const bigint_kernel_2& rhs
 )
 {
 // if there are other references to our data
 if (data->references != 1)
 {
 data_record* temp = new data_record(std::max(data->digits_used,rhs.data->digits_used)+slack);
 data->references -= 1; 
 long_add(data,rhs.data,temp);
 data = temp;
 }
 // if data is not big enough for the result
 else if (data->size <= std::max(data->digits_used,rhs.data->digits_used))
 {
 data_record* temp = new data_record(std::max(data->digits_used,rhs.data->digits_used)+slack);
 long_add(data,rhs.data,temp);
 delete data;
 data = temp;
 }
 // there is enough size and no references
 else
 {
 long_add(data,rhs.data,data);
 }
 return *this;
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 bigint_kernel_2::
 operator- (
 const bigint_kernel_2& rhs
 ) const
 {
 data_record* temp = new data_record (
 data->digits_used + slack
 );
 long_sub(data,rhs.data,temp);
 return bigint_kernel_2(temp,0);
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator-= (
 const bigint_kernel_2& rhs
 )
 {
 // if there are other references to this data 
 if (data->references != 1)
 {
 data_record* temp = new data_record(data->digits_used+slack);
 data->references -= 1;
 long_sub(data,rhs.data,temp);
 data = temp;
 }
 else
 {
 long_sub(data,rhs.data,data);
 }
 return *this;
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 bigint_kernel_2::
 operator* (
 const bigint_kernel_2& rhs
 ) const
 {
 data_record* temp = new data_record (
 data->digits_used + rhs.data->digits_used + slack
 );
 long_mul(data,rhs.data,temp);
 return bigint_kernel_2(temp,0);
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator*= (
 const bigint_kernel_2& rhs
 )
 {
 // create a data_record to store the result of the multiplication in
 data_record* temp = new data_record(rhs.data->digits_used+data->digits_used+slack); 
 long_mul(data,rhs.data,temp);
 // if there are other references to data
 if (data->references != 1)
 {
 data->references -= 1;
 }
 else
 {
 delete data;
 }
 data = temp;
 return *this;
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 bigint_kernel_2::
 operator/ (
 const bigint_kernel_2& rhs
 ) const
 {
 data_record* temp = new data_record(data->digits_used+slack);
 data_record* remainder;
 try {
 remainder = new data_record(data->digits_used+slack); 
 } catch (...) { delete temp; throw; }
 long_div(data,rhs.data,temp,remainder);
 delete remainder;
 
 return bigint_kernel_2(temp,0);
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator/= (
 const bigint_kernel_2& rhs
 )
 {
 data_record* temp = new data_record(data->digits_used+slack);
 data_record* remainder;
 try {
 remainder = new data_record(data->digits_used+slack);
 } catch (...) { delete temp; throw; } 
 long_div(data,rhs.data,temp,remainder);
 // check if there are other references to data
 if (data->references != 1)
 {
 data->references -= 1;
 }
 // if there are no references to data then it must be deleted
 else
 {
 delete data;
 }
 data = temp;
 delete remainder;
 
 return *this;
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 bigint_kernel_2::
 operator% (
 const bigint_kernel_2& rhs
 ) const
 {
 data_record* temp = new data_record(data->digits_used+slack);
 data_record* remainder;
 try {
 remainder = new data_record(data->digits_used+slack);
 } catch (...) { delete temp; throw; }
 long_div(data,rhs.data,temp,remainder);
 delete temp;
 return bigint_kernel_2(remainder,0);
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator%= (
 const bigint_kernel_2& rhs
 )
 {
 data_record* temp = new data_record(data->digits_used+slack);
 data_record* remainder;
 try {
 remainder = new data_record(data->digits_used+slack);
 } catch (...) { delete temp; throw; }
 long_div(data,rhs.data,temp,remainder);
 // check if there are other references to data
 if (data->references != 1)
 {
 data->references -= 1;
 }
 // if there are no references to data then it must be deleted
 else
 {
 delete data;
 }
 data = remainder;
 delete temp;
 return *this;
 }
// ----------------------------------------------------------------------------------------
 bool bigint_kernel_2::
 operator < (
 const bigint_kernel_2& rhs
 ) const
 {
 return is_less_than(data,rhs.data);
 }
// ----------------------------------------------------------------------------------------
 bool bigint_kernel_2::
 operator == (
 const bigint_kernel_2& rhs
 ) const
 {
 return is_equal_to(data,rhs.data);
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator= (
 const bigint_kernel_2& rhs
 )
 {
 if (this == &rhs)
 return *this;
 // if we have the only reference to our data then delete it
 if (data->references == 1)
 {
 delete data;
 data = rhs.data;
 data->references += 1;
 }
 else
 {
 data->references -= 1;
 data = rhs.data;
 data->references += 1;
 }
 return *this;
 }
// ----------------------------------------------------------------------------------------
 std::ostream& operator<< (
 std::ostream& out_,
 const bigint_kernel_2& rhs
 )
 {
 std::ostream out(out_.rdbuf());
 typedef bigint_kernel_2 bigint;
 
 bigint::data_record* temp = new bigint::data_record(*rhs.data,0);
 // get a char array big enough to hold the number in ascii format
 char* str;
 try {
 str = new char[(rhs.data->digits_used)*5+10];
 } catch (...) { delete temp; throw; }
 char* str_start = str;
 str += (rhs.data->digits_used)*5+9;
 *str = 0; --str;
 uint16 remainder;
 rhs.short_div(temp,10000,temp,remainder);
 // pull the digits out of remainder
 char a = remainder % 10 + '0';
 remainder /= 10;
 char b = remainder % 10 + '0';
 remainder /= 10;
 char c = remainder % 10 + '0';
 remainder /= 10;
 char d = remainder % 10 + '0';
 remainder /= 10;
 *str = a; --str;
 *str = b; --str;
 *str = c; --str;
 *str = d; --str;
 // keep looping until temp represents zero
 while (temp->digits_used != 1 || *(temp->number) != 0)
 { 
 rhs.short_div(temp,10000,temp,remainder);
 // pull the digits out of remainder
 char a = remainder % 10 + '0';
 remainder /= 10;
 char b = remainder % 10 + '0';
 remainder /= 10;
 char c = remainder % 10 + '0';
 remainder /= 10;
 char d = remainder % 10 + '0';
 remainder /= 10;
 *str = a; --str;
 *str = b; --str;
 *str = c; --str;
 *str = d; --str; 
 }
 // throw away and extra leading zeros
 ++str;
 if (*str == '0')
 ++str;
 if (*str == '0')
 ++str;
 if (*str == '0')
 ++str;
 
 out << str;
 delete [] str_start;
 delete temp;
 return out_;
 }
// ----------------------------------------------------------------------------------------
 std::istream& operator>> (
 std::istream& in_,
 bigint_kernel_2& rhs
 )
 {
 std::istream in(in_.rdbuf());
 // ignore any leading whitespaces
 while (in.peek() == ' ' || in.peek() == '\t' || in.peek() == '\n')
 {
 in.get();
 }
 // if the first digit is not an integer then this is an error
 if ( !(in.peek() >= '0' && in.peek() <= '9'))
 {
 in_.clear(std::ios::failbit);
 return in_;
 }
 int num_read;
 bigint_kernel_2 temp;
 do
 {
 // try to get 4 chars from in
 num_read = 1;
 char a = 0;
 char b = 0; 
 char c = 0;
 char d = 0;
 if (in.peek() >= '0' && in.peek() <= '9')
 {
 num_read *= 10;
 a = in.get();
 }
 if (in.peek() >= '0' && in.peek() <= '9')
 {
 num_read *= 10;
 b = in.get();
 }
 if (in.peek() >= '0' && in.peek() <= '9')
 {
 num_read *= 10;
 c = in.get();
 }
 if (in.peek() >= '0' && in.peek() <= '9')
 {
 num_read *= 10;
 d = in.get();
 }
 
 // merge the for digits into an uint16
 uint16 num = 0;
 if (a != 0)
 {
 num = a - '0';
 }
 if (b != 0)
 {
 num *= 10;
 num += b - '0';
 }
 if (c != 0)
 {
 num *= 10;
 num += c - '0';
 }
 if (d != 0)
 {
 num *= 10;
 num += d - '0';
 }
 if (num_read != 1)
 {
 // shift the digits in temp left by the number of new digits we just read
 temp *= num_read;
 // add in new digits
 temp += num;
 }
 } while (num_read == 10000);
 rhs = temp;
 return in_;
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 operator+ (
 uint16 lhs,
 const bigint_kernel_2& rhs
 )
 {
 typedef bigint_kernel_2 bigint;
 bigint::data_record* temp = new bigint::data_record
 (rhs.data->digits_used+rhs.slack);
 rhs.short_add(rhs.data,lhs,temp);
 return bigint_kernel_2(temp,0);
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 operator+ (
 const bigint_kernel_2& lhs,
 uint16 rhs
 )
 {
 typedef bigint_kernel_2 bigint;
 bigint::data_record* temp = new bigint::data_record
 (lhs.data->digits_used+lhs.slack);
 lhs.short_add(lhs.data,rhs,temp);
 return bigint_kernel_2(temp,0);
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator+= (
 uint16 rhs
 )
 {
 // if there are other references to this data 
 if (data->references != 1)
 {
 data_record* temp = new data_record(data->digits_used+slack);
 data->references -= 1; 
 short_add(data,rhs,temp);
 data = temp;
 }
 // or if we need to enlarge data then do so
 else if (data->digits_used == data->size)
 {
 data_record* temp = new data_record(data->digits_used+slack);
 short_add(data,rhs,temp);
 delete data;
 data = temp;
 }
 // or if there is plenty of space and no references
 else
 {
 short_add(data,rhs,data);
 }
 return *this;
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 operator- (
 uint16 lhs,
 const bigint_kernel_2& rhs
 )
 {
 typedef bigint_kernel_2 bigint;
 bigint::data_record* temp = new bigint::data_record(rhs.slack);
 *(temp->number) = lhs - *(rhs.data->number);
 return bigint_kernel_2(temp,0);
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 operator- (
 const bigint_kernel_2& lhs,
 uint16 rhs
 )
 {
 typedef bigint_kernel_2 bigint;
 bigint::data_record* temp = new bigint::data_record
 (lhs.data->digits_used+lhs.slack);
 lhs.short_sub(lhs.data,rhs,temp);
 return bigint_kernel_2(temp,0);
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator-= (
 uint16 rhs
 )
 {
 // if there are other references to this data 
 if (data->references != 1)
 {
 data_record* temp = new data_record(data->digits_used+slack);
 data->references -= 1;
 short_sub(data,rhs,temp);
 data = temp;
 }
 else
 {
 short_sub(data,rhs,data);
 }
 return *this;
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 operator* (
 uint16 lhs,
 const bigint_kernel_2& rhs
 )
 {
 typedef bigint_kernel_2 bigint;
 bigint::data_record* temp = new bigint::data_record
 (rhs.data->digits_used+rhs.slack);
 rhs.short_mul(rhs.data,lhs,temp);
 return bigint_kernel_2(temp,0); 
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 operator* (
 const bigint_kernel_2& lhs,
 uint16 rhs
 )
 {
 typedef bigint_kernel_2 bigint;
 bigint::data_record* temp = new bigint::data_record
 (lhs.data->digits_used+lhs.slack);
 lhs.short_mul(lhs.data,rhs,temp);
 return bigint_kernel_2(temp,0); 
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator*= (
 uint16 rhs
 )
 {
 // if there are other references to this data 
 if (data->references != 1)
 {
 data_record* temp = new data_record(data->digits_used+slack);
 data->references -= 1;
 short_mul(data,rhs,temp);
 data = temp;
 }
 // or if we need to enlarge data
 else if (data->digits_used == data->size)
 {
 data_record* temp = new data_record(data->digits_used+slack);
 short_mul(data,rhs,temp);
 delete data;
 data = temp;
 }
 else
 {
 short_mul(data,rhs,data);
 }
 return *this;
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 operator/ (
 uint16 lhs,
 const bigint_kernel_2& rhs
 )
 {
 typedef bigint_kernel_2 bigint;
 bigint::data_record* temp = new bigint::data_record(rhs.slack);
 // if rhs might not be bigger than lhs
 if (rhs.data->digits_used == 1)
 {
 *(temp->number) = lhs/ *(rhs.data->number);
 }
 
 return bigint_kernel_2(temp,0); 
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 operator/ (
 const bigint_kernel_2& lhs,
 uint16 rhs
 )
 {
 typedef bigint_kernel_2 bigint;
 bigint::data_record* temp = new bigint::data_record
 (lhs.data->digits_used+lhs.slack);
 uint16 remainder;
 lhs.short_div(lhs.data,rhs,temp,remainder);
 return bigint_kernel_2(temp,0); 
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator/= (
 uint16 rhs
 )
 {
 uint16 remainder;
 // if there are other references to this data
 if (data->references != 1)
 {
 data_record* temp = new data_record(data->digits_used+slack);
 data->references -= 1;
 short_div(data,rhs,temp,remainder); 
 data = temp;
 }
 else
 {
 short_div(data,rhs,data,remainder);
 }
 return *this;
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 operator% (
 uint16 lhs,
 const bigint_kernel_2& rhs
 )
 {
 typedef bigint_kernel_2 bigint;
 // temp is zero by default
 bigint::data_record* temp = new bigint::data_record(rhs.slack);
 if (rhs.data->digits_used == 1)
 {
 // if rhs is just an uint16 inside then perform the modulus
 *(temp->number) = lhs % *(rhs.data->number);
 }
 else
 {
 // if rhs is bigger than lhs then the answer is lhs
 *(temp->number) = lhs;
 }
 
 return bigint_kernel_2(temp,0); 
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 operator% (
 const bigint_kernel_2& lhs,
 uint16 rhs
 )
 {
 typedef bigint_kernel_2 bigint;
 bigint::data_record* temp = new bigint::data_record(lhs.data->digits_used+lhs.slack);
 uint16 remainder;
 lhs.short_div(lhs.data,rhs,temp,remainder);
 temp->digits_used = 1;
 *(temp->number) = remainder;
 return bigint_kernel_2(temp,0); 
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator%= (
 uint16 rhs
 )
 {
 uint16 remainder;
 // if there are other references to this data
 if (data->references != 1)
 {
 data_record* temp = new data_record(data->digits_used+slack);
 data->references -= 1;
 short_div(data,rhs,temp,remainder);
 data = temp;
 }
 else
 {
 short_div(data,rhs,data,remainder);
 }
 data->digits_used = 1;
 *(data->number) = remainder;
 return *this;
 }
// ----------------------------------------------------------------------------------------
 bool operator < (
 uint16 lhs,
 const bigint_kernel_2& rhs
 )
 {
 return (rhs.data->digits_used > 1 || lhs < *(rhs.data->number) );
 }
// ----------------------------------------------------------------------------------------
 bool operator < (
 const bigint_kernel_2& lhs,
 uint16 rhs
 )
 {
 return (lhs.data->digits_used == 1 && *(lhs.data->number) < rhs);
 }
// ----------------------------------------------------------------------------------------
 bool operator == (
 const bigint_kernel_2& lhs,
 uint16 rhs
 )
 {
 return (lhs.data->digits_used == 1 && *(lhs.data->number) == rhs);
 }
// ----------------------------------------------------------------------------------------
 bool operator == (
 uint16 lhs,
 const bigint_kernel_2& rhs
 )
 {
 return (rhs.data->digits_used == 1 && *(rhs.data->number) == lhs);
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator= (
 uint16 rhs
 )
 {
 // check if there are other references to our data
 if (data->references != 1)
 {
 data->references -= 1;
 try {
 data = new data_record(slack);
 } catch (...) { data->references += 1; throw; }
 }
 else
 {
 data->digits_used = 1;
 }
 
 *(data->number) = rhs;
 return *this;
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator++ (
 )
 {
 // if there are other references to this data then make a copy of it
 if (data->references != 1)
 {
 data_record* temp = new data_record(data->digits_used+slack);
 data->references -= 1;
 increment(data,temp);
 data = temp;
 }
 // or if we need to enlarge data then do so
 else if (data->digits_used == data->size)
 {
 data_record* temp = new data_record(data->digits_used+slack);
 increment(data,temp);
 delete data;
 data = temp;
 }
 else
 {
 increment(data,data);
 }
 return *this;
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 bigint_kernel_2::
 operator++ (
 int
 )
 {
 data_record* temp; // this is the copy of temp we will return in the end
 
 data_record* temp2 = new data_record(data->digits_used+slack);
 increment(data,temp2);
 
 temp = data;
 data = temp2;
 return bigint_kernel_2(temp,0);
 }
// ----------------------------------------------------------------------------------------
 bigint_kernel_2& bigint_kernel_2::
 operator-- (
 )
 {
 // if there are other references to this data 
 if (data->references != 1)
 {
 data_record* temp = new data_record(data->digits_used+slack);
 data->references -= 1;
 decrement(data,temp);
 data = temp;
 }
 else
 {
 decrement(data,data);
 }
 return *this;
 }
// ----------------------------------------------------------------------------------------
 const bigint_kernel_2 bigint_kernel_2::
 operator-- (
 int
 )
 {
 data_record* temp; // this is the copy of temp we will return in the end
 
 data_record* temp2 = new data_record(data->digits_used+slack);
 decrement(data,temp2);
 
 temp = data;
 data = temp2;
 return bigint_kernel_2(temp,0);
 }
// ----------------------------------------------------------------------------------------
// ----------------------------------------------------------------------------------------
 // private member function definitions
// ----------------------------------------------------------------------------------------
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 short_add (
 const data_record* data, 
 uint16 value,
 data_record* result
 ) const
 {
 // put value into the carry part of temp
 uint32 temp = value;
 temp <<= 16;
 
 const uint16* number = data->number;
 const uint16* end = number + data->digits_used; // one past the end of number
 uint16* r = result->number;
 while (number != end)
 {
 // add *number and the current carry
 temp = *number + (temp>>16);
 // put the low word of temp into *r
 *r = static_cast<uint16>(temp & 0xFFFF);
 ++number;
 ++r;
 }
 // if there is a final carry
 if ((temp>>16) != 0)
 {
 result->digits_used = data->digits_used + 1;
 // store the carry in the most significant digit of the result
 *r = static_cast<uint16>(temp>>16); 
 }
 else
 {
 result->digits_used = data->digits_used;
 }
 }
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 short_sub (
 const data_record* data, 
 uint16 value,
 data_record* result
 ) const
 {
 
 const uint16* number = data->number;
 const uint16* end = number + data->digits_used - 1;
 uint16* r = result->number;
 uint32 temp = *number - value;
 // put the low word of temp into *data 
 *r = static_cast<uint16>(temp & 0xFFFF);
 
 while (number != end)
 {
 ++number;
 ++r;
 // subtract the carry from *number 
 temp = *number - (temp>>31);
 // put the low word of temp into *r 
 *r = static_cast<uint16>(temp & 0xFFFF);
 }
 // if we lost a digit in the subtraction
 if (*r == 0)
 {
 if (data->digits_used == 1)
 result->digits_used = 1;
 else
 result->digits_used = data->digits_used - 1;
 }
 else
 {
 result->digits_used = data->digits_used;
 }
 }
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 short_mul (
 const data_record* data, 
 uint16 value,
 data_record* result
 ) const
 {
 
 uint32 temp = 0;
 const uint16* number = data->number; 
 uint16* r = result->number;
 const uint16* end = r + data->digits_used;
 while ( r != end)
 {
 // multiply *data and value and add in the carry
 temp = *number*(uint32)value + (temp>>16);
 // put the low word of temp into *data
 *r = static_cast<uint16>(temp & 0xFFFF);
 ++number;
 ++r;
 }
 // if there is a final carry
 if ((temp>>16) != 0)
 {
 result->digits_used = data->digits_used + 1;
 // put the final carry into the most significant digit of the result
 *r = static_cast<uint16>(temp>>16);
 }
 else
 {
 result->digits_used = data->digits_used;
 }
 }
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 short_div (
 const data_record* data, 
 uint16 value, 
 data_record* result,
 uint16& rem
 ) const
 {
 
 uint16 remainder = 0;
 uint32 temp;
 
 const uint16* number = data->number + data->digits_used - 1;
 const uint16* end = number - data->digits_used;
 uint16* r = result->number + data->digits_used - 1;
 // if we are losing a digit in this division
 if (*number < value)
 {
 if (data->digits_used == 1)
 result->digits_used = 1;
 else
 result->digits_used = data->digits_used - 1;
 }
 else
 {
 result->digits_used = data->digits_used;
 }
 // perform the actual division
 while (number != end)
 {
 
 temp = *number + (((uint32)remainder)<<16);
 *r = static_cast<uint16>(temp/value);
 remainder = static_cast<uint16>(temp%value);
 --number;
 --r;
 }
 rem = remainder;
 }
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 long_add (
 const data_record* lhs,
 const data_record* rhs,
 data_record* result
 ) const
 {
 // put value into the carry part of temp
 uint32 temp=0; 
 uint16* min_num; // the number with the least digits used
 uint16* max_num; // the number with the most digits used
 uint16* min_end; // one past the end of min_num
 uint16* max_end; // one past the end of max_num
 uint16* r = result->number;
 uint32 max_digits_used;
 if (lhs->digits_used < rhs->digits_used)
 {
 max_digits_used = rhs->digits_used;
 min_num = lhs->number;
 max_num = rhs->number;
 min_end = min_num + lhs->digits_used;
 max_end = max_num + rhs->digits_used;
 }
 else
 {
 max_digits_used = lhs->digits_used;
 min_num = rhs->number;
 max_num = lhs->number;
 min_end = min_num + rhs->digits_used;
 max_end = max_num + lhs->digits_used;
 }
 
 while (min_num != min_end)
 {
 // add *min_num, *max_num and the current carry
 temp = *min_num + *max_num + (temp>>16);
 // put the low word of temp into *r
 *r = static_cast<uint16>(temp & 0xFFFF);
 ++min_num;
 ++max_num;
 ++r;
 }
 while (max_num != max_end)
 {
 // add *max_num and the current carry
 temp = *max_num + (temp>>16);
 // put the low word of temp into *r
 *r = static_cast<uint16>(temp & 0xFFFF);
 ++max_num;
 ++r;
 } 
 // check if there was a final carry
 if ((temp>>16) != 0)
 {
 result->digits_used = max_digits_used + 1;
 // put the carry into the most significant digit in the result
 *r = static_cast<uint16>(temp>>16);
 }
 else
 {
 result->digits_used = max_digits_used;
 }
 
 }
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 long_sub (
 const data_record* lhs,
 const data_record* rhs,
 data_record* result
 ) const
 {
 const uint16* number1 = lhs->number;
 const uint16* number2 = rhs->number;
 const uint16* end = number2 + rhs->digits_used;
 uint16* r = result->number;
 uint32 temp =0;
 
 while (number2 != end)
 {
 // subtract *number2 from *number1 and then subtract any carry
 temp = *number1 - *number2 - (temp>>31);
 // put the low word of temp into *r 
 *r = static_cast<uint16>(temp & 0xFFFF);
 ++number1;
 ++number2;
 ++r;
 }
 end = lhs->number + lhs->digits_used;
 while (number1 != end)
 {
 // subtract the carry from *number1 
 temp = *number1 - (temp>>31);
 // put the low word of temp into *r 
 *r = static_cast<uint16>(temp & 0xFFFF);
 ++number1;
 ++r;
 }
 result->digits_used = lhs->digits_used;
 // adjust the number of digits used appropriately 
 --r;
 while (*r == 0 && result->digits_used > 1)
 {
 --r;
 --result->digits_used;
 }
 }
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 long_div (
 const data_record* lhs,
 const data_record* rhs,
 data_record* result,
 data_record* remainder
 ) const
 {
 // zero result
 result->digits_used = 1;
 *(result->number) = 0;
 uint16* a;
 uint16* b;
 uint16* end;
 // copy lhs into remainder
 remainder->digits_used = lhs->digits_used;
 a = remainder->number;
 end = a + remainder->digits_used;
 b = lhs->number;
 while (a != end)
 {
 *a = *b;
 ++a;
 ++b;
 }
 // if rhs is bigger than lhs then result == 0 and remainder == lhs
 // so then we can quit right now
 if (is_less_than(lhs,rhs))
 {
 return; 
 }
 // make a temporary number
 data_record temp(lhs->digits_used + slack);
 // shift rhs left until it is one shift away from being larger than lhs and
 // put the number of left shifts necessary into shifts
 uint32 shifts; 
 shifts = (lhs->digits_used - rhs->digits_used) * 16;
 shift_left(rhs,&temp,shifts);
 // while (lhs > temp)
 while (is_less_than(&temp,lhs))
 {
 shift_left(&temp,&temp,1);
 ++shifts;
 }
 // make sure lhs isn't smaller than temp
 while (is_less_than(lhs,&temp))
 {
 shift_right(&temp,&temp);
 --shifts;
 }
 
 
 // we want to execute the loop shifts +1 times
 ++shifts;
 while (shifts != 0)
 {
 shift_left(result,result,1);
 // if (temp <= remainder)
 if (!is_less_than(remainder,&temp))
 {
 long_sub(remainder,&temp,remainder);
 
 // increment result
 uint16* r = result->number;
 uint16* end = r + result->digits_used;
 while (true)
 {
 ++(*r);
 // if there was no carry then we are done
 if (*r != 0)
 break;
 ++r;
 // if we hit the end of r and there is still a carry then
 // the next digit of r is 1 and there is one more digit used
 if (r == end)
 {
 *r = 1;
 ++(result->digits_used);
 break;
 }
 }
 }
 shift_right(&temp,&temp);
 --shifts;
 }
 
 
 }
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 long_mul (
 const data_record* lhs,
 const data_record* rhs,
 data_record* result
 ) const
 {
 // if one of the numbers is small then use this simple but O(n^2) algorithm
 if (std::min(lhs->digits_used, rhs->digits_used) < 10)
 {
 // make result be zero
 result->digits_used = 1;
 *(result->number) = 0;
 const data_record* aa;
 const data_record* bb;
 if (lhs->digits_used < rhs->digits_used)
 {
 // make copies of lhs and rhs and give them an appropriate amount of
 // extra memory so there won't be any overflows
 aa = lhs;
 bb = rhs;
 }
 else
 {
 // make copies of lhs and rhs and give them an appropriate amount of
 // extra memory so there won't be any overflows
 aa = rhs;
 bb = lhs;
 }
 // copy the larger(approximately) of lhs and rhs into b
 data_record b(*bb,aa->digits_used+slack); 
 uint32 shift_value = 0;
 uint16* anum = aa->number;
 uint16* end = anum + aa->digits_used;
 while (anum != end )
 {
 uint16 bit = 0x0001;
 for (int i = 0; i < 16; ++i)
 {
 // if the specified bit of a is 1
 if ((*anum & bit) != 0)
 {
 shift_left(&b,&b,shift_value);
 shift_value = 0;
 long_add(&b,result,result);
 }
 ++shift_value;
 bit <<= 1;
 }
 ++anum; 
 }
 }
 else // else if both lhs and rhs are large then use the more complex 
 // O(n*logn) algorithm
 {
 uint32 size = 1;
 // make size a power of 2
 while (size < (lhs->digits_used + rhs->digits_used)*2) 
 {
 size *= 2;
 }
 // allocate some temporary space so we can do the FFT
 ct* a = new ct[size];
 ct* b; try {b = new ct[size]; } catch (...) { delete [] a; throw; }
 // load lhs into the a array. We are breaking the input number into 
 // 8bit chunks for the purpose of using this fft algorithm. The reason 
 // for this is so that we have smaller numbers coming out of the final 
 // ifft. This helps avoid overflow.
 for (uint32 i = 0; i < lhs->digits_used; ++i)
 {
 a[i*2] = ct((t)(lhs->number[i]&0xFF),0);
 a[i*2+1] = ct((t)(lhs->number[i]>>8),0);
 }
 for (uint32 i = lhs->digits_used*2; i < size; ++i)
 {
 a[i] = 0;
 }
 // load rhs into the b array
 for (uint32 i = 0; i < rhs->digits_used; ++i)
 {
 b[i*2] = ct((t)(rhs->number[i]&0xFF),0);
 b[i*2+1] = ct((t)(rhs->number[i]>>8),0);
 }
 for (uint32 i = rhs->digits_used*2; i < size; ++i)
 {
 b[i] = 0;
 }
 // perform the forward fft of a and b
 fft(a,size);
 fft(b,size);
 const double l = 1.0/size;
 // do the pointwise multiply of a and b and also apply the scale
 // factor in this loop too.
 for (unsigned long i = 0; i < size; ++i)
 {
 a[i] = l*a[i]*b[i];
 }
 // Now compute the inverse fft of the pointwise multiplication of a and b. 
 // This is basically the result. We just have to take care of any carries 
 // that should happen.
 ifft(a,size);
 // loop over the result and propagate any carries that need to take place. 
 // We will also be moving the resulting numbers into result->number at 
 // the same time.
 uint64 carry = 0;
 result->digits_used = 0;
 int zeros = 0;
 const uint32 len = lhs->digits_used + rhs->digits_used;
 for (unsigned long i = 0; i < len; ++i)
 {
 uint64 num1 = static_cast<uint64>(std::round(a[i*2].real()));
 num1 += carry;
 carry = 0;
 if (num1 > 255)
 {
 carry = num1 >> 8;
 num1 = (num1&0xFF);
 }
 uint64 num2 = static_cast<uint64>(std::round(a[i*2+1].real()));
 num2 += carry;
 carry = 0;
 if (num2 > 255)
 {
 carry = num2 >> 8;
 num2 = (num2&0xFF);
 }
 // put the new number into its final place
 num1 = (num2<<8) | num1;
 result->number[i] = static_cast<uint16>(num1);
 // keep track of the number of leading zeros
 if (num1 == 0)
 ++zeros;
 else
 zeros = 0;
 ++(result->digits_used);
 }
 // adjust digits_used so that it reflects the actual number
 // of non-zero digits in our representation.
 result->digits_used -= zeros;
 // if the result was zero then adjust the result accordingly
 if (result->digits_used == 0)
 {
 // make result be zero
 result->digits_used = 1;
 *(result->number) = 0;
 }
 // free all the temporary buffers
 delete [] a;
 delete [] b;
 }
 }
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 shift_left (
 const data_record* data,
 data_record* result,
 uint32 shift_amount
 ) const
 {
 uint32 offset = shift_amount/16;
 shift_amount &= 0xf; // same as shift_amount %= 16;
 uint16* r = result->number + data->digits_used + offset; // result
 uint16* end = data->number;
 uint16* s = end + data->digits_used; // source
 const uint32 temp = 16 - shift_amount;
 *r = (*(--s) >> temp);
 // set the number of digits used in the result
 // if the upper bits from *s were zero then don't count this first word
 if (*r == 0)
 {
 result->digits_used = data->digits_used + offset;
 }
 else
 {
 result->digits_used = data->digits_used + offset + 1;
 }
 --r;
 while (s != end)
 {
 *r = ((*s << shift_amount) | ( *(s-1) >> temp));
 --r;
 --s;
 }
 *r = *s << shift_amount;
 // now zero the rest of the result
 end = result->number;
 while (r != end)
 *(--r) = 0;
 }
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 shift_right (
 const data_record* data,
 data_record* result
 ) const
 {
 uint16* r = result->number; // result
 uint16* s = data->number; // source
 uint16* end = s + data->digits_used - 1;
 while (s != end)
 {
 *r = (*s >> 1) | (*(s+1) << 15);
 ++r;
 ++s;
 }
 *r = *s >> 1;
 // calculate the new number for digits_used
 if (*r == 0)
 {
 if (data->digits_used != 1)
 result->digits_used = data->digits_used - 1;
 else
 result->digits_used = 1;
 }
 else
 {
 result->digits_used = data->digits_used;
 }
 
 }
// ----------------------------------------------------------------------------------------
 bool bigint_kernel_2::
 is_less_than (
 const data_record* lhs,
 const data_record* rhs
 ) const
 {
 uint32 lhs_digits_used = lhs->digits_used;
 uint32 rhs_digits_used = rhs->digits_used;
 // if lhs is definitely less than rhs
 if (lhs_digits_used < rhs_digits_used )
 return true;
 // if lhs is definitely greater than rhs
 else if (lhs_digits_used > rhs_digits_used)
 return false;
 else 
 {
 uint16* end = lhs->number;
 uint16* l = end + lhs_digits_used;
 uint16* r = rhs->number + rhs_digits_used;
 
 while (l != end)
 {
 --l;
 --r;
 if (*l < *r)
 return true;
 else if (*l > *r)
 return false;
 }
 // at this point we know that they are equal
 return false;
 }
 }
// ----------------------------------------------------------------------------------------
 bool bigint_kernel_2::
 is_equal_to (
 const data_record* lhs,
 const data_record* rhs
 ) const
 {
 // if lhs and rhs are definitely not equal
 if (lhs->digits_used != rhs->digits_used )
 {
 return false;
 }
 else 
 { 
 uint16* l = lhs->number;
 uint16* r = rhs->number;
 uint16* end = l + lhs->digits_used;
 
 while (l != end)
 {
 if (*l != *r)
 return false;
 ++l;
 ++r;
 }
 // at this point we know that they are equal
 return true;
 }
 }
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 increment (
 const data_record* source,
 data_record* dest
 ) const
 {
 uint16* s = source->number;
 uint16* d = dest->number;
 uint16* end = s + source->digits_used;
 while (true)
 {
 *d = *s + 1;
 // if there was no carry then break out of the loop
 if (*d != 0)
 {
 dest->digits_used = source->digits_used;
 // copy the rest of the digits over to d
 ++d; ++s;
 while (s != end)
 {
 *d = *s;
 ++d;
 ++s;
 }
 break;
 }
 
 ++s; 
 // if we have hit the end of s and there was a carry up to this point
 // then just make the next digit 1 and add one to the digits used
 if (s == end)
 {
 ++d;
 dest->digits_used = source->digits_used + 1;
 *d = 1;
 break;
 }
 ++d;
 }
 }
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 decrement (
 const data_record* source,
 data_record* dest
 ) const
 {
 uint16* s = source->number;
 uint16* d = dest->number;
 uint16* end = s + source->digits_used;
 while (true)
 {
 *d = *s - 1;
 // if there was no carry then break out of the loop
 if (*d != 0xFFFF)
 {
 // if we lost a digit in the subtraction 
 if (*d == 0 && s+1 == end)
 {
 if (source->digits_used == 1)
 dest->digits_used = 1;
 else
 dest->digits_used = source->digits_used - 1;
 }
 else
 {
 dest->digits_used = source->digits_used;
 }
 break;
 }
 else
 {
 ++d;
 ++s;
 }
 }
 // copy the rest of the digits over to d
 ++d;
 ++s;
 while (s != end)
 {
 *d = *s;
 ++d;
 ++s;
 }
 }
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 fft (
 ct* data, 
 unsigned long len
 ) const 
 {
 const t pi2 = -2.0*3.1415926535897932384626433832795028841971693993751;
 const unsigned long half = len/2;
 std::vector<ct> twiddle_factors;
 twiddle_factors.resize(half);
 // compute the complex root of unity w
 const t temp = pi2/len;
 ct w = ct(std::cos(temp),std::sin(temp));
 ct w_pow = 1;
 // compute the twiddle factors
 for (std::vector<ct>::size_type j = 0; j < twiddle_factors.size(); ++j)
 {
 twiddle_factors[j] = w_pow; 
 w_pow *= w;
 }
 ct a, b;
 // now compute the decimation in frequency. This first
 // outer loop loops log2(len) number of times
 unsigned long skip = 1;
 for (unsigned long step = half; step != 0; step >>= 1)
 {
 // do blocks of butterflies in this loop
 for (unsigned long j = 0; j < len; j += step*2)
 {
 // do step butterflies
 for (unsigned long k = 0; k < step; ++k)
 {
 const unsigned long a_idx = j+k;
 const unsigned long b_idx = j+k+step;
 a = data[a_idx] + data[b_idx];
 b = (data[a_idx] - data[b_idx])*twiddle_factors[k*skip];
 data[a_idx] = a;
 data[b_idx] = b;
 }
 }
 skip *= 2;
 }
 }
// ----------------------------------------------------------------------------------------
 void bigint_kernel_2::
 ifft(
 ct* data, 
 unsigned long len
 ) const 
 {
 const t pi2 = 2.0*3.1415926535897932384626433832795028841971693993751;
 const unsigned long half = len/2;
 std::vector<ct> twiddle_factors;
 twiddle_factors.resize(half);
 // compute the complex root of unity w
 const t temp = pi2/len;
 ct w = ct(std::cos(temp),std::sin(temp));
 ct w_pow = 1;
 // compute the twiddle factors
 for (std::vector<ct>::size_type j = 0; j < twiddle_factors.size(); ++j)
 {
 twiddle_factors[j] = w_pow; 
 w_pow *= w;
 }
 ct a, b;
 // now compute the inverse decimation in frequency. This first
 // outer loop loops log2(len) number of times
 unsigned long skip = half;
 for (unsigned long step = 1; step <= half; step <<= 1)
 {
 // do blocks of butterflies in this loop
 for (unsigned long j = 0; j < len; j += step*2)
 {
 // do step butterflies
 for (unsigned long k = 0; k < step; ++k)
 {
 const unsigned long a_idx = j+k;
 const unsigned long b_idx = j+k+step;
 data[b_idx] *= twiddle_factors[k*skip];
 a = data[a_idx] + data[b_idx];
 b = data[a_idx] - data[b_idx];
 data[a_idx] = a;
 data[b_idx] = b;
 }
 }
 skip /= 2;
 }
 }
// ----------------------------------------------------------------------------------------
}
#endif // DLIB_BIGINT_KERNEL_2_CPp_

AltStyle によって変換されたページ (->オリジナル) /