「CmdOI2019」算力训练
2022/5/5 23:13:53
本文主要是介绍「CmdOI2019」算力训练，对大家解决编程问题具有一定的参考价值，需要的程序猿们随着小编来一起学习吧！
题目
点这里看题目。
分析
在此之前，约定 \(x=(x_{m-1}x_{m-2}x_{m-3}\dots x_0)_k\)，后者为 \(x\) 的 \(k\) 进制表示（从高位到低位），且根据题意至多 \(m\) 位。
约定 \(\omega_k\) 为 \(k\) 次单位根。
确实是大大地考验了我的算力。
显然，我们需要计算如下形式的卷积：
\[\bigoplus_{j=1}^n(1+x^{A_j}) \]在这里，我们不妨认为 \(\bigoplus\) 表示的是对于指数的 K 进制异或卷积。
如果你做过类似的题目，比如「UNR #2」黎明前的巧克力，那么应该记得，这种题目有一个普适的方法，就是直接计算答案序列 FWT 之后的结果，然后一遍 IFWT 回去。
这道题也适用于这个思路。首先考察一个序列的 FWT：
\[\operatorname{FWT}(1+x^{a})[x]=1+\prod_{j=0}^{m-1}\omega_k^{a_jx_j} \]所以，一个序列 FWT 之后的结果应当只有 \(k\) 种，形如 \(1+\omega_k^j\)。显然，我们可以尝试计算在第 \(x\) 位上，\(1+\omega_k^j\) 这样的因子出现了多少次，则最后可以直接算幂次得到结果。所以，不妨设 \(c_{x,j}\) 表示在第 \(x\) 位上，形如 \(1+\omega_k^j\) 因子出现的次数。
注意到，我们已然知道 \(\sum_j c_{x,j}=n\)，这暗示我们可以通过解方程得到各个 \(c\)。尝试利用 FWT 本身来建立一个方程：注意到，根据 FWT 运算的线性性，我们很容易得到 \(\sum_{j} c_{x,j}\omega_k^j\)，只需要考察 \(\operatorname{FWT}(\sum_{j=1}^nx^{A_j})[x]\) 即可。
再看一眼已有的两个方程：
\[\begin{cases} \sum_{j=0}^{k-1}c_{x,j}&=n\\ \sum_{j=0}^{k-1}c_{x,j}\omega_k^j&=\dots \end{cases} \]注意到，方程的系数已经有点像 DFT 矩阵的系数了。因此，我们尝试补全这个矩阵。由于 FWT 实际上就是高维 DFT，所以我们不难修改它，从而得到其余的几条方程。实际上，只需要改一下插入的值——无非就是将 \(\omega_k\) 改成 \(\omega_k^j\)——即可生成一条新的方程。
下面先考虑如何实现这个算法。由于 998244353 性质（在这里）不够好，或者 \(k=5,6\) 实在是两个精心构造的值，反正我们不能直接用原根的幂次替代 \(\omega_k\)。但是，在这种背景下，绝对不能忘记扩域这样一个选项——直接引入 \(\omega_k\) 及其幂次即可。
但是，又会有新的问题出现：由于单位根的某些性质，比如 \(\sum_{j=0}^{k-1}\omega_k^j=0\)，一个整数在扩张后的域内不能唯一表示。当然，放到几何背景下，这其实就说明某些（我们认为是基向量的）向量可以被其它向量线性表示。解决办法也比较简单，消去多余的元（在复平面上画图观察）即可。
最终的结果是：在 \(k=5\) 的时候，我们只需要留下 \(\omega_5^0,\omega_5^1,\omega_5^2,\omega_5^3\)；而在 \(k=6\) 的时候，我们只需要留下 \(\omega_k^0,\omega_k^2\)。
最后是一个观察——当我们试图生成下式：
\[\sum_{j=0}^{k-1}c_{x,j}\omega_k^{vj} \]的值的时候，我们实际上是在计算：
\[\sum_{i=1}^n\prod_{j=0}^{m-1}(\omega_k^v)^{x_j(a_i)_j} \]如果每次都换一下插入的值实在是太浪费了。注意到 \((\omega_k^v)^{x_j}=\omega_k^{vx_j}\)，所以我们可以只需要插入 \(\omega_k\) 计算一次，就可以通过下标移位找到需要的值，从而降一个 \(k\)。
注意，这里不能把插入 \(\omega_k^v\) 看作是“复合”上一个新的元，因为外部的运算并不会满足真正的“复合”的逻辑。
最终复杂度大约为 \(O(k^{m+2}m)\)。
小结：
-
扩域的想法，或者说是解决问题的手段，要积累一下。
-
简单序列 FWT 的结果可以直接推导，这个是一个常见的思路；
考察不同的项的指数，通过解方程得到指数，也算是一个常用的思路，当然关键在于如何建立方程。
-
最后一步的观察，思路比较朴素，也就是要充分利用计算出来了的信息。支撑这个方法的还有一点，就是整个下标集在按位运算下是封闭的。
代码
#include <cmath>
#include <cstdio>
#define rep( i, a, b ) for( int i = (a) ; i <= (b) ; i ++ )
#define per( i, a, b ) for( int i = (a) ; i >= (b) ; i -- )
const int mod = 998244353, g = 3, phi = mod - 1;
const int MAXN = 1e5 + 5, MAXRT = 1e3 + 5;
template<typename _T>
void read( _T &x ) {
x = 0; char s = getchar(); bool f = false;
while( ! ( '0' <= s && s <= '9' ) ) { f = s == '-', s = getchar(); }
while( '0' <= s && s <= '9' ) { x = ( x << 3 ) + ( x << 1 ) + ( s - '0' ), s = getchar(); }
if( f ) x = -x;
}
template<typename _T>
void write( _T x ) {
if( x < 0 ) putchar( '-' ), x = -x;
if( 9 < x ) write( x / 10 );
putchar( x % 10 + '0' );
}
template<typename _T>
struct SwiftPow {
_T sml[MAXRT], lrg[MAXRT];
int B;
SwiftPow(): B( 1 ) {}
inline void Init( const _T &x, const int &n ) {
B = ceil( sqrt( n ) );
lrg[0] = sml[0] = _T( 1 );
rep( i, 1, B ) sml[i] = sml[i - 1] * x;
lrg[1] = sml[B];
rep( i, 2, B ) lrg[i] = lrg[i - 1] * lrg[1];
}
inline _T operator () ( const int &n ) const {
return lrg[n / B] * sml[n % B];
}
};
int app[MAXN];
int N, K, M;
inline int Qkpow( int, int );
inline int Inv( const int &a ) { return Qkpow( a, mod - 2 ); }
inline int Mul( int x, const int &v ) { return 1ll * x * v % mod; }
inline int Sub( int x, const int &v ) { return ( x -= v ) < 0 ? x + mod : x; }
inline int Add( int x, const int &v ) { return ( x += v ) >= mod ? x - mod : x; }
inline int& MulEq( int &x, const int &v ) { return x = 1ll * x * v % mod; }
inline int& SubEq( int &x, const int &v ) { return ( x -= v ) < 0 ? ( x += mod ) : x; }
inline int& AddEq( int &x, const int &v ) { return ( x += v ) >= mod ? ( x -= mod ) : x; }
inline int Qkpow( int base, int indx ) {
int ret = 1;
while( indx ) {
if( indx & 1 ) MulEq( ret, base );
MulEq( base, base ), indx >>= 1;
}
return ret;
}
inline int readK() {
static char buf[20];
int ret = 0;
scanf( "%s", buf );
for( int i = 0 ; buf[i] ; i ++ )
ret = ret * K + ( int ) ( buf[i] - '0' );
return ret;
}
namespace KEq5 {
struct Complex {
int a[4];
Complex(): a{} {}
Complex( int R ): a{} { a[0] = R; }
inline int& operator [] ( const int &idx ) { return a[idx]; }
inline const int operator [] ( const int &idx ) const { return a[idx]; }
inline Complex& operator += ( const Complex &b ) {
rep( i, 0, 3 ) AddEq( a[i], b[i] );
return *this;
}
inline Complex& operator -= ( const Complex &b ) {
rep( i, 0, 3 ) SubEq( a[i], b[i] );
return *this;
}
inline Complex& operator /= ( int b ) {
b = Inv( b );
rep( i, 0, 3 ) MulEq( a[i], b );
return *this;
}
inline Complex& operator *= ( const Complex &b ) {
return *this = *this * b;
}
inline Complex operator + ( const Complex &b ) const {
Complex ret;
rep( i, 0, 3 ) ret[i] = Add( a[i], b[i] );
return ret;
}
inline Complex operator - ( const Complex &b ) const {
Complex ret;
rep( i, 0, 3 ) ret[i] = Sub( a[i], b[i] );
return ret;
}
inline Complex operator * ( const Complex &b ) const {
Complex ret; int out = 0;
rep( i, 0, 3 ) if( a[i] )
rep( j, 0, 3 ) if( b[j] ) {
if( ( i + j ) % 5 == 4 )
AddEq( out, Mul( a[i], b[j] ) );
else
AddEq( ret[( i + j ) % 5], Mul( a[i], b[j] ) );
}
rep( i, 0, 3 ) SubEq( ret[i], out );
return ret;
}
inline Complex operator / ( int b ) const {
Complex ret; b = Inv( b );
rep( i, 0, 3 ) ret[i] = Mul( a[i], b );
return ret;
}
};
Complex coe[MAXN], nw[MAXN], w[5], res[5];
SwiftPow<Complex> pw[5];
int lim;
Complex Work( const int x ) {
static int val[5] = {};
static int inv5 = Inv( 5 );
Complex f = 0, ret = 1;
rep( i, 0, K - 1 ) {
val[i] = 0;
for( int j = 0, s = 1 ; j < M ; j ++, s *= K )
val[i] += s * ( ( x / s % 5 ) * i % 5 );
}
rep( i, 0, K - 1 ) {
f = 0;
rep( j, 0, K - 1 )
f += coe[val[j]] * w[( K - i * j % K ) % K];
ret *= pw[i]( Mul( inv5, f.a[0] ) );
}
return ret;
}
void Init() {
lim = 1;
rep( i, 1, M ) lim *= K;
rep( i, 0, lim - 1 ) coe[i] = app[i];
rep( i, 0, 3 ) w[i].a[i] = 1, w[4].a[i] = mod - 1;
rep( i, 0, 4 ) pw[i].Init( w[i] + 1, N );
}
void Solve() {
Init(); Complex tmp;
for( int s = 1 ; s < lim ; s *= K )
for( int i = 0 ; i < lim ; i += s * K )
for( int j = i ; j < i + s ; j ++ ) {
for( int k = 0 ; k < 5 ; k ++ ) {
tmp = 0;
for( int t = 0 ; t < 5 ; t ++ )
tmp += coe[j + t * s] * w[k * t % 5];
res[k] = tmp;
}
for( int k = 0 ; k < 5 ; k ++ )
coe[j + k * s] = res[k];
}
rep( i, 0, lim - 1 ) nw[i] = Work( i );
rep( i, 0, lim - 1 ) coe[i] = nw[i];
for( int s = 1 ; s < lim ; s *= K )
for( int i = 0 ; i < lim ; i += s * K )
for( int j = i ; j < i + s ; j ++ ) {
for( int k = 0 ; k < 5 ; k ++ ) {
tmp = 0;
for( int t = 0 ; t < 5 ; t ++ )
tmp += coe[j + t * s] * w[( 5 - k * t % 5 ) % 5];
res[k] = tmp / 5;
}
for( int k = 0 ; k < 5 ; k ++ )
coe[j + k * s] = res[k];
}
rep( i, 0, lim - 1 ) write( coe[i].a[0] ), putchar( '\n' );
}
}
namespace KEq6 {
struct Complex {
int w0, w2;
Complex(): w0( 0 ), w2( 0 ) {}
Complex( int R ): w0( R ), w2( 0 ) {}
Complex( int W0, int W2 ): w0( W0 ), w2( W2 ) {}
inline Complex& operator += ( const Complex &b ) {
AddEq( w0, b.w0 );
AddEq( w2, b.w2 );
return *this;
}
inline Complex& operator -= ( const Complex &b ) {
SubEq( w0, b.w0 );
SubEq( w2, b.w2 );
return *this;
}
inline Complex& operator /= ( int b ) {
b = Inv( b );
MulEq( w0, b );
MulEq( w2, b );
return *this;
}
inline Complex& operator *= ( const Complex &b ) {
return *this = *this * b;
}
inline Complex operator + ( const Complex &b ) const {
return Complex( Add( w0, b.w0 ), Add( w2, b.w2 ) );
}
inline Complex operator - ( const Complex &b ) const {
return Complex( Sub( w0, b.w0 ), Sub( w2, b.w2 ) );
}
inline Complex operator * ( const Complex &b ) const {
return Complex( Sub( Mul( w0, b.w0 ), Mul( w2, b.w2 ) ),
Sub( Add( Mul( w2, b.w0 ), Mul( w0, b.w2 ) ), Mul( w2, b.w2 ) ) );
}
inline Complex operator / ( int b ) const {
b = Inv( b );
return Complex( Mul( w0, b ), Mul( w2, b ) );
}
};
Complex nw[MAXN], coe[MAXN], w[6], res[6];
SwiftPow<Complex> pw[6];
int lim;
Complex Work( const int x ) {
static int val[6] = {};
static int inv6 = Inv( 6 );
Complex f = 0, ret = 1;
rep( i, 0, K - 1 ) {
val[i] = 0;
for( int j = 0, s = 1 ; j < M ; j ++, s *= K )
val[i] += s * ( ( x / s % K ) * i % K );
}
rep( i, 0, K - 1 ) {
f = 0;
rep( j, 0, K - 1 )
f += coe[val[j]] * w[( K - i * j % K ) % K];
ret *= pw[i]( Mul( inv6, f.w0 ) );
}
return ret;
}
void Init() {
lim = 1;
rep( i, 1, M ) lim *= K;
rep( i, 0, lim - 1 ) coe[i] = app[i];
w[0] = Complex( 1, 0 );
w[1] = Complex( 1, 1 );
w[2] = Complex( 0, 1 );
w[3] = Complex( mod - 1, 0 );
w[4] = Complex( mod - 1, mod - 1 );
w[5] = Complex( 0, mod - 1 );
rep( i, 0, 5 ) pw[i].Init( w[i] + 1, N );
}
void Solve() {
Init(); Complex tmp;
for( int s = 1 ; s < lim ; s *= K )
for( int i = 0 ; i < lim ; i += s * K )
for( int j = i ; j < i + s ; j ++ ) {
for( int k = 0 ; k < 6 ; k ++ ) {
tmp = 0;
for( int t = 0 ; t < 6 ; t ++ )
tmp += coe[j + t * s] * w[k * t % 6];
res[k] = tmp;
}
for( int k = 0 ; k < 6 ; k ++ )
coe[j + k * s] = res[k];
}
rep( i, 0, lim - 1 ) nw[i] = Work( i );
rep( i, 0, lim - 1 ) coe[i] = nw[i];
for( int s = 1 ; s < lim ; s *= K )
for( int i = 0 ; i < lim ; i += s * K )
for( int j = i ; j < i + s ; j ++ ) {
for( int k = 0 ; k < 6 ; k ++ ) {
tmp = 0;
for( int t = 0 ; t < 6 ; t ++ )
tmp += coe[j + t * s] * w[( 6 - k * t % 6 ) % 6];
res[k] = tmp / 6;
}
for( int k = 0 ; k < 6 ; k ++ )
coe[j + k * s] = res[k];
}
rep( i, 0, lim - 1 ) write( coe[i].w0 ), putchar( '\n' );
}
}
int main() {
read( N ), read( K ), read( M );
rep( i, 1, N ) AddEq( app[readK()], 1 );
if( K == 5 ) KEq5 :: Solve();
if( K == 6 ) KEq6 :: Solve();
return 0;
}
这篇关于「CmdOI2019」算力训练的文章就介绍到这儿，希望我们推荐的文章对大家有所帮助，也希望大家多多支持为之网！
- 2024-06-06Package Easy（基于 NSIS 的打包exe安装包工具）使用方法-icode9专业技术文章分享
- 2024-06-06基于 casdoor 的 ELK 开源登录认证解决方案: elk-auth-casdoor-icode9专业技术文章分享
- 2024-05-29Elasticsearch慢查询日志配置
- 2024-05-29揭秘华为如此多成功项目的产品关键——Charter模板
- 2024-05-29海外IDC业务拓展的7大挑战
- 2024-05-29InLine Chat功能优化对标Github Copilot，CodeGeeX带来更高效、更直观的编程体验！
- 2024-05-29CodeGeeX 智能编程助手 6 项功能升级，在Visual Studio插件市场霸榜2周！
- 2024-05-29AutoMQ 生态集成 Apache Doris
- 2024-05-292024年IDC行业的深度挖掘：机遇、挑战与未来展望
- 2024-05-29五款扩展组件齐发 —— Volcano、Keda、Crane-scheduler 等，邀你体验
