WC 2020 【取数游戏】
首先根据原问题把图建出来如果(exist m,a_i^m=a_jpmod p),则(i)到(j)有边。
可以发现这个图有一个性质：若(i)到(j)存在边,j到k存在边，则i到k存在边。
然后我们将其拓扑排序（若有一个团，则任意排）。
则第(i)个点对答案的贡献为(2^{n-1-cnt})，cnt为能到达它的点数。
那么问题转化成：
如何快速判断两个数是否可以转换？
首先可以发现(P)一定存在原根。
设(g)为(P)的原根。
一个重要的结论：
({g^0,g^1...g^{phi(p)-1}}={x|gcd(p,x)=0})
(proof:)
若(g^x=k*b+l*(p/k)*k),则可以发现(g^x)是(k)的倍数，其中k是(p)的因子，由于(g,p)互质，所以不可能。
所以将和(p)互质的用BSGS求算离散对数，其它的暴力连边就可以了。
/*
{
######################
# Author #
# Gary #
# 2021 #
######################
*/
#include<bits/stdc++.h>
#define rb(a,b,c) for(int a=b;a<=c;++a)
#define rl(a,b,c) for(int a=b;a>=c;--a)
#define LL long long
#define IT iterator
#define PB push_back
#define II(a,b) make_pair(a,b)
#define FIR first
#define SEC second
#define FREO freopen("check.out","w",stdout)
#define rep(a,b) for(int a=0;a<b;++a)
#define SRAND mt19937 rng(chrono::steady_clock::now().time_since_epoch().count())
#define random(a) rng()%a
#define ALL(a) a.begin(),a.end()
#define POB pop_back
#define ff fflush(stdout)
#define fastio ios::sync_with_stdio(false)
#define check_min(a,b) a=min(a,b)
#define check_max(a,b) a=max(a,b)
using namespace std;
//inline int read(){
// int x=0;
// char ch=getchar();
// while(ch<'0'||ch>'9'){
// ch=getchar();
// }
// while(ch>='0'&&ch<='9'){
// x=(x<<1)+(x<<3)+(ch^48);
// ch=getchar();
// }
// return x;
//}
const int INF=0x3f3f3f3f;
typedef pair<int,int> mp;
/*}
*/
int n,p;
int a[5005];
int pw2[5005];
bool e[5005][5005];
int cnt[5005];
int ind[5005];
vector<int> prime_div(int x){
vector<int> ret;
rb(i,2,x){
if(x%i==0){
ret.PB(i);
while(x%i==0) x/=i;
}
}
return ret;
}
int phi(int x){
vector<int> tmp=prime_div(x);
for(auto it:tmp) x/=it;
for(auto it:tmp) x*=(it-1);
return x;
}
LL quick(LL A,LL B){
if(B==0) return 1;
LL tmp=quick(A,B>>1);
tmp*=tmp;
tmp%=p;
if(B&1)
tmp*=A,tmp%=p;
return tmp;
}
int gr(int x){
int Phi=phi(x);
vector<int > pf=prime_div(Phi);
rb(i,2,x){
if(__gcd(i,Phi)!=1) continue;
bool ok=1;
for(auto it:pf){
if(quick(i,Phi/it)==1){
ok=0;
break;
}
}
if(ok){
return i;
}
}
return 0;
}
int g;
unordered_map<int,int> hs;
int t,m;
int bsgs(int B,int mod){
int now;
B=quick(B,mod-2);
int tmp=quick(g,t);
now=1ll*tmp*B%mod;
rb(i,1,m){
// cout<<i<<' '<<now<<' '<<tmp<<' '<<t<<endl;
if(hs.find(now)!=hs.end()){
return i*t-hs[now];
}
now=1ll*now*tmp%mod;
}
cout<<"!"<<endl;
return -1;
}
bool ok[5005];
unordered_map<int,int> um;
int main(){
scanf("%d%d",&n,&p);
rb(i,1,n) scanf("%d",&a[i]);
pw2[0]=1;
rb(i,1,n) pw2[i]=pw2[i-1]<<1,pw2[i]%=998244353;
g=gr(p);
int now=1;
rb(i,1,n) ok[i]=(__gcd(p,a[i])==1);
rb(i,1,n) if(!ok[i]) um[a[i]]=i;
t=sqrt(double(p))+3;
m=p/t+3;
rb(i,0,t-1){
hs[now]=i;
now=1ll*now*g%p;
}
rb(i,1,n){
if(ok[i])
a[i]=bsgs(a[i],p);
}
int phi_=phi(p);
rb(i,1,n){
if(!ok[i]) continue;
int gap=__gcd(a[i],phi_);
rb(j,1,n) if(ok[j]&&i!=j) if(a[j]%gap==0) e[i][j]=1;
}
int k=0;
vector<int> pf=prime_div(p);
int p0=pf[0];
int savep=p;
while(p!=1){
k++;
p/=p0;
}
p=savep;
rb(i,1,n){
if(!ok[i]){
int tmp=a[i];
rb(T,1,k-1){
if(um.find(tmp)!=um.end()){
int id=um[tmp];
if(id!=i){
e[i][id]=1;
}
}
// cout<<i<<' '<<tmp<<endl;
tmp=1ll*tmp*a[i]%p;
}
}
}
// rb(i,1,n)
// rb(j,1,n){
// cout<<e[i][j];
// if(j==n){
// cout<<endl;
// }
// }
int rest=0;
rb(i,1,n){
rb(j,1,n){
if(e[j][i])
{
if(e[i][j]&&j>i) continue;
cnt[i]++;
}
}
}
rb(i,1,n) rest+=pw2[n-1-cnt[i]],rest%=998244353;
cout<<rest<<endl;
return 0;
}