C/C++高精度(加减乘除)算法的实现
作者:CodeOfCC
前言
C/C++基本类型做算术运算时长度有限,但也基本满足大部分场景的使用,有时需要计算大长度数据就有点无能为力了,比如1000的阶乘结果就有2000多位,用基本类型是无法计算的。
高精度的算法,一般的方式是用一个很长的数组去记录数据,数组的每一位记录固定位数的数字,记录顺序是低位到高位。计算方式则通常采用模拟立竖式计算。计算方式有一些优化方法,比如FFT快速傅里叶变换优化乘法,牛顿迭代优化除法。
本方法采用的,是上述的一般方式,用数组从低到高位记录数据,计算方式采用模拟立竖式计算。
一、必要的参数
1、计算位数
需要设定数组每个元素记录的位数即计算的进制(十进制、百进制、千进制)。计算位数越大,计算性能越好,在32/64位程序中,位数的范围可以是[1,9]。
2、计算进制
与计算位数直接对应,比如1位等于10进制,2位等于百进制。
3、数组类型
根据计算位数定义能装载数据的整型,尽量不浪费空间又不溢出,比如2位以内可以用int8,4位可以用int16等。
代码实现如下:
#include<stdint.h> #include<stdio.h> //计算位数,范围[1,9] #define NUM_DIGIT 9 //数组类型 #if NUM_DIGIT<=0 static_assert(1, "NUM_DIGIT must > 0 and <=9"); #elif NUM_DIGIT <=2 #define _NUM_T int8_t #elif NUM_DIGIT <=4 #define _NUM_T int16_t #elif NUM_DIGIT <=9 //NUM_DIGIT在[5,9]时用int32就可以装载数据,但是实测发现用int32在32位程序中,对除法性能有影响,用int64反而正常,这是比较奇怪的现象。 #define _NUM_T int64_t #else static_assert(1, "NUM_DIGIT must > 0 and <=9"); #endif //计算进制 static const size_t _hex = NUM_DIGIT == 1 ? 10 : NUM_DIGIT == 2 ? 100 : NUM_DIGIT == 3 ? 1000 : NUM_DIGIT == 4 ? 10000 : NUM_DIGIT == 5 ? 100000 : NUM_DIGIT == 6 ? 1000000 : NUM_DIGIT == 7 ? 10000000 : NUM_DIGIT == 8 ? 100000000 : NUM_DIGIT == 9 ? 1000000000 : -1; #define _MIN_NUM_SIZE 30/NUM_DIGIT
上述代码只要设置NUM_DIGIT(计算位数)就可以,其他值(计算进制、数组类型)都会自动生成。
二、辅助函数
1、初始化函数
由整型、字符串初始化高精度数组的方法。
2、比较函数
比较两个高精度数的大小等于。
3、输出函数
将高精度数组打印输出。
代码实现如下:
/// <summary> /// 通过字符串初始化 /// </summary> /// <param name="a">[in]高精度数组</param> /// <param name="value">[in]字符串首地址</param> /// <param name="len">[in]字符串长度</param> /// <returns>数据长度</returns> static size_t initByStr(_NUM_T* a, const char* value, size_t len) { size_t shift = 10; size_t i, j, k = 0, n; n = len / NUM_DIGIT; for (i = 0; i < n; i++) { a[i] = (value[len - k - 1] - '0'); for (j = 1; j < NUM_DIGIT; j++) { a[i] += (value[len - k - j - 1] - '0') * shift; shift *= 10; } shift = 10; k += NUM_DIGIT; } if (n = len - n * NUM_DIGIT) { a[i] = (value[len - k - 1] - '0'); for (j = 1; j < n; j++) { a[i] += (value[len - k - j - 1] - '0') * shift; shift *= 10; } i++; } return i; } /// <summary> /// 通过无符号32整型初始化 /// </summary> /// <param name="a">[in]高精度数组</param> /// <param name="value">[in]整型值</param> /// <returns>数据长度</returns> static size_t initByInt32(_NUM_T* a, uint32_t value) { size_t i; for (i = 0; i < 8096; i++) { a[i] = value % _hex; value /= _hex; if (!value) { i++; break; } } return i; } /// <summary> /// 通过无符号64整型初始化 /// </summary> /// <param name="a">[in]高精度数组</param> /// <param name="value">[in]整型值</param> /// <returns>数据长度</returns> static size_t initByInt64(_NUM_T* a, uint64_t value) { size_t i; for (i = 0; i < 8096; i++) { a[i] = value % _hex; value /= _hex; if (!value) { i++; break; } } return i; } /// <summary> /// 比较两个高精度数的大小 /// </summary> /// <param name="a">[in]第一个数</param> /// <param name="aLen">[in]第一个数的长度</param> /// <param name="b">[in]第二个数</param> /// <param name="bLen">[in]第二个数的长度</param> /// <returns>1是a>b,0是a==b,-1是a<b</returns> static int compare(const _NUM_T* a, size_t aLen, const _NUM_T* b, size_t bLen) { if (aLen > bLen) { return 1; } else if (aLen < bLen) { return -1; } else { for (int i = aLen - 1; i > -1; i--) { if (a[i] > b[i]) { return 1; } else if (a[i] < b[i]) { return -1; } } } return 0; } /// <summary> /// 打印输出结果 /// </summary> static void print(const _NUM_T* a, size_t aLen) { int i = aLen - 1, j = 0, n = _hex; char format[32]; sprintf(format, "%%0%dd", NUM_DIGIT); printf("%d", a[i--]); for (; i > -1; i--) printf(format, a[i]); }
三、实现加减乘除
1、加法
由于数组存储数据的顺序是由低位到高位的,所以只需要两个数从数组0位开始相加,结果除以计算位数,余数保存在第一个数数组当前位,商累加到第一个数数组的下一位,然后下标加1重复前面的计算,直到全部元素计算完。
代码实现如下:
/// <summary> /// 加法(累加) ///结果会保存在a中 /// </summary> /// <param name="a">[in]被加数</param> /// <param name="aLen">[in]被加数长度</param> /// <param name="b">[in]加数</param> /// <param name="bLen">[in]加数长度</param> /// <returns>结果长度</returns> static size_t acc(_NUM_T* a, size_t aLen, const _NUM_T* b, size_t bLen) { size_t i, n = bLen; const size_t max = _hex - 1; if (aLen < bLen) { memset(a + aLen, 0, (bLen - aLen + 1) * sizeof(_NUM_T)); } else { a[aLen] = 0; } for (i = 0; i < bLen; i++) { uint64_t temp = (uint64_t)a[i] + (uint64_t)b[i]; a[i] = temp % _hex; a[i + 1] += temp / _hex; } for (; i < aLen; i++) { uint64_t temp = a[i]; if (temp <= max) break; a[i] = temp % _hex; a[i + 1] += temp / _hex; } n = aLen > bLen ? aLen : bLen; if (a[n]) return n + 1; return n; }
2、减法
减法的原理和加法是类似的,两个数从数组0位开始相减,结果大于等于0直接保存到第一个数数组当前位,否则结果小于零结果+计算进制保存到第一个数数组当前位,同时下一位减等于1,然后下标加1重复前面的计算,直到全部元素计算完。
代码实现如下:
/// <summary> /// 减法(累减) ///结果会保存在a中 /// </summary> /// <param name="a">[in]被减数,被减数必须大于等于减数</param> /// <param name="aLen">[in]被减数长度</param> /// <param name="b">[in]减数</param> /// <param name="bLen">[in]减数长度</param> /// <returns>结果长度</returns> static size_t subc(_NUM_T* a, size_t aLen, const _NUM_T* b, size_t bLen) { size_t m, n, i; if (aLen < bLen) { m = bLen; n = aLen; } else { n = bLen; m = aLen; } for (i = 0; i < n; i++) { int64_t temp = (int64_t)a[i] - (int64_t)b[i]; a[i] = temp; if (temp < 0) { a[i + 1] -= 1; a[i] += _hex; } } for (; i < m; i++) { _NUM_T temp = a[i]; if (temp < 0) { a[i + 1] -= 1; a[i] += _hex; } else { break; } } for (size_t i = aLen - 1; i != SIZE_MAX; i--) if (a[i]) return i + 1; return 1; }
3、乘法
依然是参考立竖式计算方法,两层循环,遍历一个数的元素去乘以另一个数的每个元素,结果记录到新的数组中。
代码实现如下:
/// <summary> /// 乘法 /// </summary> /// <param name="a">[in]被乘数</param> /// <param name="aLen">[in]被乘数长度</param> /// <param name="b">[in]乘数</param> /// <param name="bLen">[in]乘数长度</param> /// <param name="c">[out]结果,数组长度必须大于等于aLen+bLen+1,且全部元素为0</param> /// <returns>结果长度</returns> static size_t multi(const _NUM_T* a, size_t aLen, const _NUM_T* b, size_t bLen, _NUM_T c[]) { _NUM_T d; size_t j; c[aLen + bLen] = 0; for (size_t i = 0; i < aLen; i++) { d = 0; for (j = 0; j < bLen; j++) { uint64_t temp = (uint64_t)a[i] * (uint64_t)b[j] + c[j + i] + d; d = temp / _hex; c[j + i] = temp % _hex; } if (d) { c[j + i] = d; } } for (size_t k = aLen + bLen; k != SIZE_MAX; k--) if (c[k]) return k + 1; return 1; }
4、除法
基本的除法直接利用减法就可以实现,本方法使用的是模拟立竖式计算的方式,将移动除数使其与被除数高位对其,乘以一定倍数让除数与被除数值接近再进行减法运算。
代码实现如下:
/// <summary> /// 除法 /// </summary> /// <param name="a">[in]被除数,被除数必须大于除数。小于商为0、余数为除数,等于商为1、余数为0,不需要在函数内实现,在外部通过compare就可以做到。</param> /// <param name="aLen">[in]被除数长度</param> /// <param name="b">[in]除数</param> /// <param name="bLen">[in]除数长度</param> /// <param name="c">[out]商,数组长度大于等于aLen,且全部元素为0</param> /// <param name="mod">[out]余数,数组长度大于等于aLen</param> /// <param name="modLen">[out]余数长度</param> /// <param name="temp">[in]临时缓冲区,由外部提供以提高性能,数组长度大于等于aLen+bLen+1</param> /// <returns>商长度</returns> static size_t divide(const _NUM_T* a, size_t aLen, const _NUM_T* b, size_t bLen, _NUM_T* c, _NUM_T* mod, size_t* modLen, _NUM_T* temp) { intptr_t n, l; int equal; memcpy(mod, a, (aLen) * sizeof(_NUM_T)); n = aLen - bLen; l = aLen; while (n > -1) { equal = compare(mod + n, l - n, b, bLen); if (equal == -1) { n--; if (!mod[l - 1]) l--; continue; } if (equal == 0) { c[n]++; n -= bLen; l -= bLen; continue; } uint64_t x; if ((l - n) > bLen) { x = ((mod[l - 1]) * _hex) / ((uint64_t)b[bLen - 1] + 1); if (x == 0) { x = (mod[l - 2]) / ((uint64_t)b[bLen - 1] + 1); } } else { x = (mod[l - 1]) / ((uint64_t)b[bLen - 1] + 1); } if (x == 0) x = 1; c[n] += x; if (x == 1) { l = n + subc(mod + n, l - n, b, bLen); } else { _NUM_T num[_MIN_NUM_SIZE]; size_t len = initByInt32(num, x); memset(temp, 0, (aLen + bLen + 1) * sizeof(_NUM_T)); len = multi(b, bLen, num, len, temp); l = n + subc(mod + n, l - n, temp, len); } } intptr_t i = l - 1; for (; i > -1; i--) if (mod[i]) break; if (i == -1) { mod[0] = 0; *modLen = 1; } else { *modLen = i + 1; } if (c[aLen - bLen] != 0) return aLen - bLen + 1; else return aLen - bLen; }
四、使用例子
1、加法例子
求 n 累加和(ja)。用高精度方法,求 s=1+2+3+……+n 的精确值(n 以一般整数输入)。
int main(){ int64_t n; size_t len, len2; _NUM_T num[1024]; _NUM_T num2[1024]; std::cin >> n; len = initByInt32(num, 0); for (int64_t i = 1; i <= n; i++) { len2 = initByInt64(num2, i); len = acc(num, len, num2, len2); } print(num, len); return 0; }
2、减法例子
两个任意十一位数的减法;(小于二十位)
int main() { _NUM_T a1[8096], a2[8096]; std::string s1, s2; std::cin >> s1 >> s2; size_t len1,len2; len1 =initByStr(a1, s1.c_str(), s1.size()); len2 = initByStr(a2, s2.c_str(), s2.size()); len1 =subc(a1, len1, a2, len2); print(a1,len1); return 0; }
3、乘法例子
求 N!的值(ni)。用高精度方法,求 N!的精确值(N 以一般整数输入)。
int main(){ int64_t n; size_t len, len2; _NUM_T num[8096]; _NUM_T num2[8096]; _NUM_T num3[8096]; _NUM_T* p1 = num; _NUM_T* p2 = num3; std::cin >> n; len = initByInt32(num,1); for (int64_t i = 1; i <= n; i++) { len2 = initByInt64(num2, i); memset(p2, 0, 8096* sizeof(_NUM_T)); len = multi(p1, len, num2, len2, p2); _NUM_T* temp = p1; p1 = p2; p2 = temp; } print(p1, len); return 0; }
4、除法例子
给定两个非负整数A,B,请你计算 A / B的商和余数。
int main() { _NUM_T a1[8096], a2[8096],c[8096],mod[8096], temp[8096]; std::string s1, s2; std::cin >> s1 >> s2; size_t len1,len2,ml; len1 =initByStr(a1, s1.c_str(), s1.size()); len2 = initByStr(a2, s2.c_str(), s2.size()); memset(c, 0, 8096 * sizeof(_NUM_T)); len1 =divide(a1, len1, a2, len2,c, mod,&ml, temp); print(c,len1); std::cout<< std::endl ; print(mod, ml); return 0; }
以上就是C/C++高精度(加减乘除)算法的实现的详细内容,更多关于C++高精度算法的资料请关注脚本之家其它相关文章!