在windows 10 下使用qemu 虛擬化 arm/amd cpu 架構
阿新 • • 發佈：2021-11-22
emm……推導過程什麼的左轉去看各位神犇的部落格吧QwQ
這裡只貼程式碼（以後有時間的話可能會寫寫推導過程）。
首先是喜聞樂見的 NTT 多項式乘法 及 求逆 板子:
namespace NTT{ ll lim, len; inline ll qpow(ll a, ll b){ ll res = 1; while(b){ if(b & 1) res = res * a % mod; a = a * a % mod, b >>= 1; } return res; } inline void get_rev(cl n){ lim = 1, len = 0; while(lim < n) lim <<= 1, ++len; for(int i = 0; i < lim; ++i) p[i] = (p[i >> 1] >> 1) | ((i & 1) << (len - 1)); } inline void ntt(ll A[], cl lim, cl type){ for(int i = 0; i < lim; ++i) if(i < p[i]) swap(A[i], A[p[i]]); for(int mid = 1; mid < lim; mid <<= 1){ ll Wn = qpow(type == 1 ? G : Gi, (mod - 1) / (mid << 1)); for(int i = 0; i < lim; i += (mid << 1)){ ll w = 1; for(int j = 0; j < mid; ++j, w = w * Wn % mod){ ll x = A[i + j], y = w * A[i + j + mid] % mod; A[i + j] = (x + y) % mod; A[i + j + mid] = (x - y + mod) % mod; } } } if(type == 1) return; ll inv = qpow(lim, mod - 2); for(int i = 0; i < lim; ++i) A[i] = A[i] * inv % mod; } inline void Mul(cl n, cl m, ll a[], ll b[]){ get_rev(n + m); ntt(a, lim, 1), ntt(b, lim, 1); for(int i = 0; i < lim; ++i) a[i] = a[i] * b[i] % mod; ntt(a, lim, -1); } inline void Inv(ll n, ll a[], ll b[]){ if(n == 1) return b[0] = qpow(a[0], mod - 2), void(); Inv((n + 1) >> 1, a, b); get_rev(n << 1); for(int i = 0; i < n; ++i) c[i] = a[i]; for(int i = n; i < lim; ++i) c[i] = 0; ntt(c, lim, 1), ntt(b, lim, 1); for(int i = 0; i < lim; ++i) b[i] = (2ll - c[i] * b[i] % mod + mod) * b[i] % mod; ntt(b, lim, -1); for(int i = n; i < lim; ++i) b[i] = 0; } } using namespace NTT;
然後是多項式開根
inline void Sqrt(cl n, ll a[], ll b[]){ if(n == 1) return b[0] = 1, void(); Sqrt((n + 1) >> 1, a, b); get_rev(n << 1); memset(f, 0, sizeof(f)); Inv(n, b, f);//f(x) 是 b(x) 的逆多項式 for(int i = 0; i < n; ++i) d[i] = a[i]; for(int i = n; i < lim; ++i) d[i] = 0; ntt(d, lim, 1), ntt(b, lim, 1), ntt(f, lim, 1); for(int i = 0; i < lim; ++i) b[i] = (b[i] + d[i] * f[i] % mod) * inv2 % mod; ntt(b, lim, -1); for(int i = n; i < lim; ++i) b[i] = 0; }
接著自然是是求 ln
namespace Ln{ inline void Diff(cl n, ll a[], ll b[]){//微分求導 for(int i = 1; i < n; ++i) b[i - 1] = i * a[i]; b[n - 1] = 0; } inline void Integral(cl n, ll a[], ll b[]){//積分 for(int i = 1; i < n; ++i) b[i] = a[i - 1] * qpow(i, mod - 2) % mod; b[0] = 0; } inline void ln(cl n, ll a[], ll b[]){ memset(f1, 0, sizeof(f1)); memset(f2, 0, sizeof(f2)); Diff(n, a, f1), Inv(n, a, f2);//f1(x) = a'(x)，f2(x) = a(x)^(-1) Mul(n, n, f1, f2);//f1(x) = f1(x) * f2(x) Integral(n, f1, b);//g(x) = f(x) 的積分 } } using namespace Ln;
exp 當然也不能少
inline void Exp(cl n, ll a[], ll b[]){
if(n == 1) return b[0] = 1, void();
Exp((n + 1) >> 1, a, b);
get_rev(n << 1);
ln(n, b, d);// d(x) = ln(b(x))
for(int i = 0; i < n; ++i) e[i] = a[i];
for(int i = n; i < lim; ++i) e[i] = 0;
ntt(d, lim, 1), ntt(e, lim, 1), ntt(b, lim, 1);
for(int i = 0; i < lim; ++i) b[i] = (1ll - d[i] + e[i] + mod) * b[i] % mod;
ntt(b, lim, -1);
for(int i = n; i < lim; ++i) b[i] = 0;
}
最後是多項式除法
除法由於其優秀的邊界故不和上面放到一起，並且這裡給出完整程式碼。
#include <bits/stdc++.h>
#define ll long long
#define cl const long long
using namespace std;
namespace IO{
inline ll read(){
ll x = 0;
char ch = getchar();
while(!isdigit(ch)) ch = getchar();
while(isdigit(ch)) x = (x << 3) + (x << 1) + ch - '0', ch = getchar();
return x;
}
template <typename T> inline void write(T x){
if(x > 9) write(x / 10);
putchar(x % 10 + '0');
}
}
using namespace IO;
const ll N = 3e5 + 10;
const ll mod = 998244353;
const ll G = 3, Gi = 332748118;
ll n, m;
ll f[N], g[N], fr[N], gr[N];
ll gi[N], c[N], q[N], r[N], p[N];
namespace NTT{
int lim, len;
inline ll qpow(ll a, ll gi){
ll res = 1;
while(gi){
if(gi & 1) res = res * a % mod;
gi >>= 1, a = a * a % mod;
}
return res;
}
inline void get_rev(cl n){
lim = 1, len = 0;
while(lim <= n) lim <<= 1, ++len;
for(int i = 0; i <= lim; ++i) p[i] = (p[i >> 1] >> 1) | ((i & 1) << (len - 1));
}
inline void ntt(ll A[], cl lim, cl type){
for(ll i = 0; i <= lim; ++i)
if(i < p[i]) swap(A[i], A[p[i]]);
for(ll mid = 1; mid < lim; mid <<= 1){
ll Wn = qpow(type == 1 ? G : Gi, (mod - 1) / (mid << 1));
for(ll i = 0; i < lim; i += (mid << 1)){
ll w = 1;
for(ll j = 0; j < mid; ++j, w = w * Wn % mod){
ll x = A[i + j], y = w * A[i + j + mid] % mod;
A[i + j] = (x + y) % mod;
A[i + j + mid] = (x - y + mod) % mod;
}
}
}
if(type == 1) return;
ll inv = qpow(lim, mod - 2);
for(int i = 0; i <= lim; ++i) A[i] = A[i] * inv % mod;
}
inline void Mul(cl n, cl m, ll a[], ll b[]){
get_rev(n + m);
ntt(a, lim, 1), ntt(b, lim, 1);
for(ll i = 0; i <= lim; ++i)
a[i] = a[i] * b[i] % mod;
ntt(a, lim, -1);
}
inline void Inv(cl n, ll a[], ll b[]){
if(n == 0) return b[0] = qpow(a[0], mod - 2), void();
Inv(n >> 1, a, b);
get_rev(n << 1);
for(ll i = 0; i <= n; ++i) c[i] = a[i];
for(ll i = n + 1; i <= lim; ++i) c[i] = 0;
ntt(c, lim, 1), ntt(b, lim, 1);
for(ll i = 0; i <= lim; ++i)
b[i] = (2 - c[i] * b[i] % mod + mod) * b[i] % mod;
ntt(b, lim, -1);
for(ll i = n + 1; i <= lim; ++i) b[i] = 0;
}
}
using namespace NTT;
inline void solve(){
for(ll i = n - m + 1; i <= m; ++i) gr[i] = 0;
Inv(n - m, gr, gi);
Mul(n, n - m, fr, gi);
for(ll i = 0; i <= n - m; ++i) write(q[i] = fr[n - m - i]), putchar(' ');
puts("");
Mul(m, n - m, g, q);
for(ll i = 0; i < m; ++i) write((f[i] - g[i] + mod) % mod), putchar(' ');
puts("");
}
signed main(){
n = read(), m = read();
for(ll i = 0; i <= n; ++i) f[i] = read(), fr[n - i] = f[i];
for(ll i = 0; i <= m; ++i) g[i] = read(), gr[m - i] = g[i];
solve();
return 0;
}
本文來自部落格園，作者：xixike，轉載請註明原文連結：https://www.cnblogs.com/xixike/p/15614455.html