/* $Id$ */

/*
 * Copyright (C) 2005-2006 Richard Braun
 * 
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or
 * (at your option) any later version.
 * 
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 * 
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 */

#include <string>
#include <string.h>
#include <iostream>

#include "mtx.h"
#include "mtxmatrix.h"
#include "mtxinteger.h"
#include "mtxrational.h"
#include "mtxexception.h"
#include "mtxdivbyzeroexception.h"

using namespace std;

unsigned int MtxInteger::base = MTX_DEFAULT_BASE;

MtxInteger::MtxInteger()
{
  mtx_int_t word;

  words = NULL;
  word = 0;
  setWords(&word, 1);
}

MtxInteger::MtxInteger(int integer)
{
  mtx_int_t word;

  words = NULL;
  word = integer;
  setWords(&word, 1);
}

MtxInteger::MtxInteger(const char *digits)
{
  size_t i, length;
  unsigned int n;
  bool negative;

  {
    mtx_int_t word;

    words = NULL;
    word = 0;
    setWords(&word, 1);
  }

  while (isspace(*digits))
    digits++;

  if (*digits == '-')
    {
      negative = true;
      digits++;
    }

  else
    {
      negative = false;

      if (*digits == '+')
        digits++;
    }

  if (*digits == '\0')
    throw MtxException("invalid string representation");

  length = strlen(digits);

  if (isPowerOfTwo(base))
    {
      unsigned int bits_per_digit;

      bits_per_digit = bitsPerDigit(base);

      for (i = length; i > 0; i--)
        {
          n = charToInt(digits[i - 1]);

          if (n >= base)
            throw MtxException("digit not allowed with this base");

          *this += MtxInteger(n) << (bits_per_digit * (length - i));
        }
    }

  else
    {
      MtxInteger tmp_base;

      tmp_base = base;

      for (i = length; i > 0; i--)
        {
          n = charToInt(digits[i - 1]);

          if (n >= base)
            throw MtxException("digit not allowed with this base");

          *this += tmp_base.power(length - i) * n;
        }
    }

  if (negative)
    *this = -*this;
}

MtxInteger::MtxInteger(const MtxInteger &integer)
{
  words = NULL;
  setWords(integer.words, integer.words_count);
}

MtxInteger::~MtxInteger()
{
  /*
   * words is always allocated.
   */
  delete[] words;
}

bool
MtxInteger::isNegative() const
{
  if (words[words_count - 1] & MTX_INT_HALF)
    return true;

  else
    return false;
}

void
MtxInteger::normalize()
{
  if (words_count == 1)
    return;

  else if (words[words_count - 1] != 0
      && words[words_count - 1] != ~((mtx_int_t)0))
    return;

  else
    {
      size_t i, new_words_count;
      mtx_int_t *tmp_words;
      bool negative;

      for (i = words_count - 2; i > 0; i--)
        if (words[i] != 0 && words[i] != ~((mtx_int_t)0))
          break;

      negative = isNegative();

      if ((!negative && (words[i] & MTX_INT_HALF) == 0)
          || (negative && (words[i] & MTX_INT_HALF) != 0))
        new_words_count = i + 1;

      else
        {
          new_words_count = i + 2;

          if (words_count == new_words_count)
            return;
        }

      tmp_words = new mtx_int_t[new_words_count];
      memcpy(tmp_words, words, sizeof(mtx_int_t) * new_words_count);
      delete[] words;
      words = tmp_words;
      words_count = new_words_count;
    }
}

void
MtxInteger::setWords(const mtx_int_t *words, size_t words_count)
{
  if (this->words != NULL)
    delete[] this->words;

  this->words = new mtx_int_t[words_count];
  memcpy(this->words, words, words_count * sizeof(mtx_int_t));
  this->words_count = words_count;
  normalize();
}

void
MtxInteger::unsigned_divide(const MtxInteger &divisor,
                            MtxInteger &quotient,
                            MtxInteger &remainder) const
{
  MtxInteger abs_dividend, abs_divisor;

  abs_dividend = abs();
  abs_divisor = divisor.abs();
  quotient = 0;
  remainder = 0;

  if (abs_divisor == 0)
    throw MtxDivByZeroException();

  else if (abs_dividend == 0)
    return;

  else if (abs_divisor == 1)
    quotient = abs_dividend;

  else if (abs_divisor == abs_dividend)
    quotient = 1;

  else if (abs_divisor > abs_dividend)
    remainder = abs_dividend;

  else
    {
      MtxInteger inf, sup, tmp;

      inf = 0;
      sup = abs_dividend;

      do
        {
          quotient = inf + ((sup - inf) >> 1);
          tmp = quotient * abs_divisor;

          if (tmp > abs_dividend)
            sup = quotient;

          else
            {
              remainder = abs_dividend - tmp;

              if (remainder < abs_divisor)
                return;

              inf = quotient;
            }
        }
      while(true);
    }
}

void
MtxInteger::divide(const MtxInteger &divisor,
                   MtxInteger &quotient,
                   MtxInteger &remainder) const
{
  bool negative;

  unsigned_divide(divisor, quotient, remainder);
  negative = isNegative();

  if (negative != divisor.isNegative())
    quotient = -quotient;

  if (negative)
    remainder = -remainder;
}

bool
MtxInteger::isPowerOfTwo(unsigned int n)
{
  while ((n & 1) == 0)
    n >>= 1;

  if (n == 1)
    return true;

  else
    return false;
}

unsigned int
MtxInteger::bitsPerDigit(unsigned int base)
{
  unsigned int bits, mask;

  if (base < 2)
    base = 2;

  mask = 1 << ((sizeof(mask) << 3) - 1);

  while ((base & mask) == 0)
    mask >>= 1;

  base &= mask;
  bits = 0;

  while (base != 1)
    {
      bits++;
      base >>= 1;
    }

  return bits;
}

char
MtxInteger::intToChar(unsigned int n)
{
  switch (n)
    {
      case 0 :
        return '0';
      case 1 :
        return '1';
      case 2 :
        return '2';
      case 3 :
        return '3';
      case 4 :
        return '4';
      case 5 :
        return '5';
      case 6 :
        return '6';
      case 7 :
        return '7';
      case 8 :
        return '8';
      case 9 :
        return '9';
      case 10 :
        return 'a';
      case 11 :
        return 'b';
      case 12 :
        return 'c';
      case 13 :
        return 'd';
      case 14 :
        return 'e';
      case 15 :
        return 'f';
      default :
        throw MtxException("value out of range");
    }
}

unsigned int
MtxInteger::charToInt(char c)
{
  switch (c)
    {
      case '0' :
        return 0;
      case '1' :
        return 1;
      case '2' :
        return 2;
      case '3' :
        return 3;
      case '4' :
        return 4;
      case '5' :
        return 5;
      case '6' :
        return 6;
      case '7' :
        return 7;
      case '8' :
        return 8;
      case '9' :
        return 9;
      case 'a' : case 'A' :
        return 10;
      case 'b' : case 'B' :
        return 11;
      case 'c' : case 'C' :
        return 12;
      case 'd' : case 'D' :
        return 13;
      case 'e' : case 'E' :
        return 14;
      case 'f' : case 'F' :
        return 15;
      default :
        throw MtxException("character out of range");
    }
}

MtxInteger &
MtxInteger::operator=(const MtxInteger &integer)
{
  if (this != &integer)
    setWords(integer.words, integer.words_count);

  return *this;
}

MtxInteger
MtxInteger::operator+() const
{
  return *this;
}

MtxInteger
MtxInteger::operator-() const
{
  MtxInteger new_integer(*this);
  size_t i;

  for (i = 0; i < new_integer.words_count; i++)
    new_integer.words[i] = ~new_integer.words[i];

  new_integer++;
  return new_integer;
}

bool
MtxInteger::operator==(const MtxInteger &integer) const
{
  if (this == &integer)
    return true;

  else if (words_count != integer.words_count)
    return false;

  else
    {
      size_t i;

      for (i = 0; i < words_count; i++)
        if (words[i] != integer.words[i])
          return false;

      return true;
    }
}

bool
MtxInteger::operator!=(const MtxInteger &integer) const
{
  return !(*this == integer);
}

bool
MtxInteger::operator>(const MtxInteger &integer) const
{
  if (this == &integer)
    return false;

  else if (*this == integer)
    return false;

  else
    {
      bool negative;

      negative = isNegative();

      if (negative != integer.isNegative())
        {
          if (negative)
            return false;

          else
            return true;
        }

      else
        {
          if (words_count > integer.words_count)
            {
              if (negative)
                return false;

              else
                return true;
            }

          else if (words_count < integer.words_count)
            {
              if (negative)
                return true;

              else
                return false;
            }

          else
            {
              size_t i;

              for (i = words_count; i > 0; i--)
                {
                  if (words[i - 1] > integer.words[i - 1])
                    return true;

                  else if (words[i - 1] < integer.words[i - 1])
                    return false;
                }

              /*
               * Normally never reached as equality is cheched at the
               * beginning of the method.
               */
              cerr << "warning: running unreachable code." << endl;
              return false;
            }
        }
    }
}

bool
MtxInteger::operator>=(const MtxInteger &integer) const
{
  return (*this > integer || *this == integer);
}

bool
MtxInteger::operator<(const MtxInteger &integer) const
{
  return integer > *this;
}

bool
MtxInteger::operator<=(const MtxInteger &integer) const
{
  return (*this < integer || *this == integer);
}

MtxInteger
MtxInteger::operator>>(unsigned int bits) const
{
  MtxInteger new_integer;
  unsigned int i;
  size_t j;

  new_integer = *this;

  for (i = 0; i < bits; i++)
    {
      new_integer.words[0] >>= 1;

      for (j = 1; j < new_integer.words_count; j++)
        {
          new_integer.words[j - 1] |= new_integer.words[j]
                                   << (MTX_INT_BITS - 1);
          new_integer.words[j] >>= 1;
        }
    }

  new_integer.normalize();
  return new_integer;
}

MtxInteger &
MtxInteger::operator>>=(unsigned int bits)
{
  *this = *this >> bits;
  return *this;
}

MtxInteger
MtxInteger::operator<<(unsigned int bits) const
{
  bool overflow_previous, overflow_next;
  MtxInteger new_integer;
  unsigned int i;
  size_t j;

  new_integer = *this;

  for (i = 0; i < bits; i++)
    {
      overflow_previous = false;

      for (j = 0; j < new_integer.words_count; j++)
        {
          overflow_next = new_integer.words[j] & MTX_INT_HALF;
          new_integer.words[j] <<= 1;

          if (overflow_previous)
            new_integer.words[j] |= 1;

          overflow_previous = overflow_next;
        }

      if (overflow_previous)
        {
          mtx_int_t *tmp_words;

          tmp_words = new mtx_int_t[new_integer.words_count + 1];
          memcpy(tmp_words, new_integer.words,
                 sizeof(mtx_int_t) * new_integer.words_count);
          delete[] new_integer.words;
          new_integer.words = tmp_words;
          new_integer.words[new_integer.words_count] = 1;
          new_integer.words_count++;
        }
    }

  /*
   * Sign must be preserved for positive integers only.
   */
  if (!isNegative() && new_integer.isNegative())
    {
      mtx_int_t *tmp_words;

      tmp_words = new mtx_int_t[new_integer.words_count + 1];
      memcpy(tmp_words, new_integer.words,
             sizeof(mtx_int_t) * new_integer.words_count);
      delete[] new_integer.words;
      new_integer.words = tmp_words;
      new_integer.words[new_integer.words_count] = 0;
      new_integer.words_count++;
    }

  return new_integer;
}

MtxInteger &
MtxInteger::operator<<=(unsigned int bits)
{
  *this = *this << bits;
  return *this;
}

MtxInteger
MtxInteger::operator+(const MtxInteger &integer) const
{
  bool overflow, negative, integer_negative;
  mtx_int_t word = 0, word1, word2;
  MtxInteger new_integer;
  size_t i, words_count;

  words_count = (this->words_count > integer.words_count)
              ? this->words_count
              : integer.words_count;
  delete[] new_integer.words;
  new_integer.words = new mtx_int_t[words_count];
  overflow = false;
  negative = isNegative();
  integer_negative = integer.isNegative();

  for (i = 0; i < words_count; i++)
    {
      if (i < this->words_count)
        word1 = words[i];

      else
        word1 = negative ? ~((mtx_int_t)0) : 0;

      if (i < integer.words_count)
        word2 = integer.words[i];

      else
        word2 = integer_negative ? ~((mtx_int_t)0) : 0;

      word = word1 + word2;

      if (overflow)
        word++;

      overflow = (word < word1) ? true : false;
      new_integer.words[i] = word;
    }

  /*
   * Handling the most significant word is tricky.
   */
  if ((negative == integer_negative)
      && (((word & MTX_INT_HALF) && !isNegative())
          || (!(word & MTX_INT_HALF) && overflow)))
    {
      mtx_int_t *tmp_words;

      tmp_words = new mtx_int_t[words_count + 1];
      memcpy(tmp_words, new_integer.words,
             sizeof(mtx_int_t) * words_count);
      tmp_words[words_count] = overflow ? ~((mtx_int_t)0) : 0;
      delete[] new_integer.words;
      new_integer.words = tmp_words;
      words_count++;
    }

  new_integer.words_count = words_count;
  new_integer.normalize();
  return new_integer;
}

MtxInteger &
MtxInteger::operator++()
{
  *this = *this + 1;
  return *this;
}

MtxInteger &
MtxInteger::operator++(int a)
{
  *this = *this + 1;
  return *this;
}

MtxInteger &
MtxInteger::operator+=(const MtxInteger &integer)
{
  *this = *this + integer;
  return *this;
}

MtxInteger
MtxInteger::operator-(const MtxInteger &integer) const
{
  return *this + (-integer);
}

MtxInteger &
MtxInteger::operator--()
{
  *this = *this - 1;
  return *this;
}

MtxInteger &
MtxInteger::operator--(int a)
{
  *this = *this - 1;
  return *this;
}

MtxInteger &
MtxInteger::operator-=(const MtxInteger &integer)
{
  *this = *this - integer;
  return *this;
}

MtxInteger
MtxInteger::operator*(const MtxInteger &integer) const
{
  MtxInteger a, b, new_integer, tmp;
  size_t i, bits_count;

  a = abs();
  b = integer.abs();
  new_integer = 0;
  bits_count = b.words_count * MTX_INT_BITS;

  for (i = 0; i < bits_count; i++)
    {
      if (((b >> i).words[0] & 1) == 0)
        tmp = 0;

      else
        tmp = a << i;

      new_integer += tmp;
    }

  if (isNegative() == integer.isNegative())
    return new_integer;

  else
    return (-new_integer);
}

MtxInteger &
MtxInteger::operator*=(const MtxInteger &integer)
{
  *this = *this * integer;
  return *this;
}

MtxInteger
MtxInteger::operator/(const MtxInteger &integer) const
{
  MtxInteger quotient, remainder;

  divide(integer, quotient, remainder);
  return quotient;
}

MtxInteger &
MtxInteger::operator/=(const MtxInteger &integer)
{
  *this = *this / integer;
  return *this;
}

MtxInteger
MtxInteger::operator%(const MtxInteger &integer) const
{
  MtxInteger quotient, remainder;

  divide(integer, quotient, remainder);
  return remainder;
}

MtxInteger &
MtxInteger::operator%=(const MtxInteger &integer)
{
  *this = *this % integer;
  return *this;
}

ostream &
operator<<(ostream &cout, const MtxInteger &integer)
{
  cout << integer.toString();
  return cout;
}

MtxInteger
MtxInteger::power(unsigned int power) const
{
  MtxInteger new_integer;
  unsigned int i;

  new_integer = 1;

  for (i = 0; i < power; i++)
    new_integer *= *this;

  return new_integer;
}

MtxInteger
MtxInteger::abs() const
{
  if (isNegative())
    return -(*this);

  else
    return *this;
}

string
MtxInteger::toString() const
{
  return toString(base);
}

string
MtxInteger::toString(unsigned int base) const
{
  unsigned int bits_per_digit;
  size_t i, sign, length;
  MtxInteger integer;
  bool negative;
  string str;
  char *s;

  if (base < 2 || base > MTX_MAX_BASE)
    base = MtxInteger::base;

  bits_per_digit = bitsPerDigit(base);
  negative = isNegative();
  sign = negative ? 1 : 0;
  integer = abs();
  length = sign + ((words_count * MTX_INT_BITS) / bits_per_digit);
  s = new char[length + 1];
  s[length] = '\0';
  i = length - 1;

  if (isPowerOfTwo(base))
    {
      unsigned int n, tmp;
      mtx_int_t mask;

      tmp = bits_per_digit - 1;
      mask = 1;

      while (tmp > 0)
        {
          mask = (mask << 1) | 1;
          tmp--;
        }

      do
        {
          n = integer.words[0] & mask;
          s[i] = intToChar(n);
          i--;
          integer >>= bits_per_digit;
        }
      while (integer > 0);
    }

  else
    {
      MtxInteger n, tmp;

      do
        {
          integer.unsigned_divide(base, tmp, n);
          s[i] = intToChar((unsigned int)n.words[0]);
          i--;
          integer = tmp;
        }
      while (integer > 0);
    }

  if (negative)
    {
      s[i] = '-';
      str = s + i;
    }

  else
    str = s + i + 1;

  delete[] s;
  return str;
}

MTXVALUE_CLONE_METHOD(MtxInteger);

#define ADD_SUBTRACT_MULTIPLY_DIVIDE(name, op) \
MtxValue * \
MtxInteger::name(const MtxValue *value) const \
{ \
  const MtxInteger *integer; \
  integer = dynamic_cast<const MtxInteger *>(value); \
  if (integer == NULL) \
    { \
      const MtxRational this_rational(*this), *rational; \
      rational = dynamic_cast<const MtxRational *>(value); \
      if (rational == NULL) \
        { \
          const MtxMatrix *matrix; \
          matrix = dynamic_cast<const MtxMatrix *>(value); \
          if (matrix == NULL) \
            throw MtxException("internal error: operation on unknown type"); \
          return (this_rational op *matrix).clone(); \
        } \
      return (this_rational op *rational).clone(); \
    } \
  else \
    return (*this op *integer).clone(); \
}

ADD_SUBTRACT_MULTIPLY_DIVIDE(add, +)
ADD_SUBTRACT_MULTIPLY_DIVIDE(subtract, -)
ADD_SUBTRACT_MULTIPLY_DIVIDE(multiply, *)
ADD_SUBTRACT_MULTIPLY_DIVIDE(divide, /)

MtxValue *
MtxInteger::opposite() const
{
  return (-*this).clone();
}

void
MtxInteger::setBase(unsigned int base)
{
  if (base < 2 || base > MTX_MAX_BASE)
    return;

  MtxInteger::base = base;
}

unsigned int
MtxInteger::getBase()
{
  return base;
}

MtxInteger
MtxInteger::gcd(MtxInteger a, MtxInteger b)
{
  MtxInteger tmp;

  a = a.abs();
  b = b.abs();

  while (b != 0)
    {
      tmp = b;
      b = a % b;
      a = tmp;
    }

  return a;
}