约定:在博主的校内 OJ 上, long double 类型支持阶码 $15$ 位、有效值 $63$ 位(即科学记数法保留 $64$ 位有效数字)。以下讨论以此为前提。
设非负整数 $a, b, p$ 满足 $a, b<p<2^{63},\ p \ne 0,$ 求 $ab \bmod p.$
$$ab \bmod p=ab-p\left\lfloor {ab \over p} \right\rfloor$$
我们知道 $\left\lfloor {ab \over p} \right\rfloor \leqslant {ab \over p} < p,$ 因此 $\left\lfloor {ab \over p} \right\rfloor$ 可以使用 $64$ 位无符号整数 unsigned long long 存储。
在计算过程中,我们可能会产生上溢出,且 unsigned long long 上溢出的结果是对 $2^{64}$ 取模,这会导致我们无法还原 $\left\lfloor {ab \over p} \right\rfloor.$
但是,注意到浮点数的精度溢出会在二进制下舍入,规则与“四舍六入五成双”类似,如果多余的尾数最高位为 $0,$ 就舍去;如果多余的尾数有且仅有最高位为 $1,$ 就向使结果末尾为 $0$ 的方向舍入;否则就进位。
我们考虑使用 long double 来计算 ${ab \over p},$ 然后转化为 unsigned long long,通过向 $0$ 取整规则得到 $\left\lfloor {ab \over p} \right\rfloor.$
通过分析,我们得到用 long double 计算 $\left\lfloor {ab \over p} \right\rfloor$ 的误差不超过 $\pm1.$
最后我们在模 $2^{64}$ 意义下计算 $ab-p\left\lfloor {ab \over p} \right\rfloor$ 得到 $ab \bmod p.$
于是代码就是这么写的:
1 typedef unsigned long long QWORD; 2 QWORD multi(QWORD x, QWORD y, QWORD mod) 3 { 4 return (x*y-(QWORD)((long double)x*y/mod)*mod+mod)%mod; 5 }
我们来理论验证一下以上细节:
- 由于 $2^{-64}<{1 \over p} \leqslant {ab \over p}<p<2^{64}$ 并且 $2^e \leqslant {ab \over p} < 2^{e+1},$ 阶码 $e$ 的范围只有 $[-64, 64)$,完全足够使用。
- 引理 $1.$ 如果 $1 \leqslant p, q<2^{64},$ 那么计算 $\left\lfloor q \over p \right\rfloor$ 没有误差。
- 产生误差的唯一原因是 $q \over p$ 能够表示为 $2^e(1.\underset{k}{\underbrace{111\cdots1}}0\cdots)_{(2)}\ (k \geqslant 64),$ 即 $2^ep\left( 2^{k+2}-2 \right) \leqslant 2^{k+1}q < 2^ep\left( 2^{k+2}-1 \right).$
- 当 $e \geqslant 0$ 时,因为 $2^{e+1}p, q$ 是整数,要想 $p\left( 2^{e+1}-2^{e-k} \right) \leqslant q < p\left( 2^{e+1}-2^{e-k-1} \right),$ 必然 $2^{e-k}p \geqslant 1.$ $q \geqslant p\left( 2^{e+1}-2^{e-k} \right)>2^ep \geqslant 2^k \geqslant 2^{64},$ 矛盾。
- 当 $e<0$ 时,因为 $4p, 2^{-e}q$ 是整数,要想 $p\left(4-2^{-k}\right) \leqslant 2^{-e}q < p\left(4-2^{-1-k}\right),$ 必然 $2^{-k}p \geqslant 1.\ p \geqslant 2^k \geqslant 2^{64},$ 矛盾。
- 综上,整数下取整除法没有误差。
- 当 $ab<2^{64}$ 时,计算 $ab$ 没有误差,根据上述引理,计算 $\left\lfloor {ab \over p} \right\rfloor$ 也没有误差。
- 当 $ab \geqslant 2^{64}$ 时,计算 $ab$ 的误差可以这么估计:$2^ex-2^{e-1} \leqslant ab \leqslant 2^ex+2^{e-1},$ 其中 $2^ex$ 是舍入结果,$2^{63} \leqslant x<2^{64}.$ 所以 $\left(2^{63}-1\right)p>p^2>ab \geqslant \left(2x-1\right)2^{e-1} \geqslant \left(2^{64}-1\right)2^{e-1}>\left(2^{63}-1\right)2^e,$ 因此 $\left\lvert 2^ex-ab \right\rvert < p.$ 根据引理 $1,$ 由于 $2^ex$ 和 $x$ 只有阶码上的差别,$\large \left\lfloor {2^ex \over p} \right\rfloor$ 同样没有误差。由此,所有误差产生于 $\large\frac{ab}p$ 和 $\large\frac{2^ex}p$ 之间的差异,这不超过 $\pm 1.$
- 博主正在思考还有哪些情况没有误差,如果大家知道烦请告诉我;也可以给出导致误差的反例。
- 博主的数学很菜,恳求大家检查并指出上述证明中的错误。
原文地址:https://www.cnblogs.com/nealchen/p/10346414.html