Advertisement

中国剩余定理

阅读量:

问题引入

在《孙子算经》中有这样一个问题:“今有物不知其数,三三数之剩二(除以3余2),五五数之剩三(除以5余3),七七数之剩二(除以7余2),问物几何?”这个问题称为“孙子问题”,该问题的一般解法国际上称为“中国剩余定理”。

我们把上述问题抽象出来,就是求解同余方程组:
\begin{cases} x\equiv b_1(\operatorname{mod}a_1) \\ x\equiv b_2(\operatorname{mod}a_2) \\ x\equiv b_3(\operatorname{mod}a_3) \\ \cdots \\ x\equiv b_n(\operatorname{mod}a_n) \end{cases}

其中a_1,a_2,\cdots,a_n互质。

解法

M=\prod_{i=1}^na_i,\frac{M}{a_i}|n_in_i\equiv b_i(\operatorname{mod}a_i).

那么方程的一个解为:

\sum_{i=1}^nn_i

虽然解法难以想到,但原理却容易验证;我们设x=\sum_{i=1}^nn_i,则:

x\operatorname{mod}a_i=\sum_{j=1}^n(n_j\operatorname{mod}a_i)\operatorname{mod}a_i

由于对k\not=i都有n_k=\frac{M}{a_k}=a_1a_2\cdots a_{k-1}a_{k+1}\cdots a_n,其包含因子a_i,则n_j\operatorname{mod}a_i=0.

故:

x\operatorname{mod}a_i=\sum_{j=1}^n(n_j\operatorname{mod}a_i)\operatorname{mod}a_i=n_i\operatorname{mod}a_i=b_i

那么我们如何寻找n_i呢?我们设n_i=t\frac{M}{a_i},则:

t\frac{M}{a_i}\equiv b_i(\operatorname{mod}a_i)

由欧拉定理:

\left(\frac{M}{a_i}\right)^{\phi(a_i)}\equiv1(\operatorname{mod}a_i)

故:

t=\left(\frac{M}{a_i}\right)^{\phi(a_i)-1}

即:

n_i=\left(\frac{M}{a_i}\right)^{\phi(a_i)}

x的一个解为:

\sum_{i=1}^nn_i

其最小正整数解为:

\sum_{i=1}^nn_i\operatorname{mod} M

这是因为x-a_it\equiv b_i(\operatorname{mod}a_i).

例题

https://www.luogu.com.cn/problem/P1495

复制代码
    #include<cstdio>
    using namespace std;
    	int n,cnt;
    	int prime[500001];
    	long long phi[2000001];
    	bool vis[2000001];
    	long long m;
    	long long a[11];
    	long long b[11];
    long long multi(long long a,long long b,long long p)
    {
    return ((a * b - (long long)(((long double)a * b + 0.5) / p) * p) % p + p) % p;
    }
    long long quickpower(long long a,long long b,long long mod)
    {
    	long long t = 1;
    	while (b)
    	{
    		if (b & 1)
    			t = multi(t,a,mod);
    		a = multi(a,a,mod);
    		b >>= 1;
    	}
    	return t;
    }
    int main()
    {
    	phi[1] = 1;
    	for (int i = 2;i <= 2000000;i ++)
    	{
    		if (!vis[i])
    		{
    			prime[++cnt] = i;
    			phi[i] = i - 1;
    		}
    		for (int j = 1;j <= cnt && prime[j] * i <= 2000000;j ++)
    		{
    			vis[i * prime[j]] = true;
    			if (i % prime[j] == 0)
    			{
    				phi[i * prime[j]] = phi[i] * prime[j];
    				break;
    			}
    			phi[i * prime[j]] = phi[i] * phi[prime[j]];
    		}
    	}
    	scanf("%d",&n);
    	for (int i = 1;i <= n;i ++)
    		scanf("%lld%lld",&a[i],&b[i]);
    	m = 1;
    	for (int i = 1;i <= n;i ++)
    		m *= a[i];
    	long long ans = 0;
    	for (int i = 1;i <= n;i ++)
    	{
    		long long now = m / a[i];
    		ans += multi(multi(b[i],now,m),quickpower(now % a[i],phi[a[i]] - 1,a[i]),m);
    		ans %= m;
    	}
    	printf("%lld",ans);
    	return 0;
    }

全部评论 (0)

还没有任何评论哟~