当前位置:网站首页>[cf1054h] epic Revolution -- number theory, convolution, arbitrary modulus NTT
[cf1054h] epic Revolution -- number theory, convolution, arbitrary modulus NTT
2022-07-29 07:00:00 【Ouye xjx】
CF1054H Epic Convolution
Answer key
This question is called Epic Convolution, Which translates as “ Epic convolution ”( Also called “ Shit shit convolution ”)
Because the modulus is very small and a prime number , So by Euler's theorem , We can turn the problem into a quest
f x = ∑ i = 0 n − 1 ∑ j = 0 m − 1 A i B j [ i 2 j 3 m o d 490018 = x ] f_x=\sum_{i=0}^{n-1}\sum_{j=0}^{m-1}A_iB_j[i^2j^3\bmod 490018=x] fx=i=0∑n−1j=0∑m−1AiBj[i2j3mod490018=x] Then in order to solve more efficiently , We need to consider transforming the above multiplicative convolution into additive convolution .
can 490018 490018 490018 Obviously, there is no original root ? How to convert ?
At this time, the idea of dividing the range into Multiplication rings is obviously unrealistic ( That's why I'm a mathematical shit ), We need more abstract ideas . Consider putting 490018 490018 490018 Make a unique decomposition :
490018 = 2 × 491 × 499 490018=2\times 491 \times 499 490018=2×491×499 And when we get i 2 j 3 i^2j^3 i2j3 modulus 2 , 491 , 499 2,491,499 2,491,499 After the value of , It can be used CRT combined , So we might as well change the formula , Change to seek
f x , y , z = ∑ i = 0 n − 1 ∑ j = 0 m − 1 A i B j [ i 2 j 3 m o d 2 = x ] [ i 2 j 3 m o d 491 = y ] [ i 2 j 3 m o d 499 = z ] f_{x,y,z}=\sum_{i=0}^{n-1}\sum_{j=0}^{m-1}A_iB_j[i^2j^3\bmod 2=x][i^2j^3\bmod 491=y][i^2j^3\bmod 499=z] fx,y,z=i=0∑n−1j=0∑m−1AiBj[i2j3mod2=x][i2j3mod491=y][i2j3mod499=z]
The rest is simple , We just need to convolute each dimension in turn , The first dimension can be violently convoluted , The latter two dimensions can be transformed into additive convolution by taking logarithm of the original root .
Notice that 0 Cases where logarithms cannot be taken , stay 0 Place the convolution of violence .
Code
What I write is the convolution of violence classification discussion , In fact, a little bit of exclusion can significantly reduce the number of operations .
#include<bits/stdc++.h>//JZM yyds!!
#define ll long long
#define lll __int128
#define uns unsigned
#define fi first
#define se second
#define IF (it->fi)
#define IS (it->se)
#define lowbit(x) ((x)&-(x))
#define END putchar('\n')
#define inline jzmyyds
using namespace std;
const int MAXN=5e5+5;
const ll INF=1e17;
ll read(){
ll x=0;bool f=1;char s=getchar();
while((s<'0'||s>'9')&&s>0){
if(s=='-')f^=1;s=getchar();}
while(s>='0'&&s<='9')x=(x<<1)+(x<<3)+(s^48),s=getchar();
return f?x:-x;
}
int ptf[50],lpt;
void print(ll x,char c='\n'){
if(x<0)putchar('-'),x=-x;
ptf[lpt=1]=x%10;
while(x>9)x/=10,ptf[++lpt]=x%10;
while(lpt)putchar(ptf[lpt--]^48);
if(c>0)putchar(c);
}
//490018=2*491*499,g_491=2,g_499=7
const ll MOD=490019,NM=1000000000245104641ll;// One pass large modulus NTT
const lll E=1;
ll ksml(ll a,ll b,ll mo=NM){
ll res=1;
for(;b;b>>=1,a=E*a*a%mo)if(b&1)res=E*res*a%mo;
return res;
}
ll ksm(ll a,ll b,ll mo=MOD){
ll res=1;
for(;b;b>>=1,a=a*a%mo)if(b&1)res=res*a%mo;
return res;
}
#define g 7ll
int rev[1025];
int initrev(){
for(int i=0;i<(1<<10);i++)rev[i]=(rev[i>>1]>>1)|((i&1)<<9);
return 114514;
}
ll omg[1025],cbddl=initrev(),IVN=ksml(1024,NM-2);
void NTT(ll*a,int inv,int n=1024){
for(int i=0;i<n;i++)if(i<rev[i])swap(a[i],a[rev[i]]);
ll x,y,tmp=233,mi=(NM-1)>>1;omg[0]=1;
for(int m=1;m<n;m<<=1,mi>>=1){
if(m>1)tmp=ksml(g,inv>0?mi:NM-1-mi);
for(int i=1;i<m;i++)omg[i]=E*omg[i-1]*tmp%NM;
for(int i=0,om=0;i<n;i+=(m<<1),om=0)
for(int j=i;j<i+m;j++,om++){
x=a[j],y=E*a[j+m]*omg[om]%NM,a[j]=x+y,a[j+m]=x-y;
if(a[j]>=NM)a[j]-=NM;
if(a[j+m]<0)a[j+m]+=NM;
}
}if(inv<0)for(int i=0;i<n;i++)a[i]=E*a[i]*IVN%NM;
}
#undef g
ll f[1025][1025],g[1025][1025];
void polymul(){
for(int i=0;i<500;i++)NTT(f[i],1),NTT(g[i],1);
for(int i=1;i<1024;i++)
for(int j=0;j<i;j++)swap(f[i][j],f[j][i]),swap(g[i][j],g[j][i]);
for(int i=0;i<1024;i++){
NTT(f[i],1),NTT(g[i],1);
for(int j=0;j<1024;j++)f[i][j]=E*f[i][j]*g[i][j]%NM;
NTT(f[i],-1);
}
for(int i=1;i<1024;i++)for(int j=0;j<i;j++)swap(f[i][j],f[j][i]);
for(int i=0;i<1000;i++)NTT(f[i],-1);
}
int ex1[810],lg1[810],ex2[810],lg2[810];
int initexp(){
ex1[0]=ex2[0]=1,lg1[1]=lg2[1]=0;
for(int i=1;i<490;i++)ex1[i]=(ex1[i-1]<<1)%491,lg1[ex1[i]]=i;
for(int i=1;i<498;i++)ex2[i]=ex2[i-1]*7%499,lg2[ex2[i]]=i;
return 1919810;
}
int sh_t=initexp(),n,m;
ll a1[499],b1[499],ans,mi[MAXN],A[MAXN];
void solve(ll(*a)[499],ll(*b)[499],ll(*c)[499],bool CL){
if(CL){
memset(f,0,sizeof(f));
memset(g,0,sizeof(g));
memset(a1,0,sizeof(a1));
memset(b1,0,sizeof(b1));
}
for(int i=1;i<491;i++)
for(int j=1;j<499;j++)
(f[lg1[i]][lg2[j]]+=a[i][j])%=MOD,(g[lg1[i]][lg2[j]]+=b[i][j])%=MOD;
for(int i=1;i<491;i++)
for(int j=1;j<499;j++)(a1[i]+=a[i][j])%=MOD;
for(int i=1;i<491;i++)
for(int j=0;j<499;j++)(b1[i]+=b[i][j])%=MOD;
polymul();
for(int i=0,mi=0;i<1000;i++,mi=i%490)
for(int j=0;j<1000;j++)
(c[ex1[mi]][ex2[j%498]]+=f[i][j])%=MOD;
for(int i=1;i<491;i++)
for(int j=1,to=i;j<491;j++,to=(to+i)%491)
(c[to][0]+=a[i][0]*b1[j]+a1[i]*b[j][0])%=MOD;
memset(a1,0,sizeof(a1));
memset(b1,0,sizeof(b1));
for(int i=1;i<491;i++)
for(int j=0;j<499;j++)(a1[j]+=a[i][j])%=MOD;
for(int i=0;i<491;i++)
for(int j=0;j<499;j++)(b1[j]+=b[i][j])%=MOD;
for(int i=0;i<499;i++)
for(int j=0,to=0;j<499;j++,to=(to+i)%499)
(c[0][to]+=a[0][i]*b1[j]+a1[i]*b[0][j])%=MOD;
}
int crt(int a,int b,int c){
static const ll iv1=492>>1,iv2=ksm(491<<1,499-2,499);
int k1=(b-a+491)*iv1%491,k2=(c-a-k1-k1+4990)*iv2%499;
return (a+(k1<<1)+(k2*491<<1))%(MOD-1);
}
ll a[2][491][499],b[2][491][499],b_[491][499],as[2][491][499],c;
const bool KFC=1;
int main()
{
n=read(),m=read(),c=read(),mi[0]=1;
for(int i=1;i<MOD;i++)mi[i]=mi[i-1]*c%MOD;
if(1ll*n*m<=1000000ll&&KFC){
for(int i=0;i<n;i++)A[i]=read();
for(int j=0;j<m;j++){
ll B=read();
for(int i=0;i<n;i++)
(ans+=A[i]*B*mi[1ll*i*i*j%(MOD-1)*j*j%(MOD-1)])%=MOD;
}return print(ans),0;
}
for(int i=0,k=0;i<n;i++,k=1ll*i*i%(MOD-1))
(a[k&1][k%491][k%499]+=read())%=MOD;
for(int i=0,k=0,d;i<m;i++,k=1ll*i*i*i%(MOD-1))
d=read(),(b[k&1][k%491][k%499]+=d)%=MOD,(b_[k%491][k%499]+=d)%=MOD;
solve(a[1],b[1],as[1],0);
solve(a[0],b_,as[0],1);
solve(a[1],b[0],as[0],1);
for(int i=0;i<2;i++)for(int j=0;j<491;j++)
for(int k=0;k<499;k++)(ans+=as[i][j][k]*mi[crt(i,j,k)])%=MOD;
print(ans);
return 0;
}
边栏推荐
- 联邦学习后门攻击总结(2019-2022)
- 【干货备忘】50种Matplotlib科研论文绘图合集,含代码实现
- 吴恩达老师机器学习课程笔记 00 写在前面
- Idea cannot find a database solution
- Overview of database system
- MySql基础知识(高频面试题)
- Teacher Wu Enda's machine learning course notes 00 are written in the front
- Simulation volume leetcode [general] 150. evaluation of inverse Polish expression
- Shallow reading of shared lock source code of abstractqueuedsynchronizer (AQS)
- Junda technology | applicable to "riyueyuan" brand ups wechat cloud monitoring card
猜你喜欢

Is online legend software testing training really so black hearted? Are they all scams?

JVM之垃圾回收机制(GC)

5g service interface and reference point

Idea cannot find a database solution

线程同步—— 生产者与消费者、龟兔赛跑、双线程打印

Unity free element special effect recommendation

Windows 上 php 7.4 连接 oracle 配置

Actual combat! Talk about how to solve the deep paging problem of MySQL

二次元卡通渲染——进阶技巧

【冷冻电镜】Relion4.0——subtomogram教程
随机推荐
矩阵分解与梯度下降
说一下 TCP/IP 协议?以及每层的作用?
Security in quantum machine learning
Basic knowledge of MySQL (high frequency interview questions)
Sword finger offer II 115: reconstruction sequence
vscode通过remotessh结合xdebug远程调试php解决方案
模拟卷Leetcode【普通】061. 旋转链表
【冷冻电镜|论文阅读】emClarity:用于高分辨率冷冻电子断层扫描和子断层平均的软件
Share some tips for better code, smooth coding and improve efficiency
Teacher wangshuyao's notes on operations research 04 fundamentals of linear algebra
好文佳句摘录
Difference between CNAME record and a record
Teacher wangshuyao's notes on operations research 05 linear programming and simplex method (concept, modeling, standard type)
Invalid access control
线程同步—— 生产者与消费者、龟兔赛跑、双线程打印
剑指 Offer II 115:重建序列
如何优雅的写 Controller 层代码?
MySql基础知识(高频面试题)
【经验】通过跳板机远程连接内网服务器的相关配置
Ali gave several SQL messages and asked how many tree search operations need to be performed?