大整数&C

Question

我正在编写一个程序，生成大整数，将它们保存在数组中，并执行一些基本操作，如乘或加。

我真的很担心实际代码的性能，并希望得到一些提示或改进，以使其更快。任何建议都是受欢迎的，即使它改变了我的整个程序或数据类型。

我将在下面添加一些代码，以便您能够看到我正在使用的结构以及我如何处理这个B.I.N。

unsigned int seed;


void initCharArray( char *L, unsigned N )
{
  for ( int i=0; i< N; i++ ) 
  {
    L[i] = i%50;
  }
}


char Addition( char *Vin1, char *Vin2, char *Vout, unsigned N )
{
  char CARRY = 0;
  for ( int i=0; i< N; i++ ) 
  {
    char R = Vin1[i] + Vin2[i] + CARRY;
    if ( R <= 9 )
    {
      Vout[i] = R; CARRY = 0;
    }
    else
    {
      Vout[i] = R-10; CARRY = 1;
    }
  }
  return CARRY;
}



int main(int argc, char **argv)
{
int N=10000;
unsigned char *V1=new int[N];
unsigned char *V2=new int[N];
unsigned char *V3=new int[N];

initCharArray(V1,N); initCharArray(V2,N);

Addition(V1,V2,V3,N);
}

Answer 1

既然现代拥有人在处理固定位长数字时效率很高，你为什么没有一个数组呢？

假设您使用unsigned long long。它们应该是64位宽，所以最大可能的unsigned long long应该是2^64 - 1。让我们将任何数字表示为数字集合，如下所示：

-big_num = ( n_s, n_0, n_1, ...)

-n_s只需0和1表示+和-符号

-n_0代表0到10^a -1之间的数字(指数a)

-n_1表示10^a到10^(a+1) -1之间的数字，以此类推.

确定a：

所有n_都必须以10^a-1为界。因此，当添加两个big_num时，这意味着我们需要添加n_，如下所示：

// A + B = ( (to be determent later),
//          bound(n_A_1 + n_B_1) and carry to next,
//          bound(n_A_2 + n_B_2 + carry) and carry to next,
//        ...)

可以按以下方式进行边界：

bound(n_A_i + n_B_i + carry) = (n_A_i + n_B_i + carry)%(10^a)

因此，i+1的进位被确定为：

// carry (to be used in i+1) = (n_A_i + n_B_i + carry)/(10^a) 
// (division of unsigned in c++ will floor the result by construction)

这告诉我们，最坏的情况是carry = 10^a -1，因此最糟糕的加法(n_A_i + n_B_i + carry)是：(worst case) (10^a-1) + (10^a-1) + (10^a-1) = 3*(10^a-1)，因为类型是无符号的long，如果我们不想在这个加法上有溢出，我们必须对指数a进行界，这样：

//    3*(10^a-1) <= 2^64 - 1, and a an positive integer
// => a <= floor( Log10((2^64 - 1)/3 + 1) )
// => a <= 18

因此，现在修正的是最大可能的a=18，因此用unsigned long long表示的最大可能的n_是10^18 -1 = 999,999,999,999,999,999。有了这个基本的设置，现在让我们来了解一些实际的代码。现在，我将使用std::vector来保存我们讨论过的big_num，但这可能会改变：

// Example code with unsigned long long
#include <cstdlib>
#include <vector>
//
// FOR NOW BigNum WILL BE REPRESENTED
// BY std::vector. YOU CAN CHANGE THIS LATTER
// DEPENDING ON WHAT OPTIMIZATIONS YOU WANT
//
using BigNum = std::vector<unsigned long long>;

// suffix ULL garanties number be interpeted as unsigned long long
#define MAX_BASE_10 999999999999999999ULL

// random generate big number
void randomize_BigNum(BigNum &a){
  // assuming MyRandom() returns a random number
  // of type unsigned long long
  for(size_t i=1; i<a.size(); i++)
    a[i] = MyRandom()%(MAX_NUM_BASE_10+1); // cap the numbers
}

// wrapper functions
void add(const BigNum &a, const BigNum &b, BigNum &c); // c = a + b
void add(const BigNum &a, BigNum &b);         // b = a + b

// actual work done here
void add_equal_size(const BigNum &a, const BigNum &b, BigNum &c, size_t &N);
void add_equal_size(const BigNum &a, const BigNum &b, size_t &N);
void blindly_add_one(BigNum &c);
// Missing cases
// void add_equal_size(BigNum &a, BigNum &b, BigNum &c, size_t &Na, size_t &Nb);
// void add_equal_size(BigNum &a, BigNum &b, size_t &Na, size_t &Nb);

int main(){
  size_t n=10;
  BigNum a(n), b(n), c(n);
  randomize_BigNum(a);
  randomize_BigNum(b);
  add(a,b,c);
  return;
}

包装器函数应该如下所示。它们将安全地防止数组调用的不正确大小：

// To do: add support for when size of a,b,c not equal

// c = a + b
void add(const BigNum &a, const BigNum &b, BigNum &c){

  c.resize(std::max(a.size(),b.size()));

  if(a.size()==b.size())
    add_equal_size(a,b,c,a.size());
  else
    // To do: add_unequal_size(a,b,c,a.size(),b.size());

  return;
};
// b = a + b
void add(const BigNum &a, const BigNum &b){

  if(a.size()==b.size())
    add_equal_size(a,b,a.size());
  else{
    b.resize(a.size());
    // To do: add_unequal_size(a,b,a.size());
  }

  return;
};

这项工作的主要内容将在这里完成(如果您知道自己在做什么，可以直接调用它并跳过函数调用)：

// If a,b,c same size array
// c = a + b
void add_equal_size(const BigNum &a, const BigNum &b, BigNum &c, const size_t &N){
  // start with sign of c is sign of a
  // Specific details follow on whether I need to flip the
  // sign or not
  c[0] = a[0];
  unsigned long long carry=0;

  // DISTINGUISH TWO GRAND CASES:
  //
  // a and b have the same sign
  // a and b have oposite sign
  // no need to check which has which sign (details follow)
  //
  if(a[0]==b[0]){// if a and b have the same sign
    //
    // This means that either +a+b or -a-b=-(a+b)
    // In both cases I just need to add two numbers a and b
    // and I already got the sign of the result c correct form the
    // start
    //
    for(size_t i=1; i<N;i++){
      c[i]  = (a[i] + b[i] + carry)%(MAX_BASE_10+1);
      carry = c[i]/(MAX_BASE_10+1);      
    }
    if(carry){// if carry>0 then I need to extend my array to fit the final carry
      c.resize(N+1);
      c[N]=carry;
    }
  }
  else{// if a and b have opposite sign
    //
    // If I have opposite sign then I am subtracting the
    // numbers. The following is inspired by how
    // you can subtract two numbers with bitwise operations.
    for(size_t i=1; i<N;i++){
      c[i]  = (a[i] + (MAX_BASE_10 - b[i]) + carry)%(MAX_BASE_10+1);
      carry = c[i]/(MAX_BASE_10+1);      
    }
    if(carry){ // I carried? then I got the sign right from the start
      // just add 1 and I am done
      blindly_add_one(c);
    }
    else{ // I didn't carry? then I got the sign wrong from the start
      // flip the sign
      c[0] ^= 1ULL;
      // and take the compliment
      for(size_t i=1; i;<N;i++)
    c[i] = MAX_BASE_10 - c[i];
    }
  }
  return;
};

关于// if a and b have opposite sign案例的一些细节如下:让我们在基数10中工作。假设我们正在减去a - b，让我们将其转换为一个加法。定义以下操作：

让我们命名一个数字di的10位基数。那么任何数字都是n = d1 + 10*d2 + 10*10*d3...，一个数字的补充现在将被定义为：

     `compliment(d1) = 9-d1`

那么，对一个数字n的称赞是：

   compliment(n) =         compliment(d1)
                   +    10*compliment(d2)
                   + 10*10*compliment(d3)
                   ...

考虑两种情况，a>b和a<b

a>b的例子:唯恐a=830和b=126。执行下面的830 - 126 -> 830 + compliment(126) = 830 + 873 = 1703确定，如果是a>b，我会删除1，然后添加1，结果是704！

a<b的例子:唯恐a=126和b=830。做下面的126 - 830 -> 126 + compliment(830) = 126 + 169 = 295 .？如果我赞美它呢？compliment(295) = 704！所以如果a<b我已经有结果了..。有相反的标志。

在我们的例子中，由于数组中的每个数字都有MAX_BASE_10的限制，所以对我们的数字的补充是

compliment(n) = MAX_BASE_10 - n

因此，使用这种恭维来将减法转换为加法，我只需要注意如果我在加法结束时多带了1( a>b的情况)。现在的算法是

iteration):

• na_i - nb_i + carry(i-1)

• convert -> na_i +na_i(Nb_i)+ carry(i-1)

• bound结果-> (na_i + -> (Nb_i)+ carry(i-1))%MAX_BASE_10

• find )进位-> (na_i +->(Nb_i)+进位(I-1))))/MAX_BASE_10

• keep关于将数组numbers...

• At添加到数组的末尾(如果我携带忘记进位，加1。否则接受结果

的恭维。

这个“并添加一个”是由另一个函数完成的：

// Just add 1, no matter the sign of c
void blindly_add_one(BigNum &c){
  unsigned long long carry=1;
  for(size_t i=1; i<N;i++){
    c[i]  = carry%(MAX_BASE_10+1);
    carry = c[i]/(MAX_BASE_10+1);
  }
  if(carry){ // if carry>0 then I need to extend my basis to fit the number
    c.resize(N+1);
    c[N]=carry;
  }
};

好到这为止。具体来说，在这段代码中，不要忘记在函数开始时，我们将c的符号设置为a符号。因此，如果我在最后进行，这意味着我有|a|>|b|，我做了+a-b>0或-a+b=-(a-b)<0。无论是哪种情况，将结果c符号设置为a符号都是正确的。如果我不带，我就用|a|<|b|和+a-b<0或者-a+b=-(a-b)>0。在这两种情况下，将结果c符号设置为a符号是不正确的，所以如果我不携带，我需要翻转这个标志。

以下函数的操作方式与上面的相同，而不是c = a + b，而是b = a + b

// same logic as above, only b = a + b
void add_equal_size(BigNum &a, BigNum &b, size_t &N){

  unsigned long long carry=0;
  if(a[0]==b[0]){// if a and b have the same sign
    for(size_t i=1; i<N;i++){
      b[i]  = (a[i] + b[i] + carry)%(MAX_BASE_10+1);
      carry = b[i]/(MAX_BASE_10+1);      
    }
    if(carry){// if carry>0 then I need to extend my basis to fit the number
      b.resize(N+1);
      b[N]=carry;
    }
  }
  else{ // if a and b have oposite sign
    b[0] = a[0];
    for(size_t i=1; i<N;i++){
      b[i]  = (a[i] + (MAX_BASE_10 - b[i]) + carry)%(MAX_BASE_10+1);
      carry = b[i]/(MAX_BASE_10+1);      
    }
    if(carry){
      add_one(b);
    }
    else{
      b[0] ^= 1ULL;
      for(size_t i=1; i;<N;i++)
    b[i] = MAX_BASE_10 - b[i];
    }
  }
  return;
};

这是关于如何在数组中使用无符号数字来表示非常大的整数的基本设置。

从这里到哪里

从现在开始，为了优化代码，他们有很多事情要做，我将提到我可以想到的几点：

-Try和用可能的BLAS调用替换数组

-Make确定您正在利用矢量化的优势。取决于您如何编写循环，您可能生成或不生成矢量化代码。如果您的数组变得很大，您可能会从中受益。

-In的精神，上面的精神，确保您有适当的对齐数组在内存中，以实际利用矢量化。根据我的理解，std::vector不能保证对齐。这两种药物都不是盲malloc。我认为boost库有一个矢量版本，您可以声明一个固定的对齐方式，在这种情况下，您可以为unsigned long long数组请求一个64位对齐数组。另一种选择是让您自己的类管理原始指针和剂量对齐分配，并使用自定义分配器。借用aligned_malloc和aligned_free中的https://embeddedartistry.com/blog/2017/02/22/generating-aligned-memory/，您可以使用这样的类来替换std：：https://embeddedartistry.com/blog/2017/02/22/generating-aligned-memory/：

// aligned_malloc and aligned_free from:
// https://embeddedartistry.com/blog/2017/02/22/generating-aligned-memory/

// wrapping in absolutly minimal class to handle
// memory allocation and freeing 
class BigNum{
private:
  unsigned long long *ptr;
  size_t size;
public:  
  BigNum() : ptr(nullptr)
       , size(0)
  {};  

  BigNum(const size_t &N) : ptr(nullptr)
              , size(N)
  {
    resize(N);
  }  
  // Defining destructor this will now delete copy and move constructor and assignment. Make your own if you need them
  ~BigNum(){
    aligned_free(ptr);
  }  

  // Access an object in aligned storage
  const unsigned long long& operator[](std::size_t pos) const{
    return *reinterpret_cast<const unsigned long long*>(&ptr[pos]);
  }
  // return my size
  void size(){    
    return size;
  }
  // resize memory allocation
  void resize(const size_t &N){
    size = N;
    if(N){
      void* temp = aligned_malloc(ptr,N+1); // N+1, always keep first entry for the sign of BigNum 
      if(temp!=nullptr)
    ptr = static_cast<unsigned long long>(temp);
      else
    throw std::bad_alloc();
    }
    else{
      aligned_free(ptr);
    }
  }
};