Здравствуйте!

Понадобилось тут немного арифметики, решил сам наваять, дело вообщем-то нехитрое

Использовал эту

Автор: Lorenzo_LAMAS
Дата: 12.04.12

тему.

Умножение чисел произвольной длинны взял из книги "Алгоритмические трюки для программистов" (Уоррен). Тут пришлось поковыряться немного, в моем скане куча ошибок в исходниках, да и стиль у автора не очень внятный.

#include <limits>
#include <algorithm>

inline bool unsigned_addition( UINT16 a, UINT16 b, UINT16 *c = 0 )
{
    UINT16 r = a + b; if (c) *c = r; return !(r>=a && r>=b);
}

inline bool unsigned_addition( UINT32 a, UINT32 b, UINT32 *c = 0 )
{
    UINT32 r = a + b; if (c) *c = r; return !(r>=a && r>=b);
}

inline bool unsigned_addition( UINT64 a, UINT64 b, UINT64 *c = 0 )
{
    UINT64 r = a + b; if (c) *c = r; return !(r>=a && r>=b);
}

inline bool unsigned_subtraction( UINT16 a, UINT16 b, UINT16 *c = 0 )
{
    UINT16 r = a - b; if (c) *c = r; return a<b;
}

inline bool unsigned_subtraction( UINT32 a, UINT32 b, UINT32 *c = 0 )
{
    UINT32 r = a - b; if (c) *c = r; return a<b;
}

inline bool unsigned_subtraction( UINT64 a, UINT64 b, UINT64 *c = 0 )
{
    UINT64 r = a - b; if (c) *c = r; return a<b;
}

inline bool unsigned_addition( UINT16 a, UINT16 b, UINT16 &c )
    { return unsigned_addition( a, b, &c ); }
inline bool unsigned_addition( UINT32 a, UINT32 b, UINT32 &c ) 
    { return unsigned_addition( a, b, &c ); }
inline bool unsigned_addition( UINT64 a, UINT64 b, UINT64 &c )
    { return unsigned_addition( a, b, &c ); }

inline bool unsigned_subtraction( UINT16 a, UINT16 b, UINT16 &c )
    { return unsigned_subtraction( a, b, &c ); }
inline bool unsigned_subtraction( UINT32 a, UINT32 b, UINT32 &c )
    { return unsigned_subtraction( a, b, &c ); }
inline bool unsigned_subtraction( UINT64 a, UINT64 b, UINT64 &c )
    { return unsigned_subtraction( a, b, &c ); }




// returns true if carry detected
inline bool signed_addition( INT16 a, INT16 b, INT16 *c = 0 )
{
    INT16 r = a + b; if (c) *c = r; return !(b >= 0 == r >= a);
}

inline bool signed_addition( INT32 a, INT32 b, INT32 *c = 0 )
{
    INT32 r = a + b; if (c) *c = r; return !(b >= 0 == r >= a);
}

inline bool signed_addition( INT64 a, INT64 b, INT64 *c = 0 )
{
    INT64 r = a + b; if (c) *c = r; return !(b >= 0 == r >= a);
}

inline bool signed_subtraction( INT16 a, INT16 b, INT16 *c = 0 )
    { return signed_addition( a, -b, c ); }
inline bool signed_subtraction( INT32 a, INT32 b, INT32 *c = 0 )
    { return signed_addition( a, -b, c ); }
inline bool signed_subtraction( INT64 a, INT64 b, INT64 *c = 0 )
    { return signed_addition( a, -b, c ); }


inline bool signed_addition( INT16 a, INT16 b, INT16 &c ) 
    { return signed_addition( a, b, &c ); }
inline bool signed_addition( INT32 a, INT32 b, INT32 &c )
    { return signed_addition( a, b, &c ); }
inline bool signed_addition( INT64 a, INT64 b, INT64 &c )
    { return signed_addition( a, b, &c ); }

inline bool signed_subtraction( INT16 a, INT16 b, INT16 &c )
    { return signed_subtraction( a, b, &c ); }
inline bool signed_subtraction( INT32 a, INT32 b, INT32 &c )
    { return signed_subtraction( a, b, &c ); }
inline bool signed_subtraction( INT64 a, INT64 b, INT64 &c )
    { return signed_subtraction( a, b, &c ); }


template<typename IntT>
bool integer_addition( IntT a, IntT b, IntT *c = 0)
{
    return signed_addition( (INT64)a, (INT64)b, (INT64*)c );
}

template<typename IntT>
bool integer_addition( IntT a, IntT b, IntT &c )
{
    INT64 tmp;
    bool bRes = signed_addition( (INT64)a, (INT64)b, &tmp );
    c = (IntT)tmp;
    return bRes;
}

template<> bool integer_addition<INT16> (INT16  a, INT16  b, INT16  *c) { return signed_addition( a, b, c ); }
template<> bool integer_addition<INT16> (INT16  a, INT16  b, INT16  &c) { return signed_addition( a, b, c ); }
template<> bool integer_addition<INT32> (INT32  a, INT32  b, INT32  *c) { return signed_addition( a, b, c ); }
template<> bool integer_addition<INT32> (INT32  a, INT32  b, INT32  &c) { return signed_addition( a, b, c ); }
template<> bool integer_addition<INT64> (INT64  a, INT64  b, INT64  *c) { return signed_addition( a, b, c ); }
template<> bool integer_addition<INT64> (INT64  a, INT64  b, INT64  &c) { return signed_addition( a, b, c ); }

template<> bool integer_addition<UINT16>(UINT16 a, UINT16 b, UINT16 *c) { return unsigned_addition( a, b, c ); }
template<> bool integer_addition<UINT16>(UINT16 a, UINT16 b, UINT16 &c) { return unsigned_addition( a, b, c ); }
template<> bool integer_addition<UINT32>(UINT32 a, UINT32 b, UINT32 *c) { return unsigned_addition( a, b, c ); }
template<> bool integer_addition<UINT32>(UINT32 a, UINT32 b, UINT32 &c) { return unsigned_addition( a, b, c ); }
template<> bool integer_addition<UINT64>(UINT64 a, UINT64 b, UINT64 *c) { return unsigned_addition( a, b, c ); }
template<> bool integer_addition<UINT64>(UINT64 a, UINT64 b, UINT64 &c) { return unsigned_addition( a, b, c ); }


template<typename IntT>
bool integer_subtraction( IntT a, IntT b, IntT *c = 0)
{
    return signed_subtraction( (INT64)a, (INT64)b, (INT64*)c );
}

template<typename IntT>
bool integer_subtraction( IntT a, IntT b, IntT &c )
{
    INT64 tmp;
    bool bRes = signed_subtraction( (INT64)a, (INT64)b, &tmp );
    c = (IntT)tmp;
    return bRes;
}

template<> bool integer_subtraction<INT16> (INT16  a, INT16  b, INT16  *c) { return signed_subtraction( a, b, c ); }
template<> bool integer_subtraction<INT16> (INT16  a, INT16  b, INT16  &c) { return signed_subtraction( a, b, c ); }
template<> bool integer_subtraction<INT32> (INT32  a, INT32  b, INT32  *c) { return signed_subtraction( a, b, c ); }
template<> bool integer_subtraction<INT32> (INT32  a, INT32  b, INT32  &c) { return signed_subtraction( a, b, c ); }
template<> bool integer_subtraction<INT64> (INT64  a, INT64  b, INT64  *c) { return signed_subtraction( a, b, c ); }
template<> bool integer_subtraction<INT64> (INT64  a, INT64  b, INT64  &c) { return signed_subtraction( a, b, c ); }

template<> bool integer_subtraction<UINT16>(UINT16 a, UINT16 b, UINT16 *c) { return unsigned_subtraction( a, b, c ); }
template<> bool integer_subtraction<UINT16>(UINT16 a, UINT16 b, UINT16 &c) { return unsigned_subtraction( a, b, c ); }
template<> bool integer_subtraction<UINT32>(UINT32 a, UINT32 b, UINT32 *c) { return unsigned_subtraction( a, b, c ); }
template<> bool integer_subtraction<UINT32>(UINT32 a, UINT32 b, UINT32 &c) { return unsigned_subtraction( a, b, c ); }
template<> bool integer_subtraction<UINT64>(UINT64 a, UINT64 b, UINT64 *c) { return unsigned_subtraction( a, b, c ); }
template<> bool integer_subtraction<UINT64>(UINT64 a, UINT64 b, UINT64 &c) { return unsigned_subtraction( a, b, c ); }



template<typename IntT>
IntT integer_abs( IntT i )  { return (i<0) ? -i : return i; }

template<typename IntT>
IntT integer_sign( IntT i ) { return (i<0) ? -1 : 1; }

template<typename IntT>
int sign_multiple( IntT a, IntT b )
    { return integer_sign(a)==integer_sign(b) ? 1 : -1; }



template<typename IntT>
IntT make_lo_half_mask()
{
    IntT mask = 0xFF;
    IntT resMask = mask;
    for(unsigned i =0; i!=(unsigned)(sizeof(IntT)/2); ++i)
       {
        resMask |= mask;
        mask <<= 8;
       }
    return resMask;
}

template<typename IntT> 
unsigned calc_half_shift() { return (unsigned)(sizeof(IntT)*4); }

template<typename IntT>
IntT lo_part( IntT i ) { return i & make_lo_half_mask<IntT>(); }

template<typename IntT>
IntT hi_part( IntT i ) 
    { return (i >> calc_half_shift<IntT>()) & make_lo_half_mask<IntT>(); }

/* keep sign */
template<typename IntT>
IntT lo_part_s( IntT i ) { return i & make_lo_half_mask<IntT>(); }

template<typename IntT>
IntT hi_part_s( IntT i ) 
    { return i >> calc_half_shift<IntT>(); }




template<typename UIntT, typename UIntT2>
void multiply_num_mn_unsigned_impl( UIntT *w
           , const UIntT *u, const UIntT *v
           , unsigned m, unsigned n
           )
{
    UIntT2 k, t; //, b;
    unsigned i, j; // counters

    for (i=0; i<m; i++)
        w[i] = 0;

    for(j=0; j!=n; j++)
       {
        k = 0;
        for(i=0; i!=m; i++)
           {
            t = u[i]*v[j] + w[i+j] + k;
            w[i+j] = static_cast<UIntT>(lo_part_s(t)); // t; // (Т.е., t & OxFFFF).
            k = hi_part_s(t); // t >> 16;
           }
        w[j+m] = static_cast<UIntT>(k);
       }
}

inline void 
multiply_num_mn_unsigned( UINT8 *w, const UINT8 *u, const UINT8 *v
                        , unsigned m, unsigned n )
    { multiply_num_mn_unsigned_impl<UINT8, UINT16>(w, u, v, m, n ); }

inline void 
multiply_num_mn_unsigned( UINT16 *w, const UINT16 *u, const UINT16 *v
                        , unsigned m, unsigned n )
    { multiply_num_mn_unsigned_impl<UINT16, UINT32>(w, u, v, m, n ); }

inline void 
multiply_num_mn_unsigned( UINT32 *w, const UINT32 *u, const UINT32 *v
                        , unsigned m, unsigned n )
    { multiply_num_mn_unsigned_impl<UINT32, UINT64>(w, u, v, m, n ); }




template<typename UIntT, typename SIntT, typename UIntT2>
void multiply_num_mn_signed_impl( UIntT *w
           , const UIntT *u, const UIntT *v
           , unsigned m, unsigned n
           )
{
    UIntT2 t, b; // k, 
    unsigned i, j; // counters
    multiply_num_mn_unsigned( w, u, v, m, n );

    // Теперь w[] содержит беззнаковое произведение.
    // Корректируем вычитанием v*2**16m при и < 0 и
    // вычитанием u*2**16n при v < О
    if ((SIntT)u[m-1] < 0)
       {
        b = 0; // Инициализируем заем
        for (j=0; j!=n; j++)
            {
             t = w[j+m] - v[j] - b;
             w[j+m] = static_cast<UIntT>(t);
             b = t >> (sizeof(UIntT2)*8 - 1); // 31;
            }
       }

    if ((SIntT)v[n-1] < 0)
       {
        b = 0;
        for(i=0; i!=m; i++)
           {
            t = w[i+n] - u[i] - b;
            w[i+n] = static_cast<UIntT>(t);
            b = t >> (sizeof(UIntT2)*8 - 1); // 31;
           }
       }
    return;
}

inline void 
multiply_num_mn_signed( UINT8 *w, const UINT8 *u, const UINT8 *v
                        , unsigned m, unsigned n )
    { multiply_num_mn_signed_impl<UINT8, INT8, UINT16>(w, u, v, m, n ); }

inline void 
multiply_num_mn_signed( UINT16 *w, const UINT16 *u, const UINT16 *v
                        , unsigned m, unsigned n )
    { multiply_num_mn_signed_impl<UINT16, INT16, UINT32>(w, u, v, m, n ); }

inline void 
multiply_num_mn_signed( UINT32 *w, const UINT32 *u, const UINT32 *v
                        , unsigned m, unsigned n )
    { multiply_num_mn_signed_impl<UINT32, INT32, UINT64>(w, u, v, m, n ); }


template<typename IntT, typename UIntT>
IntT multiply_impl(IntT u, IntT v, UIntT *pLowPart)
{
    UIntT u0 = lo_part_s(u); // u & 0xFF; 
    IntT  u1 = hi_part_s(u); // u >>  8;

    UIntT v0 = lo_part_s(v); // v & 0xFF; 
    IntT  v1 = hi_part_s(v); // v >>  8;

    UIntT w0 = u0*v0;
    IntT  t  = u1*v0 + hi_part_s(w0); // (w0 >> 8) ;
    IntT  w1 = lo_part_s(t); // t & 0xFF;
    IntT  w2 = hi_part_s(t); // t >>  8;
    w1       = u0*v1 + w1;
    if (pLowPart) *pLowPart = (UIntT)(w1 << calc_half_shift<UIntT>()) | (w0 & make_lo_half_mask<UIntT>() );
    return u1*v1 + w2 + hi_part_s(w1); // (w1 >> 8) ;
}



template<typename IntT, typename UIntT>
IntT multiply_high_impl(IntT u, IntT v)
{
    return multiply_impl(u, v, (UIntT*)0);
}

template<typename IntT>
IntT multiply_high( IntT a, IntT b)
{
    return multiply_high_impl<UINT64, UINT64>( a, b );
}

template<> INT16  multiply_high<INT16> (INT16  a, INT16  b) { return multiply_high_impl<INT16,  UINT16>( a, b ); }
template<> INT32  multiply_high<INT32> (INT32  a, INT32  b) { return multiply_high_impl<INT32,  UINT32>( a, b ); }
template<> INT64  multiply_high<INT64> (INT64  a, INT64  b) { return multiply_high_impl<INT64,  UINT64>( a, b ); }

template<> UINT16 multiply_high<UINT16>(UINT16 a, UINT16 b) { return multiply_high_impl<UINT16, UINT16>( a, b ); }
template<> UINT32 multiply_high<UINT32>(UINT32 a, UINT32 b) { return multiply_high_impl<UINT32, UINT32>( a, b ); }
template<> UINT64 multiply_high<UINT64>(UINT64 a, UINT64 b) { return multiply_high_impl<UINT64, UINT64>( a, b ); }



template <typename IntT>
IntT integer_pow( IntT i, unsigned p, bool *pOverflow = 0)
{
    if (pOverflow) *pOverflow = false;
    if (i==0 && p==0)
       {
        if (pOverflow) *pOverflow = true;
        return 0;
       }
    if (!i) return 0;
    if (!p) return 1;

    IntT res = 1;
    do{
       IntT tmp = res*i;
       IntT mh  = multiply_high(res,i);
       //sign_multiple( IntT a, IntT b )

       if (std::numeric_limits<IntT> :: is_signed )
          {
           if (i<0)
              {
               if (sign_multiple(tmp, res)>=0)
                  {
                   if (pOverflow) *pOverflow = true;
                  }
              }
           else
              {
               if (sign_multiple(tmp, res)<0)
                  {
                   if (pOverflow) *pOverflow = true;
                  }
              }

           int mhs  = sign_multiple(res,i);
           if (  (mhs<0 && mh!=(IntT)-1)
              || (mhs>0 && mh!=0)
              )
              {
               if (pOverflow) *pOverflow = true;
              }
          }
       else
          {
           if (mh>0)
              {
               if (pOverflow) *pOverflow = true;
              }
          }
       res = tmp;
       if (!res) { if (pOverflow) *pOverflow = true; return 0; }
       --p;
      } while(p);
    return res;
}

template <typename IntT>
IntT integer_pow( IntT i, unsigned p, bool &bOverflow )
{
    return integer_pow( i, p, &bOverflow );
}