洛谷P10446 64位整数乘法 题解 & O(1)快速乘模板

题目链接:P10446 64位整数乘法

题意:

求 \(a\) 乘 \(b\) 对 \(p\) 取模的值。

输入格式:

第一行输入整数 \(a\),第二行输入整数 \(b\),第三行输入整数 \(p\)。

输出格式:

输出一个整数,表示 a*b mod p 的值。

数据范围:

\(1 \le a,b,p \le 10^{18}\) 。

注意到 long long 类型能存储 \([-2^{63}, 2^{63}-1]\) 这一区间内的整数,也就是 \(9\times 10^{18}\) 这个数量级。

如果直接使用 long long 计算,当 \(p\) 很大时,中间 \(a\cdot b\) 最多能达到 \(10^{36}\) 数量级,显然直接就溢出了。

不过我们可以有以下几种方式来解决这个问题,首先最推荐的方案是快速乘,注意这个不是 \(O(1)\) 快速乘。

快速乘和快速幂原理相似,假如 \(b=2^k + 2^3 + 2^1\) ,那么它会将 \(a\cdot b\) 的运算转化为以下形式 \[ a\cdot b = \underset{b\, \text{个}}{\underbrace{a + a + \cdots + a}} = 2^ka + 2^3a + 2^1a \] 时间复杂度 \(O(\log M)\) ,其中 \(M\) 为 \(a,b\) 的值域。

也没啥好说的,直接看代码就知道了。

#include <bits/stdc++.h>
using namespace std;
// #define int long long
#define INF 0x3f3f3f3f3f3f3f3f
typedef long long ll;
template<typename T> void up(T &x, T y) { x < y ? x = y : 0; }
template<typename T> void down(T &x, T y) { x > y ? x = y : 0; }
#define rep(i, a, b) for(int i = (a), i##END = (b); i <= i##END; i++)
#define Rep(i, a, b) for(int i = (a), i##END = (b); i >= i##END; i--)
#define N ((int)())

// 这里 add 是快速取模的一个小技巧
void add(ll &x, ll y, const ll mod) { (x += y) >= mod ? x -= mod : 0; }
ll ksc(ll a, ll b, ll p)
{
    ll ans = 0;
    for(; b > 0; b >>= 1, add(a, a, p))
        if(b & 1) add(ans, a, p);
    return ans;
}
signed main()
{
    ios::sync_with_stdio(0);
    cin.tie(0); cout.tie(0);
    // freopen("check.in","r",stdin);
    // freopen("check.out","w",stdout);
    ll a, b, p; cin >> a >> b >> p;
    cout << ksc(a % p, b % p, p) << '\n'; // 这里别忘取模!!
    return 0;
}

除此以外,还有一个很简单但不一定什么情况都能用的。

那就是大名鼎鼎的 __int128 ,如下

cout << (long long)((__int128_t)a*b%p) << endl;
cout << (long long)((__int128)a*b%p) << endl; // 这样写也是可以的

当然上面的代码直接把 a,b,p 开成 __int128 也可以。

不推荐用 __int128,毕竟这个不一定能用,有的评测机的编译器是不支持的。

接下来就是传说中的 \(O(1)\) 快速乘了,一个非常神奇和玄学的算法,来看下代码。

#include <bits/stdc++.h>
using namespace std;
// #define int long long
#define INF 0x3f3f3f3f3f3f3f3f
typedef long long ll;
template<typename T> void up(T &x, T y) { x < y ? x = y : 0; }
template<typename T> void down(T &x, T y) { x > y ? x = y : 0; }
#define rep(i, a, b) for(int i = (a), i##END = (b); i <= i##END; i++)
#define Rep(i, a, b) for(int i = (a), i##END = (b); i >= i##END; i--)
#define N ((int)())

ll ksc(ll a, ll b, ll p)
{
    ll t = (long double)a / p * b;
    ll res = a * b - t * p;
    if(res < 0) res += p; if(res >= p) res -= p; // 取模
    return res;
}
signed main()
{
    ios::sync_with_stdio(0);
    cin.tie(0); cout.tie(0);
    // freopen("check.in","r",stdin);
    // freopen("check.out","w",stdout);
    ll a, b, p; cin >> a >> b >> p;
    cout << ksc(a, b, p) << '\n';
    return 0;
}

据说是出自 2009年全国信息学奥林匹克冬令营论文 《论程序底层优化的一些方法与技巧》骆可强

原理是余数等于 \(a\cdot b - p \cdot \left\lfloor\frac{ab}{p}\right\rfloor\) ,减号两侧的值即使溢出,他们的差值仍然是固定的。

时间复杂度 \(O(1)\) ,但是这个东西存在误差,据说当 \(p>10^{12}\) 时有可能会算错。这个方法就当了解一下吧。