简介
在实际的项目中,同事在移植一个算法时候碰到要进行64位整数的除法运算。找了一下一下,Linux内核中有支持该运算的函数do_div()
,该函数在 Linux/arch/arm/include/asm/div64.h 文件中实现。看不太懂其具体的实现方法,于是我就想能不能自己写一个大数相除的算法。下面就是算法的内容,如有不足之处,敬请指正。
注:在以下公式以及代码中,名字的含义如下:
m : 被除数
n : 除数
q(quotient) : 商
r(remainder) : 余数
以下都使用 100/7 为例进行说明。
最终的代码实现放在文章的最后。
持续减法
假设 m/n 的商为q,余数为r,那么可以得到这个公式:
m = n*q+r,进行变换得到 m-n*q = r (r < n)
于是,我们可以降 q 从0一直增大,增大到满足 r < n 的条件,那么就可以算出q和r的值。
比如100/7。
100 - 7*1 = 93 > 7
100 - 7*2 = 93 - 7 = 86 > 7
100 - 7*3 = 86 - 7 = 79 > 7
100 - 7*4 = 79 - 7 = 72 > 7
……
……
100 - 7*14 = 9 - 7 = 2 < 7
这是一个持续减法减去n的过程,直到满足 r < n 的条件,q就是商,也就是运算的次数。于是我们知道商为14,余数为2。这个算法的实现是 div1()
函数。
这个算法肯定是能够算出两个大数相除的结果,但是如果碰到 m 非常之大,n 又非常小的情况下,比如 m = 0xffffffffffffffff;n = 2,此时就悲剧了,总共要做 0x7fffffffffffffff 次减法。所以说这个算法不一定是高效的。
利用二进制性质持续减法
前面是持续尝试 7*1,7*2, 7*3, … , 7*14做,但这样子效率比较低,我们能不能尝试换一种乘法。即:
7*1,7*2, 7*4, 7*8,如果再乘下去即7*16 = 112 > 100, 所以这个时候就100 - 7*8 = 44接着下一次循环。
可以这里理解:100 = 7 * 8 + 7 * 6 + 2。这个算法的实现是 div2()
函数。
二分法
将除法操作转换为乘法操作。
7 * 100 = 700 > 100
7 * 50 = 350 > 100
7 * 25 = 175 > 100
7 * 12 = 84 < 100,此时 q = 12
7 * (12 + 6) = 126 > 100
7 * (12 + 3) = 105 > 100
7 * (12 + 1) = 91 < 100,此时 q+= 1,即 q = 13
7 * (13 + 1) = 98 < 100,此时 q+= 1,即 q = 14
7 * (14 + 1) = 105 > 100
只要满足当二分到为 1 的时候,且最后一次的结果大于 m,倒数第二次的结果小于 m,到此就要退出循环。
这个算法实际上不能称为大数的除法,因为这里面使用的是乘法操作。当 m 和 n 同时都是很大的数的话,比如 m 和 n 都是64位整数,那么两者相乘的结果肯定超过64位了。所以这个算法只适合较小的数进行相除。之所以列举在这边,是看重其运算的效率高。这个算法的实现是 div3()
函数。
实现代码
/*
三个函数 div1()/div2()/div3() 的三个参数分别为:
m : 被除数
n : 除数
*remain : 余数
函数的返回值为商。
注:测试代码中没有验证参数 n 为 0 值的异常情况。
*/
#include <stdio.h>
typedef unsigned long uint64_t;
typedef long int64_t;
uint64_t div1(uint64_t m, uint64_t n, uint64_t *remain)
{
uint64_t quot = 0;
uint64_t value = 0;
if(m < n)
{
quot = 0;
*remain = m;
}
else if(m == n)
{
quot = 1;
*remain = 0;
}
else
{
value = m-n;
quot = 1;
while(value >= n)
{
quot++;
m = value;
value = m-n;
}
}
*remain = value;
return quot;
}
uint64_t div2(uint64_t m, uint64_t n, uint64_t *remain)
{
uint64_t quot = 0;
uint64_t value = m;
uint64_t multi = 1;
if(m < n)
{
quot = 0;
*remain = m;
}
else if(m == n)
{
quot = 1;
*remain = 0;
}
else
{
while(value >= n)
{
multi = 1;
while((n*multi) <= (value>>1))
{
multi <<= 1; // 持续乘 2
}
quot += multi;
value -= (n*multi);
}
*remain = value;
}
return quot;
}
uint64_t div3(uint64_t m, uint64_t n, uint64_t *remain)
{
uint64_t quot = 0;
uint64_t value = (m>>1);
while(1)
{
if((n*(quot+value)) == m)
{
quot = quot+value;
*remain = 0;
}
else if((n*(quot+value)) > m)
{
if(value == 1)
break;
else
value >>= 1;
}
else
{
quot += value;
if(value != 1)
value >>= 1;
}
}
*remain = m-(n*quot);
return quot;
}
// 测试代码
int main(void)
{
uint64_t m[4] = {0xffffffffffffffff, 261535774000000000, 60253400000, 4LL};
uint64_t n[4] = {0xfffffffffffffffe, 6912470141400, 50698765, 2LL};
int64_t quot = 0; // 结果
uint64_t remainder = 0; // 余数
int flag = 0; // 标记符号
int i = 0;
for(i=1; i<4; i++)
{
// 判断正负符号
if(m[i] < 0)
{
m[i] = -m[i];
flag = 1;
}
if(n[i] < 0)
{
n[i] = -n[i];
flag ^= 1;
}
quot = div1(m[i], n[i], &remainder); // 修改调用的函数名验证其他的函数
if(flag)
quot = -quot;
printf("%d: %lld / %lld = %lld.....%lld\n", i, m[i], n[i], quot, remainder);
quot = 0;
remainder = 0;
}
return 0;
}