大整数 & C

Big integer numbers & C

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

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

我将在下面添加一些代码,以便您可以看到我正在使用的结构以及我如何尝试处理这个 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);
}

既然现代拥有者在处理固定位长数字时效率很高,为什么你没有它们的数组?

假设您使用 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 之间的数字 等等等等...

正在确定:

所有 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) 由于类型是 unsigned long 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 让我们将其转换为加法。定义以下操作:

让我们命名一个数字的基数 10 位 di。那么任意一个数字就是n = d1 + 10*d2 + 10*10*d3...一个数字的补码现在定义为:

     `compliment(d1) = 9-d1`

那么一个数n的恭维是:

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

考虑两种情况,a>ba<b:

a>b 的示例:免得说 a=830b=126。执行以下 830 - 126 -> 830 + compliment(126) = 830 + 873 = 1703 好的,如果 a>b,我删除 1,然后加 1 结果是 704!

a<b 的示例:免得说 a=126b=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 情况)。现在的算法是

  • 对于每个数组减法(第 i 次迭代):
  • na_i - nb_i + 进位(i-1)
  • convert -> na_i + compliment(nb_i) + carry(i-1)
  • 绑定结果 -> (na_i + 恭维(nb_i) + 进位(i-1))%MAX_BASE_10
  • 求进位-> (na_i + 恭维(nb_i) + 进位(i-1))/MAX_BASE_10

  • 继续添加数组编号...

  • 在数组的末尾如果我进位,忘了进位 并加 1。否则取结果的恭维

这个 "and add one" 是由另一个函数完成的:

// 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;
};

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

从这里到哪里去

从现在开始他们有很多事情要做来优化代码,我会提到一些我能想到的:

-尝试用可能的 BLAS 调用替换数组的添加

-确保您正在利用矢量化。根据您编写循环的方式,您可能会或可能不会生成矢量化代码。如果您的阵列变大,您可能会从中受益。

-本着上述精神,确保您在内存中正确对齐数组以实际利用矢量化。据我了解 std::vector 不保证对齐。两者都没有盲目 malloc。我认为 boost 库有一个矢量版本,您可以在其中声明固定对齐方式,在这种情况下,您可以为 unsigned long long 数组请求 64 位对齐数组。另一种选择是拥有自己的 class,它使用自定义分配器管理原始指针和剂量对齐分配。从 https://embeddedartistry.com/blog/2017/02/22/generating-aligned-memory/ 借用 aligned_mallocaligned_free 你可以用这样的 class 来替换 std::vector:

// 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);
    }
  }
};