题意
给一个n阶行列式，求行列式的值对m取余的值
思路
高斯消元，但是这个行列式的值不一定与m互质，所以如果用互质的逆元求法快速幂来做会有问题，可以将m分解质因数，然后分别求行列式模这些质因数的值，再使用中国剩余定理合成回来
代码
#include <cstdio>
#include <cmath>
#include <cstdlib>
#include <iostream>
#include <vector>
using std::
vector;
long long D[
101][
101],E[
101][
101];
long long n,m,ans,dv;
long long gcd(
long long x,
long long y)
{
if(x%y==
0)
return y;
if(y%x==
0)
return x;
if(x>y)
return gcd(y,x%y);
else return gcd(x,y%x);
}
void Gauss()
{
long long t,temp,g,t1,t2,p,q;
for(
long long i=
0,j=
0;i<n&&j<n;i++,j++)
{
t=i;
for(
long long k=i;k<n;k++)
if(D[k][j]!=
0)
{
t=k;
break;
}
if(D[t][j]==
0)
{
i--;
continue;
}
if(t!=i)
{
ans=-ans;
for(
long long k=
0;k<n;k++)
{
temp=D[i][k];
D[i][k]=D[t][k];
D[t][k]=temp;
}
}
for(
long long k=i+
1;k<n;k++)
if(D[k][j]!=
0)
{
if(D[i][j]>=
0)
p=D[i][j];
else p=-D[i][j];
if(D[k][j]>=
0)
q=D[k][j];
else q=-D[k][j];
g=gcd(p,q);
t1=D[k][j]/g;
t2=D[i][j]/g;
if(t2<
0)
{
t1=-t1;
t2=-t2;
}
t1=(t1%m+m)%m;
t2=(t2%m+m)%m;
dv=(dv*(t2%m))%m;
for(
long long l=
0;l<n;l++)
D[k][l]=(((D[k][l]*(t2%m))%m)-((D[i][l]*(t1%m))%m)+m)%m;
}
}
return;
}
long long power(
long long x,
long long y,
long long mo)
{
long long ret=
1;
while(y>
0)
{
if(y%
2)
ret=(ret*x)%mo;
x=(x*x)%mo;
y/=
2;
}
return ret;
}
vector<long long> mm;
vector<long long> re;
long long mq;
long long chinaremain()
{
long long ret=
0,cr=
1,temp;
for(
long long i=
0;i<mm.size();i++)
{
temp=re[i];
temp=(temp*power(mq/mm[i],mm[i]-
2,mq))%mq;
temp=(temp*(mq/mm[i]))%mq;
ret=(ret+temp)%mq;
}
return ret;
}
int main()
{
long long e,anss;
while(
std::
cin>>n>>m)
{
for(
long long i=
0;i<n;i++)
for(
long long j=
0;j<n;j++)
std::
cin>>D[i][j];
for(
long long i=
0;i<n;i++)
for(
long long j=
0;j<n;j++)
E[i][j]=D[i][j];
e=
sqrt(m);
mq=m;
for(
long long i=
2;i<=e;i++)
if(m%i==
0)
{
mm.push_back(i);
while(m%i==
0)
m/=i;
}
if(m!=
1)
mm.push_back(m);
anss=
1;
for(
long long i=
0;i<mm.size();i++)
{
ans=
1;
dv=
1;
m=mm[i];
for(
long long j=
0;j<n;j++)
for(
long long k=
0;k<n;k++)
D[j][k]=E[j][k];
for(
long long j=
0;j<n;j++)
for(
long long k=
0;k<n;k++)
if(D[j][k]<
0)
D[j][k]=(D[j][k]-m*(D[j][k]/m)+m)%m;
else if(D[j][k]>=m)
D[j][k]%=m;
Gauss();
for(
long long j=
0;j<n;j++)
ans=(ans*D[j][j])%m;
re.push_back((ans*power(dv,m-
2,m)%m+m)%m);
}
std::
cout<<chinaremain()<<
std::endl;
mm.clear();
re.clear();
}
return 0;
}
