数论相关
整除、同余
设\(a\),\(b\)为整数,\(a\neq 0\) 如果存在一个整数q，使得\(a\times q=b\),则\(b\)能被\(a\)整除，记为\(a\mid b\)，且称\(b\)是\(a\)的倍数，\(a\)是\(b\)的因子。
整除的几个性质：
- 传递性：如果\(a|b\)且\(b|c\)，则\(a|c\)
- a|b且a|c等价于对于任意的整数x，y，有a|(bx+cy)
- 设m不为0，则a|b等价于ma|mb
- 设整数x,y满足下式：ax+by=1，且a|n，b|n，那么(ab)|n
- 若b=q*d+c，那么d|b的充要条件是d|c.
同余的几个性质：
- 同加性：若\(a\equiv b \pmod p,\)则\(a+c \equiv b+c \pmod p\)
- 同减性：若\(a\equiv b \pmod p,\)则\(a-c \equiv b-c \pmod p\)
- 同乘性：若\(a\equiv b \pmod p,\)则\(a \times c \equiv b\times c \pmod p\)
- 同除性：若\(a\equiv b \pmod p\),且\(c \mid a,c\mid b,(c,p)=1\),则\(a/c \equiv b/c \pmod p\)
- 同幂性：若\(a\equiv b \pmod p,c>0,\)则\(a^c \equiv b^c \pmod p\)
- 若\(a\%p=x,a\%q=x\),且\(p,q\)互质，则\(a\%(p*q)=x\)
数论小常识：
- 若\(2\)能整除\(a\)的最末位，则\(2|a\)
- 若\(4\)能整除\(a\)的末两位，则\(4|a\)
- 若\(8\)能整除\(a\)的末三位，则\(8|a\)
- 若\(3\)能整除\(a\)的各位数字之和，则\(3|a\)
- 若\(9\)能整除\(a\)的各位数字之和，则\(9|a\)
- 若\(11\)能整除\(a\)的偶数位数字之和与奇数位数字之和的差，则\(11|a\)
- 能被\(7\)、\(11\)、\(13\)整除的数的特征是：这个数的末三位与末三位以前的数字所组成数之差能被\(7\)、\(11\)、\(13\)整除。
模运算
对于整数a，b，其中b不为0，求a除以b的余数，称为a模b，记为\(a\%b\).
模运算的性质：
- 分配率：模运算对加、减、乘具有分配率
\((a+b)\%c=(a\%c+b\%c)\%c\)
\((a-b)\%c=(a\%c-b\%c)\%c\)
\((a*b)\%c=(a\%c*b\%c)\%c\)
\((a^b)\%c=(a\%c)^b\%c\) - 放缩性
如果\(a\%b=c,d \neq 0\),则有\((a*d)\%(b*d)=c*d\);
如果\(a\%b=c,d\mid a,d \mid b\),则\((a/d)\%(b/d)=(c/d)\).
根据放缩性，则除法取余有这个式子：
\(a/b\%c=a\%(b*c)/b\)
快速幂取模
已知\(a,b,p\)为整数，求\(a^b\%p\)的结果。
可以用二进制倍增的思想，快速求幂。
a的b次方，可以看做是b个a相乘，b可以拆成2的幂的和，于是可以求出a的1次方，a的2次方，a的4次方,\(\dots\)，只要b的二进制的从右数第k位为1，则a的\(2^k\)方乘到结果里即可。
代码非常短小精悍：
LL ksm(LL a,LL b,LL p)
{
LL res=1;
while(b)
{
if(b&1)res=res*a%p;
a=a*a%p;
b>>=1;
}
return res;
}
类似的，还有快速乘取模，也是一样运用2进制倍增的思想。
例：求\(a*b\%p\)的结果，其中a,b,p都在长整型范围以内。保证p的两倍不超过长整型。
分析：如果两个大整数相乘，结果取幂，虽然结果在整型范围以内，但是中间结果可能超过长整型。
所以可以把乘法操作用快速加来代替。
void ksc(LL a,LL b,LL p) //求a*b%p
{
LL res=0; //此处要初始化为0.
a%=p;
while(b)
{
if(b&1)res=(res+a)%p;
a=(a+a)%p;
b>>=1;
}
return res;
}
有了快速乘以后，对于长整型的幂取模，我们也可以解决了。使用快速加来代替乘法即可。
例：求\(a^b\%p\)的结果，\(a,b,p\)的值不超过\(10^16\).
分析:如果使用乘法，则两个长整型数相乘，结果可能超出长整型。所以可以用快速加来代替乘法。
代码如下：
#include<bits/stdc++.h>
using namespace std;
#define LL long long int
LL a,b,c;
LL ans;
LL quickmul(LL a,LL b,LL p)
{
LL res=0;
while(b)
{
if(b&1)res=(res+a)%p;
a=(a+a)%p;
b>>=1;
}
return res;
}
LL quickpow(LL a,LL b,LL p)
{
LL res=1;
while(b)
{
if(b&1)res=quickmul(res,a,p);
a=quickmul(a,a,p);
b>>=1;
}
return res;
}
int main()
{
while(scanf("%lld%llld%lld",&a,&b,&c)!=-1,a||b||c)
{
ans=quickpow(a,b,c);
printf("%lld\n",ans);
}
return 0;
}
同余等价类、剩余系、缩系
在模意义下，每一个数表示的是一个同余等价类。比如在模3的意义下，1表示的是所有除以3余数为1的整数。所以\(0\)表示的是\(\{0,3,6,9,\dots,\}\) ,\(1\)表示的是\(\{1,4,7,\dots,\}\),而\(2\)表示的是\(\{2,5,8,\dots,\}\).
可以发现，在模意义下，所有的非负整数可以被分为若干同余等价类。
模p的剩余系指的是\(\{0,1,2,\dots,p-1\}\),即小于p的所有非负整数，这个集合中包含了所有模p的余数。集合中的每一个数其实代表一个同余等价类。这个集合称为模p的剩余系，记为\(Z_p\)
如p为6，则\(Z_p=\{0,1,2,3,4,5\}\).
在模p的剩余系中，与p互质的数的集合，称为模p的缩系,记为\(Z_p^*\).
如p为6，则\(Z_p^*=\{1,5\}\).
缩系又称为简化剩余系。
质数及其筛法
如果一个自然数a，只能被1和它自身整除，且不能被其他数整除，则a为质数。换句话说，质数只有两个因子。
注意：1不是质数。
100以内的质数有:2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,91,97
例1:判断一个数n是否为质数。
可以枚举\([2,n-1]\)，如果n在这个范围内有因子，则n不是质数。这个方法的时间复杂度为\(O(n)\)
更好的方法是枚举\([2 ,\sqrt(n)]\)，时间复杂度为\(O(sqrt(n))。\)因为n的因子是成对出现的，如果\(i | n\),则\(\frac n i |n\)，而\(i,\frac n i\) 必有至少一个在\([2,\sqrt(n)]\)内。
参考代码如下：
#include<bits/stdc++.h>
using namespace std;
int n;
int main()
{
scanf("%d",&n);
for(int i=2;i<=sqrt(n);i++)
if(n%i==0)
{
printf("n is not a prime\n");
return 0;
}
printf("n is a prime\n");
return 0;
}
例2：请输出\([1,n]\)内的所有质数。每行一个质数。n<=100000。
可以利用例1中判断质数的方法，将1~n的所有数都判断一下，如果是质数，就输出它。这种方法的时间复杂度为\(O(n\sqrt n)\)
代码如下：
#include<iostream>
#include<cmath>
using namespace std;
int n;
bool isprime(int n)
{
for(int i=2;i*i<=n;i++)
{
if(n%i==0)return 0;
}
return 1;
}
int main()
{
scanf("%d",&n);
for(int i=2;i<=n;i++) //将i=1排除掉。
if(isprime(i))printf("%d\n",i);
return 0;
}
这个方法还是太慢了。毕竟时间复杂度是\(O(N\sqrt N)\)。
下面介绍两种更快的质数筛法。
埃氏筛法.
时间复杂度为\(O(NlogN)\)
算法是这样的：先把n以内的2的倍数(不包含2)全部删除，再把n以内的3的倍数(不包含3)全部删除，……，这样做下去，最后剩下的数就全部是质数。
这里的删除其实不算真的删除，只是打上一个标记而已。
每个数都会被它的质因子打一次标记，而一个数的质因子个数不超过log,所以，这个时间复杂度为\(O(NlogN)\)
代码参考如下：
#include<iostream>
using namespace std;
#define MAXN 1000005
bool flgs[MAXN];//flgs[i]==0表示i为质数 不考虑i<2的情况
int n,cnt,primes[MAXN];
void getprime(int n)
{
for(int i=2;i<=n;i++)
{
if(flgs[i]==1)continue;
primes[cnt++]=i;
for(int j=i;j<=n/i;j++)
flgs[j*i]=1;
}
}
int main()
{
scanf("%d",&n);
getprime(n);
for(int i=0;i<cnt;i++)
printf("%d\n",primes[i]);
return 0;
}
欧式筛法
在上述的埃氏筛法中，一个数可能被它的各个质因子都筛了一遍，而一个数在质因子是不会超过\(logN\)的,所以时间复杂度为\(O(NlogN)\)。而欧式筛法保证每个数只被它的最小质因子筛一遍，这样，时间复杂度为\(O(N)\)。
欧式筛法的过程是这样的：
有一个质数集合Q，一开始，质数集合为空,将最小的质数2加入集合。此时\(Q=\{2\}\)
从2开始枚举依次枚举倍数a，如果倍数a为质数，即将a加入集合Q。
将质数集合中的每个质数的a倍都删掉，因为这些数都是合数。删除指的是打上删除标记，不需要真的删除掉。显然，未打删除标记的数都是质数了。
比如当前质数集合为\(Q=\{p_1,p_2,\dots,p_k\}\),当前枚举到的倍数为\(a\),如果a没有删除标记，则a为质数，先将a加入集合Q。接下来枚举质数集合中的所有质数，则将\(a*p_1,a*p_2,\dots,a*p_k\)全部打上删除标记。
在当前倍数为a，枚举质数的过程中，若发现a是某个质数\(p_i\)的倍数时，此时后续的质数就无需再枚举了，可以提前退出。为什么呢？设后续的质数为\(p_{i'},(p_{i'}>p_i)\),那么那些质数的a倍，其实也是质数\(p_i\)的倍数。即\(a*p_{i'}=a'*p_i\),因为p_i<p_{i'},所以，\(a'>a\)。因为我们要保证每一个数是被它的最小质因子给删掉，\(a*p_{i'}\)应该被\(p_i\)删掉，所以就留给倍数a变成\(a'\)再去处理了。
例：求给一个整数n，输出区间[1,n]中的所有质数。
n<=1000000
参考代码如下：
#include<bits/stdc++.h>
using namespace std;
#define MAXN 1000005
int prime[MAXN];
bool delflg[MAXN];
int n,cnt;
void sieve(int x)
{
for(int i=2;i<=x;i++)
{
if(delflg[i]==0)prime[cnt++]=i;
for(int j=0;j<cnt;j++)
{
if(prime[j]*i>x)break;
delflg[prime[j]*i]=1;
if(i%prime[j]==0)break;
}
}
}
int main()
{
scanf("%d",&n);
sieve(n);
for(int i=0;i<cnt;i++)
{
printf("%d\n",prime[i]);
}
return 0;
}
例1：给出一个整数\(n\),求\(n^n\)的最低位的值。\(n<=10^{16}\)
分析：本题即是求\(n^n\)模10的结果。
HDU1061
参考代码略。
例2：给出一个整数n，求有多少个整数k，满足k^k<n.
\(1\leq n\leq 10^{18}\)
分析：观察到k其实很小。所以可以直接枚举k，求k的快速幂。
但是因为\(n\)很大，所以k的快速幂可能会溢出。
注意：溢出不一定为负数。
本题可以打表，可以发现k不大于15.
如果不大表，一定要防止溢出。因为溢出过后，实际值一定是超过了n的。但是如何发现它溢出了呢？
可以再用除法检测一下。
参考代码：
#include<bits/stdc++.h>
using namespace std;
#define LL long long int
LL qpow(LL a,LL b)
{
LL res=1,tmp;
while(b)
{
if(b&1)
{
tmp=res;
res=res*a;
if(a==0||res/a!=tmp)return 0;
}
tmp=a;
a=a*a;
if(a/tmp!=tmp)return 0;
b>>=1;
}
return res;
}
LL n;
int main()
{
while(cin>>n)
{
LL lz=1;
while(1)
{
LL res=qpow(lz,lz);
if(res>n||res<=0)break;
lz++;
}
cout<<lz-1<<endl;
}
return 0;
}
例3：有一个n的排列，如果存在一个位置\(1<=i<=n\),使得前\(i\)个数单调递增或递减，第\(i+1\)个数到第\(n\)个数也单调递增或递减。则称它为漂亮的。问漂亮的n排列有多少个？你只需要输出答案模\(p\)的结果。\(1\leq n,p\leq10^{18}\)
source:hdu5187
分析：可以找规律，得到递推式为\(ans=2^n-2\).
这个式子其实也可以推出来。
先考虑从大到小依次选择这n个数来摆放。从第二个数开始，它只能放在已放好的那些数的左边或右边，这样形成的是一个单峰的排列。一共有\(2^{n-1}\)种。这里边包含了所有元素单调递增或所有元素单调递减的排列。
单谷排列也类似的，只需要由小到大选择这n个数即可。摆放的方法是一样的。这又有\(2^{n-1}\)种。然后其中也包含了单调递增或单调递减这两种方案，这两种在单峰排列中已经算过了。
所以答案\(ans=2*2^{n-1}-2=2^n-2\)
考虑一下特殊的几种情况：比如p=1或n=1等等。
然后因为\(n,p\)都是长整型，所以快速幂可能爆长整型，要用快速加来代替乘法。
参考代码：
#include<iostream>
#include<cstdio>
#include<cstring>
#include<cstdlib>
using namespace std;
#define LL long long int
LL n,p;
LL qmul(LL a,LL b,LL p)
{
LL res=0;
while(b)
{
if(b&1)res=(res+a)%p;
a=(a+a)%p;
b>>=1;
}
return res;
}
LL qpow(LL n,LL p)
{
LL res=1,a=2;
while(n)
{
if(n&1)res=qmul(res,a,p);
a=qmul(a,a,p);
n>>=1;
}
return res;
}
int main()
{
while(scanf("%lld%lld",&n,&p)!=-1)
{
if(p==1)printf("0\n");
else if(n==1)printf("%d\n",n);
else printf("%lld\n",(qpow(n,p)+p-2)%p);
}
return 0;
}
例4：给两个整数n和k，令\(f(i)=i^k\),求\(f(1)+f(2)+\dots+f(n)\)的值。这个值可能很大，你只需要输出结果模\(10^9+7\)的值。最多有T组数据。
\(1<=n<=10000,0<=k<=5,T<=20\)
分析：k很小，直接模拟即可。
例5：给出n对整数a和b，和一个模数m，求\((a_1^{b_1}+a_2^{b_2}+\dots+a_k^{b_k})\%m\)的值。
poj1995
最大公约数与最小公倍数
一般的，设\(a_1,a_2,a_3,\dots，a_k\)是\(k\)个正整数，如果存在一个正整数\(d\)，使得\(d| a_1,d\mid a_2,\dots d|a_k\),那么\(d\)则为\(a_1,a_2,\dots,a_k\)的公约数。在所有公约数中，最大的称为最大公约数，记为\(gcd(a_1,a_2,\dots,a_k)\)，显然它是存在的，至少为1。当\(gcd=1\)时，称这\(n\)个数是互质的或既约的。公约数一定是最大公约数的约数。
一般的，设\(a_1,a_2,a_3,……a_k\)是\(k\)个正整数，如果存在一个正整数\(d\)，使得\(a_1|d,a_2|d,a_3|d,……a_k|d\),那么\(d\)则为\(a_1,a_2,……a_k\)的公倍数，其中最小的称为最小公倍数数，记为\(lcm(a1,a2,……,ak)\)，显然它是存在的。公倍数一定是最小公倍数的倍数。
定理：\(lcm(a,b) \times gcd(a,b)=a \times b\)
证明：将\(a,b\)进行质因子分解，设\(a,b\)的质因子集合并集为\(\{p_1,p_2,\dots,p_n\}\)。
设
\(a=p_1^{k_1}p_2^{k_2}\dots p_n^{k_n}\) 其中\(0\leq k_1,0\leq k_2,0\leq k_3,\dots 0\leq k_n\)
\(b=p_1^{j_1}p_2^{j_2}\dots p_n^{j_n}\) 其中\(0\leq j_1,0\leq j_2, 0\leq j_3,\dots,0\leq j_n\)
那么显然
\(gcd(a,b)=p_1^{min(k_1,j_1)}p_2^{min(k_2,j_2)}\dots p_n^{min(k_n,j_n)}\)
\(lcm(a,b)=p_1^{max(k_1,j_1)}p_2^{max(k_2,j_2)}\dots p_n^{max(k_n,j_n)}\)
\(\because min(k_i,j_i)+max(k_i,j_i)=k_i+j_i\)
\(\therefore gcd(a,b)*lcm(a,b)=p_1^{k_1+j_1}*p_2^{k_2+j_2}\dots p_n^{k_n+j_n}=a*b\)
辗转相除法
\(gcd(a,b)=gcd(b,a\%b)\)
证明：
设\(gcd(a,b)=P,则a=kP,b=gP,gcd(b,a%\)
更相减损术的证明方法和辗转相除法差不多。
中国古代，九章算术中提到了一种求最大公约数的方法，即更相减损术。
如果要对分式\(a/b\)约分，如果a，b都是偶数，可以将a，b都折半；否则，将a和b交替地减去对方，直到最后两者相等，那么这个数可以做为原来a，b的最大公约数。
int gcd(int a,int b)
{
if(a%b==0)return b;
return gcd(b,a%b);
}
裴蜀定理
如果\(a,b\)的最大公约数为\(d\),且d能整除c，则存在正整数x和y，满足ax+by=c。
证明：
令\(a'=a/d,b'=b/d,则a'x+b'y=1\),此时，\(a',b'\)互质。即我们只需要证明：对于互质的两个整数\(a',b'\),必存在一个整数\(x\)，满足\(a'x\%b'=1\) 即可。 (1)
引理一：如果\(a,b\)为正整数，且\(a,b\)互质，则不存在小于\(b\)的正整数k，使得\(0 \equiv k*a \pmod b\)
证明显然。用反证法，如果存在\(0<k<b\),满足\(0 \equiv k*a \pmod b\)。
根据整数唯一分解定理，\(k*a\)中一定包含\(b\)的所有质因子，而\(a\)与\(b\)互质，\(a\)中不可能包含\(b\)中的质因子，所以必然k是b的倍数。而\(0<k<b\),所以，\(k\)不可能包含\(b\)中所有质因子。
推论：如果\(a,b\)为正整数，且\(a,b\)互质，\(0,a,2*a,3*a,\dots,(b-1)*a\)，这些数模b，余数互不相等。
证明：如果存在两个不同的数设为\(i*a,j*a\)，它们模\(b\)的余数相等，则\((i-j)*a\%b\)模b的余数为\(0\)。
而因为\(0 \lt i \lt n,0 \lt j \lt n\),假设\(i>j\),则\(0 \lt (i-j)\lt n\),这与引理一矛盾。
得证。
引理二：如果\(a,b\)互质，则必存在一个整数\(k\),满足\(k*a\%b=1\).
证明：根据推论，这些数模b的结果只能是集合\(\{0,1,\dots\,b-1\}\)中的数，而这些结果互不相等，且有b个，所以其中必有一个为1。
根据引理2，如果\(a,b\)互质，则必存在一个整数\(k\),满足\(k*a\%b=1\)，即\(k*a-p*b=1\)
于是，裴蜀定理得证。
威尔逊定理、费马定理、欧拉函数、欧拉定理、逆元、exgcd
威尔逊定理：\((p-1)! \equiv -1 \pmod p\)，当且仅当\(p\)为质数。
其中： \((p-1)!\)表示\(p-1\)的阶乘,即\(1*2*3*\dots*(p-1)\)。
证明：
先证充分性： 即“\(p\)为质数”\(\rightarrow\) \((p-1)! \equiv -1 \pmod p\)。
设\(0\lt i \lt p\),因为\(p\)为质数，可知\((i,p)=1\)。根据裴蜀定理中的引理二，必存在一个整数\(j(0 \lt j \lt p)\),使得\(i*j\%p=1\),称\(j\)是\(i\)的逆元，显然\(i\)的逆元具有唯一性。
而\(i\)的逆元是有可能等于\(i\)的。在哪些情况下，\(i\)的逆元等于\(i\)呢？
令\(i*i\%p=1\)，我们来求\(i\) 的值。
由$$ii%p=1$$
可得：$$ii-k*p=1$$
移项可得：
即
因为\(p\)为质数，所以：要么\((i+1)\)是\(p\)的倍数，要么\((i-1)\)是\(p\)的倍数。因为\(0 \lt i \lt p\),所以必有：
即$$\begin{cases} \begin{aligned} i_1&=p-1 \ i_2&=1 \end{aligned} \end{cases}$$
可知：只有当\(i\)为\(1\)或\(p-1\)时，\(i\)的逆元是它自身。
换句话说，在区间(1,p-1)中，每一个数都有一个不等于自己的唯一的逆元，并且逆元也在这个区间中。即这个区间中的数，可以划分成若干逆元对。
所以，\((p-1)!\%p=1*2*3*\dots*(p-1)\%p=1*(p-1)\%p=p-1\).
而p-1和-1是模p等价的。所以充分性得证。
再证必要性：即 \((p-1)! \equiv -1 \pmod p\) \(\rightarrow\) "\(p\)是质数"。
即证“p不是质数” \(\rightarrow\) \((p-1)!\not\equiv p-1\pmod p\).
如果p不是质数，则[2,p-1]中必有p的因子，设其为\(i\)，则\(i\)和\(p/i\)都是\(p\)的因子。
1.如果\(i\neq p/i\),则\(1*2*3*\dots*(p-1)\)一定是p的倍数，所以它模\(p\)一定为0.
2.如果\(i=p/i\),则则\(1*2*3*\dots*(p-1)\)必为\(i\)的倍数，所以它模\(p\)必为\(i\)的倍数。而因为\(p\)是\(i\)的倍数，且\(i>1\),所以(p-1)不可能是\(i\)的倍数。所以\((p-1)! \%p \neq (p-1)\)
得证。
费马定理: 如果\(p\)为质数，且\(a\%p\neq 0\),则有\(a^{p-1}\%p=1\)
证明：
又因为：
a与p互质，根据裴蜀定理中的推论，\((a*1)\%p,(a*2)\%p,\dots,(a*(p-1))\%p\)的结果均互不相等，且在区间[1,p-1]之间，所以：
由(1),(2)可知：
因为\(p\)是质数，所以\((p-1)!\)与\(p\)互质，所以等式(3)两边同时除以\((p-1)!\),得到：
所以，费马定理得证。
欧拉函数
对一个正整数n，它的欧拉函数记为\(\phi(n)\)，表示不超过n且与n互质的正整数的个数。
\(\phi(1)=1\)
-
如果\(n\)为质数，则\(\phi(n)=n-1\).
-
如果\(n=a^p\)，其中\(a\)为质数，则\(\phi(n)=a^p-a^{p-1}=a^{p-1}*(a-1)=a^p*(1-1/a)\)
-
令\(n=a_1^{p_1}a_2^{p_2}\dots a_k^{p_k}\)，根据积性函数的性质，
\(\begin{aligned}\phi(n)&= \phi(a_1^{p_1})*\phi(a_2^{p_2})*\dots*\phi(a_k^{p_k}) \\ &=a_1^{p_1}*(1-1/a_1)*a_2^{p_2}*(1-1/a_2)*\dots*a_k^{p_k}*(1-1/a_k)\\ &=n*(1-1/a_1)*(1-1/a_2)*\dots*(1-1/a_k) \end{aligned}\)。
例1：给出一个整数n，求小于n且与n互质的数的个数。
n<=10000000
分析：本题即为求n的欧拉函数。直接通过上述公式来求。
#include<bits/stdc++.h>
using namespace std;
int n;
int getphi(int x)
{
int res=x;
for(int i=2;i*i<=x;i++)
{
if(x%i==0)
{
res=1ll*res*(i-1)/i; //此处要用1ll。否则有可能超出整型范围。
while(x%i==0)
x=x/i;
}
}
if(x>1) //最后要注意判断x是否大于1。因为x中可能包含一个大于其平方根的质因数。
res=1ll*res*(x-1)/x;
return res;
}
int main()
{
scanf("%d",&n);
printf("%d\n",getphi(n));
return 0;
}
例2.给出一个整数n，求区间[1,n]中所有数的欧拉函数值。
n<=10000000
用例1的方法求欧拉函数，时间复杂度为\(O(\sqrt n)\).如果本题也用这种方法，肯定超时。
我们可以采用类似于质数的筛法来求欧拉函数。
算法时间复杂度为\(O(nlogn)\)
#include<bits/stdc++.h>
using namespace std;
#define MAXN 1000005
int phi[MAXN];
int n,a,b;
void getphi(int x)
{
for(int i=2;i<=x;i++)
{
if(phi[i]==0) //i一定是质数。
{
phi[i]=i-1;
for(int k=2*i;k<=x;k+=i)
{
if(phi[k]==0)phi[k]=k;
phi[k]=phi[k]/i*(i-1); //每个数都被它的质因子筛一下。
}
}
}
}
int main()
{
scanf("%d%d",&a,&b);
if(a>b)swap(a,b);
getphi(b);
for(int i=a;i<=b;i++)
printf("%d\n",phi[i]);
return 0;
}
欧拉定理
如果\(a\)与\(p\)互质，则\(a^{\phi(p)}\%p=1\)
其中\(\phi(p)\)表示\(p\)的欧拉函数。
证明：
设\(p\)的简化剩余系为\(\{p_1,p_2,\dots,p_k\}\),因为\(a\)与\(p\)互质，所以\(\{a*p_1\%p,a*p_2\%p,\dots,a*p_k\%p\}\)也构成了\(p\)的简化剩余系。
所以:$$(ap_1)(ap_2)\dots(ap_k) \equiv (p_1p_2\dots*p_k) \pmod p$$
所以$$a^\phi(p) \equiv 1 \pmod p$$
逆元
如果整数a，b满足\(a*b\%p=1\),则称\(a,b\)在模p意义下互为逆元。
只有当a与p互质时，在模p意义下a才有逆元，\(a\)的逆元记作\(a^{-1}\).
逆元的性质：若b与p互质，则\(a/b\%p=a*b^{-1}\%p\)
由此可见，通过逆元，可以将除法取余变为乘法取余。
注意此处的“除法”，不是整除，可以理解为分式。
在模意义下，分式取余是一个整数,比如\(2/3\%5\)，它是有意义的。因为在模意义下，每一个数表示的是一个同余类。在模5的意义下，2表示的是所有除以5余数为2的整数。所以2表示的是\(\{2,7,12,\dots,\}\)，于是\(2/3\%5\)和\(12/3\%5\)是一样的，结果为4.
这种分式取模，我们既可以采用同余等价类来解决，也可以采用逆元来求。采用逆元是更一般的方法。
3在模5意义下的逆元为7，因为\(3*7\%5=1\)。
于是有：\(2/3\%5=2*3^{-1}\%5=2*7\%5=4\).
区间逆元
求区间[1,n]的逆元，如果采用上述求逆元的方法做，时间复杂度为\(O(nlogn)\)
有一种递推的方法可以求区间的逆元，时间复杂度为\(O(n)\)。
设模数为m,\(f[n]\)表示\(n\)的逆元，其中1~(m-1)的逆元已经求出，则\(f[n]=(-f[m\%n]*(m/n)\%m+m)\%m\).
这个原理是什么呢？
1.设\(f[i]\)表示\(i\)在模m下的逆元，则有 \(f[i]=-f[m-i]\)。
证明：因为$$f[i]*i%m=1$$,
所以：
即:$$-f[i]*(m-i)%m=1$$
所以：$$f[m-i]=-f[i]$$
2. \(f[i]=k*f[k*i]\)
证明：
因为
即$$(f[ki]k)i=1$$
即$$f[i]=kf[k*i]$$
扩展欧几里得
扩展欧几里得是用来求不定方程如：\(ax+by=c\)。其中，\(a,b,c\)为已知整数，其中\(gcd(a,b) \mid c\)。这样，不定方程才有解。
扩展欧几里得是用递归的思想来求解的。
我们先求方程\(ax+by=gcd(a,b)。\)求出来该方程的解，则原方程的解乘上系数\((c/gcd(a,b))\)就可以了。
因为：$$ax+by=gcd(a,b)=gcd(b,a%b)=bx'+(a%b)y'$$
将a,b看做变量，移项并合并同类项,可得：
只需要求出\(x',y'\)，则可以求出\(x,y\)。
而求\(x',y'\)，可以继续递归。这样递归下去，\(gcd(a,b)\)中的参数b最终会变为0，此时gcd(a,0)=a.
于是有:\(ax+by=gcd(a,0)=a\)，可以求出\(x=1,y=0\).这是最底层的\(x,y\),然后一层层返回，就可以求出原始的\(x\)和\(y\)了。
void exgcd(int a,int b,int &d,int &x,int &y)
{
if(!b)d=a,x=1,y=0;
else
{
exgcd(b,a%b,d,y,x);
y-=x*(a/b);
}
}
中国剩余定理
类似于这样的问题：有一个数，模5余2，模7余3，模13余8，……，求这个数是多少？
实质上就是求一个模线性方程组。
\(\begin{cases} x \equiv a_1 \pmod {r_1} \\ x\equiv a_2 \pmod {r_2} \\ x\equiv a_3 \pmod{r_3} \\ \dots \\ x \equiv a_k \pmod{r_k} \end{cases}\)
其中，\(r_1,r_2,\dots,r_k\)互质。
中国剩余定理就是解决这种模数互质的模线性方程组的。
因为\(r_1,r_2,\dots,r_k\)互质，所以，对于\(r_i\),设\(A_i=\prod_{ j\neq i}r_j\),则\(A_i,r_i\)是互质的，则必存在一个整数\(c_i\),使得\(c_i*A_i \%r_i=1\),则设\(x_i=a_i*c_i*A_i\),则\(x_i\%r_i=a_i\)。因为\(A_i\)是除了\(r_i\)以外的其他\(r\)值的最小公倍数，所以，\(x_i\) 模其他的r，余数都为0，只有模\(r_i\)时，余数为\(a_i\)。
这就是说，把\(x_i\)加到满足其他方程（比如第j个方程,\(j\neq i\)）的解\(x_j\)上，则\((x_i+x_j)\) 也满足方程\(i\)。同理，如果\(x_j\) 也是这么求出来的，则\((x_i+x_j)\)也满足方程\(i\)。
那么，我们对于每一个方程，都按照这个方法求出来该方程的解。把这些解累加起来，你会发现，这个和仍然满足每一个方程。
则\(x=\sum\limits _{i=1}^k a_i*c_i*A_i\)
对于x，如果我们加上所有的r值的最小公倍数，它仍然满足所有方程。
所以，只有存在x，则意味着有无穷多组解。
通解为：
\(x=\sum\limits _{i=1}^k a_i*c_i*A_i+p*\prod_{i=1}^kr_i\)
其中，p是任意整数。
扩展中国剩余定理
中国剩余定理采用的是构造法做的，要求模数互质。但是如果模数不互质，怎么办呢？
那么这就是扩展中国剩余定理要解决的问题了。
\(\begin{cases} x \equiv a_1 \pmod {r_1} \\ x\equiv a_2 \pmod {r_2} \\ x\equiv a_3 \pmod{r_3} \\ \dots \\ x \equiv a_k \pmod{r_k} \end{cases}\)
首先取方程组的前两个方程：
\(x=k*r_1+a_1=p*r_2+a_2\)
\(k*r_1-p*r_2=(a_2-a_1)\)
这个形如\(ax+by=c\)的不定方程，用扩展欧几里得算法即可。
当\(gcd(r_1,r_2)\)不能整除\(a_2-a_1\)时，方程组无解。
否则，根据exgcd，求出\(k_0\)，满足\(k_0*r_1-p*r_2=(a_2-a_1)\)。k的通解为：\(k=k_0+b*(r_2/gcd(r_1,r_2))\)
代入到方程组\(x=k*r_1+a_1\)中，则有\(x=(k_0+b*(r_2/gcd(r_1,r_2)))*r_1+a_1=k_0*r_1+b*lcm(r_1,r_2)+a_1\)
其中b为任意整数。
即：\(x \equiv(k_0*r_1)+a_1 \pmod{lcm(r_1,r_2)}\)
这个方程相当于合并了原来的两个方程，而形式上又和原来的方程一致。
继续采用这个方法，不断地合并方程组中的两个方程，直到最后合并为一个方程，则得到了x的通解了。
注意在这个过程中，要注意求出的\(k_0\)是\(k*r_1-p*r_2=(a_2-a_1)\)的解，而不是\(k*r_1-p*r_2=gcd(r_1,r_2)\)的解。前者是后者的倍数。
一道例题：poj2891 （扩展中国剩余定理模板题）
给出一个模线性方程组，模数不互质，求方程的最小正整数解。如果无解，输出-1.
// ax+by=c
// x%a1=r1
// x%a2=r2
// x=k*a1+r1=p*a2+r2
// k*a1+q*a2=r2-r1
// a1,a2,r2-r1已知，所以可以通过exgcd求出k和q。
// 求出一个k值，设为k0.
// k=k0+a2/gcd(a1,a2)
// x=k*a1+r1=(k0+b*a2/(gcd(a1,a2)))*a1+r1=k0*a1+b*lcm(a1,a2)+r1
// x的通解 x=lcm(a1,a2)*b+k0*a1+r1
#include<iostream>
#include<cstdlib>
#include<cstdio>
#include<cmath>
using namespace std;
#define LL long long int
#define abs(a) (a>0?a:-(a))
LL m,ans;
#define MAXN 1005
LL arr[MAXN],r[MAXN];
LL a,b,d,x,y,lcm;
void exgcd(LL a,LL b,LL &d,LL &x,LL &y)
{
if(!b)
{
d=a,x=1,y=0;
}
else
{
exgcd(b,a%b,d,y,x);
y-=x*(a/b);
}
}
bool noans;
int main()
{
while(scanf("%d",&m)!=-1)
{
noans=0;
for(int i=1;i<=m;i++)
{
scanf("%lld%lld",&arr[i],&r[i]);
}
for(int i=1;i<m;i++)
{
exgcd(arr[i],arr[i+1],d,x,y);
if(abs(r[i+1]-r[i])%d)
{
noans=1;
break;
}
x=(x*((r[i+1]-r[i])/d)%(arr[i+1]/d)+arr[i+1]/d)%(arr[i+1]/d);
lcm=arr[i]*arr[i+1]/d;
arr[i+1]=lcm,r[i+1]=(x*arr[i]+r[i])%lcm;
}
if(noans==0)
ans=r[m];
else
ans=-1;
printf("%lld\n",ans);
}
return 0;
}
miller_rabbin算法
这是一个质数判断算法，它采用随机算法能够在很短的时间里判断一个数是否是质数。
首先，根据费马定理，若\(p\)为质数，取a<p,则\(a^{p-1}\%p=1\).但反过来则不成立。
因为有一类数——卡迈克尔数，它不是质数，但也满足\(a^{p-1}\%p=1\)。
如何将卡迈克尔数甄别出来。我们要用到二次探测方法。
设\(p\)为我们要判断的数，将\(p-1\)分解为\(2^r*k\)，则有\(a^{p-1}=(a^{k})^{2^r}\)
则我们可以先求出\(a^k\),然后将其不断平方，并模p。如果第一次出现模\(p\)的值为1，那么检测上一次的值是否为\(-1\),如果是,则没有毛病；否则可以断定不是质数。
证明:
令$$y^2%p=1$$
则等价于：$$y^2-1=kp$$
等价于：$$(y-1)(y+1)=k*p$$
如果\(p\)是质数，则$$\begin{cases} y-1=0 \ y+1=p \ \end{cases} $$
即y的值只能为\(1\)或\(p-1\)。
所以，如果\(y\neq 1\)且\(y\neq p-1\),则可以断定\(p\)不是质数。
而如果做r次平方，模\(p\)的值都不为1，则\(p\)显然不是质数了。
那如果在平方r次的过程中，最后的结果为1了。而且第一次结果为1时，之前的结果为-1。是不是就能保证\(p\)是质数呢?也不能百分百地肯定。但是相比于仅通过费马定理来判断，它更严格一些了。这种方法称为二次检测。
于是，我们可以多次选择\(a\)，如果选了许多\(a\)来做上述的操作，都不能判断出\(p\)是合数，则\(p\)多半就是一个质数了。
miller_rabbin算法是一个随机算法，并不能完全保证结果准确，但是出错的概率非常小。我们取a的次数越多，出错的概率越小。一般情况下，取20~30次a，正确率已经可以达到99.9999999%了。
FZU1649
给出若干整数，判断它是否是素数。
参考代码如下：
#include<iostream>
#include<cstring>
#include<ctime>
#include<cstdio>
#include<cstdlib>
#include<algorithm>
#include<cmath>
using namespace std;
#define LL unsigned long long int
LL n,k;
int t,r;
bool notprime;
LL ksc(LL a,LL b,LL p)
{
LL res=0;
while(b)
{
if(b&1)res=res+a;
a=a+a;
if(res>=p)res-=p;
if(a>=p)a-=p;
b>>=1;
}
return res;
}
LL ksm(LL a,LL b,LL p)
{
LL res=1;
while(b)
{
if(b&1)res=ksc(res,a,p);
a=ksc(a,a,p);
b>>=1;
}
return res;
}
int main()
{
srand(time(NULL));
while(scanf("%I64d",&n)!=-1)
{
if(n==2||n==3)
{
printf("It is a prime number.\n");
continue;
}
else
{
if(n==1||n&1==0)
{
printf("It is not a prime number.\n");
continue;
}
}
r=0;
k=n-1;
while(!(k&1))
{
r++;
k>>=1;
}
int t=20;
LL res=1,lastres=0;
notprime=0;
for(int id=0;id<=10;id++)
{
res=rand()%(n-1)+1;
res=ksm(res,k,n);
if(res==1)
continue;
lastres=res;
for(int i=0;i<r;i++)
{
res=ksc(res,res,n);
if(res==1)
{
if(lastres!=n-1) notprime=1;
break;
}
lastres=res;
}
if(notprime==1)break;
if(res!=1)
{
notprime=1;
break;
}
}
if(notprime==1)
printf("It is not a prime number.\n");
else
printf("It is a prime number.\n");
}
return 0;
}
高斯消元
高斯消元是解决线性方程组的，它将增广矩阵经过线性变换，使等号左边的系数矩阵变成上三角矩阵，然后从下到上依次求出各个方程的未知数。
步骤不难理解，手算麻烦，正好适合计算机做。
但难点在于如何判断无解和多解。
先看无解的情况。只要上三角矩阵已经化到了最后一列，若此时并不是最后一行，且其下面的行中，等号右边存在不为0 的值，则无解。因为此时，其下面的行中，等号左边的系数必然全部为0了，而等号右边不为0，这是矛盾的。
多解的情况：在化上三角矩阵的过程中，第i列中，从第i行开始，所有元素都为0，则意味着出现了一个自由变量。如果最后不是无解，则必为多解。
其他情况就是唯一解的情况了。
化上三角矩阵的时候，如果当前第i行第i列的系数为0，则要在下面找一个第i列系数不为0的，设为\(i'\),将第\(i'\)与第\(i\)行交换。如果第i列下面的系数均为0，则记录\(x_i\)为自由变量，然后行不变，列加1，继续做。当做到最后一列了，则化三角矩阵的过程结束。此时，查看后面的行中，等号右边是否为0，若不为0，则无解；若为0，则判断是否存在自由变量，若有自由变量，则为多解，解的个数等于自由变量的取值范围的幂。若自由变量的取值范围为\(p\),自由变量有\(m\)个，则解的个数为\(p^m\)；若无自由变量，则仍然为唯一解。
如果高斯消元为模运算，模数为2，则为异或方程组，所有的运算均可以用异或完成。
以下代码为异或方程组。
高斯消元的核心代码：
int gauss()
{
int row,col,tmprow=0;
for(row=1,col=1;col<=n;col++)
{
tmprow=0;
if(arr[row][col]==0) //当前系数为0,
{for(int j=row+1;j<=n;j++) //在下面找一个系数不为0的行,记为tmprow。
{
if(arr[j][col]==1)
{
tmprow=j;
break;
}
}
if(tmprow)for(int j=col;j<=n+1;j++) //将tmprow和row两行交换。
swap(arr[row][j],arr[tmprow][j]);
}
if(arr[row][col]==0){freecnt++;continue;} //表示第row及以下，该列的系数均为0，则出现了一个自由元。
for(int j=row+1;j<=n;j++) // 将第row及以下的行做线性变换，将第col列化为0.
{
if(arr[j][col]==1)
for(int k=col;k<=n+1;k++)
arr[j][k]^=arr[row][k];
}
row++;
}
if(row<=n) //判断无解的情况
{
for(int j=row;j<=n;j++)
if(arr[j][n+1])
return -1;
}
return freecnt; //freecnt>0表示有多解，freecnt==0表示有唯一解。
}
例题1：开关灯 poj1222
【题目大意】有一个5*6的开关灯的阵列，当改变某盏灯的状态时，其上下左右的灯也会变换。现在给出这个开关灯阵列的一个初始状态，请输出一种方案，可以使得这些灯全部关掉。每盏灯只能按最多一次开关。
【分析】：有两种方法做。
- 枚举
枚举第一列的所有操作，其他列的操作就是固定的了，只需要检测一下最后能否关掉所有灯即可。时间复杂度为\(O(2^5*30)\) - 高斯消元
设第\(i\)盏灯的开关次数为\(x_i\),考虑第i盏灯的上下左右四盏灯对i盏灯的影响，即第i盏灯的次数加上其上下左右四盏灯的操作次数的奇偶性等于第i盏灯的初始状态与目标状态异或值。
于是可以得到一个由30个方程组成的异或方程组。
代码如下：
#include<iostream>
#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<algorithm>
using namespace std;
#define MAXN 33
int arr[MAXN][MAXN];
int sou[MAXN][MAXN];
int T,n,m,dd[5][2]={0,1,1,0,0,-1,-1,0};
int ans[31];
bool freex[31];
int x[31];
int gauss()
{
int row=0,col=0;
for(row=1,col=1;col<=30;col++)
{
if(arr[row][col]==0)
for(int j=row+1;j<=30;j++)
{
if(arr[j][col]==1)
{
for(int k=col;k<=31;k++)swap(arr[j][k],arr[row][k]);
break;
}
}
if(arr[row][col]==0){freex[col]=1;continue;}
for(int j=row+1;j<=30;j++)
{
if(arr[j][col]==0)continue;
for(int k=col;k<=31;k++)
{
arr[j][k]^=arr[row][k];
}
}
row++;
}
if(row<col)
{
for(int i=row;i<=30;i++)if(arr[i][31])return -1; // 无解
}
for(int i=row-1,j=30;j>0;j--)
{
if(freex[j]==1)
{
x[j]=0;
continue;
}
x[j]=arr[i][31];
for(int k=j+1;k<=30;k++)
if(arr[i][k]&&x[k])
x[j]^=arr[i][k];
i--;
}
return 31-row;
}
int main()
{
scanf("%d",&T);
for(int t=1;t<=T;t++)
{
memset(x,0,sizeof x);
memset(freex,0,sizeof freex);
memset(arr,0,sizeof arr);
//memset(sou,0,sizeof sou);
for(int i=1;i<=5;i++)
for(int j=1;j<=6;j++)
{
scanf("%d",&sou[i][j]);
}
int id=0;
for(int i=1;i<=5;i++)
for(int j=1;j<=6;j++)
{
id++;
for(int k=0;k<4;k++)
{
int xx=i+dd[k][0],yy=j+dd[k][1];
if(xx>0&&xx<6&&yy>0&&yy<7)
{
arr[id][(xx-1)*6+yy]=1;
}
}
arr[id][(i-1)*6+j]=1;
arr[id][31]=sou[i][j];
}
int v=gauss();
printf("PUZZLE #%d\n",t);
if(v>=0)
{
for(int i=1;i<=5;i++)
{
for(int j=1;j<6;j++)
printf("%d ",x[(i-1)*6+j]);
printf("%d\n",x[i*6]);
}
}
else printf("no solution\n");
}
return 0;
}
例题2：开关问题
【题目大意】
有N个相同的开关，每个开关都与某些开关有着联系，每当你打开或者关闭某个开关的时候，其他的与此开关相关联的开关也会相应地发生变化，即这些相联系的开关的状态如果原来为开就变为关，如果为关就变为开。你的目标是经过若干次开关操作后使得最后N个开关达到一个特定的状态。对于任意一个开关，最多只能进行一次开关操作。你的任务是，计算有多少种可以达到指定状态的方法。（不计开关操作的顺序）
【分析】这个开关阵列不像上题那样是一个规则的矩形，它形成了一个有向图。因为图的复杂性和无序性，不能再用枚举的方法了。只能使用高斯消元来做。对于一个开关，将能够影响它的开关的操作次数相加，其奇偶性必然等于该开关的初始状态与结束状态的异或值的奇偶性。
于是可以得到一个方程组，即可使用高斯消元来做了。因为只需要求有多少种方法，所以先判断是否有解，若有解，则只需要统计自由变量的个数，设为cnt，则答案为\(2^{cnt}\)。
代码如下：
#include<iostream>
#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<algorithm>
using namespace std;
#define MAXN 33
int arr[MAXN][MAXN];
int n,T,freecnt,ans,st[MAXN],ed[MAXN];
bool flg[MAXN];
int gauss()
{
int row,col,tmprow=0;
for(row=1,col=1;col<=n;col++)
{
tmprow=0;
if(arr[row][col]==0)
{for(int j=row+1;j<=n;j++)
{
if(arr[j][col]==1)
{
tmprow=j;
break;
}
}
if(tmprow)for(int j=col;j<=n+1;j++)
swap(arr[row][j],arr[tmprow][j]);
}
if(arr[row][col]==0){freecnt++;continue;}
for(int j=row+1;j<=n;j++)
{
if(arr[j][col]==1)
for(int k=col;k<=n+1;k++)
arr[j][k]^=arr[row][k];
}
row++;
}
if(row<=n)
{
for(int j=row;j<=n;j++)
if(arr[j][n+1])
return -1;
}
return freecnt;
}
int main()
{
int a,b;
scanf("%d",&T);
while(T--)
{
scanf("%d",&n);
freecnt=0;
memset(arr,0,sizeof arr);
for(int i=1;i<=n;i++)scanf("%d",&st[i]);
for(int i=1;i<=n;i++)
{
scanf("%d",&ed[i]);
if(st[i]!=ed[i])arr[i][n+1]=1;
arr[i][i]=1;
}
while(scanf("%d%d",&a,&b)!=-1,a||b)
{
arr[b][a]=1;
}
int v=gauss();
if(v==-1)
printf("Oh,it's impossible~!!\n");
else
{
printf("%d\n",(1<<v));
}
}
return 0;
}
例题3：widget factory
【题目大意】给出m条记录，每条记录记载了起始时间为周几，结束时间为周几，但并不知道跨越了多少周，在这段时间里生产了每种小玩具各多少个。小玩具一共有n种，第i种玩具的生产时间为\(t_i,3\leq t_i \leq 9\)。求每种小玩具生产的时间。无解或多解输出其他的信息。
\(n\leq 300,m \leq 300\)
【分析】本题即为模7的模线性方程组，通过高斯消元直接做就行了。但\(t_i\)限制在区间[3,9]之间，这个要怎么处理呢？而区间[3,9]的长度也为7。可以先按照模7运算求解，如果解为0,1,2，给它加上7即可。
线性变换时，可以将两个系数都扩充到他们的最小公倍数，然后再消去。
因为模数为7，回代的过程需要用到逆元。
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
#define MAXN 305
int arr[MAXN][MAXN];
char day1[10],day2[10];
int n,m,k,freecnt,ans[MAXN];
int calcday(char *s)
{
if(strcmp(s,"MON")==0)return 1;
if(strcmp(s,"TUE")==0)return 2;
if(strcmp(s,"WED")==0)return 3;
if(strcmp(s,"THU")==0)return 4;
if(strcmp(s,"FRI")==0)return 5;
if(strcmp(s,"SAT")==0)return 6;
if(strcmp(s,"SUN")==0)return 7;
}
void gcd(int a,int b,int &d,int &x,int &y)
{
if(!b)d=a,x=1,y=0;
else
{
gcd(b,a%b,d,y,x);
y-=x*(a/b);
}
}
int ni(int x)
{
int res=0,d=0,tmp;
gcd(x,7,d,res,tmp);
return res%7;
}
int gauss()
{
int row,col;
freecnt=0;
for(row=1,col=1;col<=n;col++)
{
if(arr[row][col]==0)
{
for(int i=row+1;i<=m;i++)
{
if(arr[i][col])
{
for(int j=col;j<=n+1;j++)
swap(arr[i][j],arr[row][j]);
break;
}
}
}
if(arr[row][col]==0){freecnt++;continue;}
for(int i=row+1;i<=m;i++)
{
if(arr[i][col])
{
int tt1,tt2,tmp;
gcd(arr[i][col],arr[row][col],tmp,tt1,tt2);
int t1=arr[i][col];
for(int j=col;j<=n+1;j++)
{
arr[i][j]=((arr[i][j]*arr[row][col]/tmp-arr[row][j]*t1/tmp)%7+7)%7;
}
}
}
row++;
}
if(row<m+1)
{
for(int i=row;i<=m;i++)
if(arr[i][n+1]%7!=0)return -1;
}
if(freecnt>0)return freecnt;
for(int i=n,j=row-1;i>0;i--)
{
int tmp=arr[j][n+1];
for(int k=i+1;k<=n;k++)
{
if(ans[k]&&arr[j][k])tmp=(tmp-arr[j][k]*ans[k])%7;
}
ans[i]=((tmp*ni(arr[j][i])%7)+7)%7;
j--;
}
return 0;
}
int main()
{
int a;
while(scanf("%d%d",&n,&m),n||m)
{
memset(arr,0,sizeof arr);
memset(ans,0,sizeof ans);
for(int i=1;i<=m;i++)
{
scanf("%d %s %s",&k,day1,day2);
for(int j=1;j<=k;j++)
{
scanf("%d",&a);
arr[i][a]++;
}
arr[i][n+1]=calcday(day2)-calcday(day1)+8;
}
for(int i=1;i<=m;i++)
for(int j=1;j<=n+1;j++)
arr[i][j]=(arr[i][j]%7+7)%7;
int v=gauss();
if(v==-1)
{
printf("Inconsistent data.\n");
}
else if(v>0)
{
printf("Multiple solutions.\n");
}
else
{
for(int i=1;i<=n;i++)
if(ans[i]<3)ans[i]+=7;
for(int i=1;i<n;i++)
{
printf("%d ",ans[i]);
}
printf("%d\n",ans[n]);
}
}
return 0;
}
例题4：poj1681
【题目大意】一块n*n的棋盘，你每次可以对一个格子改变颜色，但是它上下左右四个格子也会改变颜色。只有白色和黄色。现在给出棋盘的初始状态，问你最少要经过多少次操作，才能将棋盘变成全部为黄色。
【分析】用高斯消元。对于自由变量，需要枚举，确定了他们的值以后，再回代到方程，求其他未知数，并找出最小的次数。自由变量可以随意取值，方程总有解。
使用高斯消元时，系数矩阵的行数必须不少于列数。即使方程的个数少于未知数个数，但系数矩阵的行数要留够，并且要列数到了最后一列才结束。
而不能看行数是否到达最后一个方程，这样写的话，在方程数少于变量个数时，会出问题。
另外出现自由变量时，回代过程中，不能从最后一个方程开始回代。而只能从n-freecnt那一行开始回代。
#include<bits/stdc++.h>
using namespace std;
#define MAXN 240
int n,T,N,minans,myfreecnt,cnt,myfreex[MAXN],ans[MAXN],arr[MAXN][MAXN];
bool myfree[MAXN];
char carr[MAXN][MAXN];
int dd[4][2]={0,1,1,0,0,-1,-1,0};
int gauss()
{
int row,col;
for(row=1,col=1;col<=N;col++)
{
if(arr[row][col]==0)
{
for(int i=row+1;i<=N;i++)
{
if(arr[i][col])
{
for(int k=col;k<=N+1;k++)
swap(arr[row][k],arr[i][k]);
break;
}
}
}
if(arr[row][col]==0)
{
myfreex[myfreecnt++]=col;
myfree[col]=1;
continue;
}
for(int i=row+1;i<=N;i++)
{
if(arr[i][col])
{
for(int j=col;j<=N+1;j++)
arr[i][j]^=arr[row][j];
}
}
row++;
}
if(row<N+1)
{
for(int i=row;i<=N;i++)
{
if(arr[i][N+1])
return -1;
}
}
return myfreecnt;
}
void reback()
{
int tot=(1<<myfreecnt);
for(int i=0;i<tot;i++)
{
cnt=0;
for(int j=0;j<myfreecnt;j++)
{
if((i>>j)&1) //tot的右数第j为1.
{
ans[myfreex[j]]=1;
}
else
ans[myfreex[j]]=0;
}
for(int j=N,i=N-myfreecnt;j>0;j--) //这里i的初始值开始写成了N，wa了好久。
{
if(myfree[j]==1)
continue;
ans[j]=arr[i][N+1];
for(int k=j+1;k<=N;k++)
if(arr[i][k]&&ans[k])ans[j]^=arr[i][k];
i--;
}
cnt=0;
for(int j=1;j<=N;j++)
if(ans[j]==1)cnt++;
if(cnt<minans)
minans=cnt;
}
}
int main()
{
scanf("%d",&T);
while(T--)
{
memset(arr,0,sizeof arr);
memset(ans,0,sizeof ans);
memset(myfree,0,sizeof myfree);
memset(myfreex,0,sizeof myfreex);
myfreecnt=0;
minans=0X3F3F3F3F;
scanf("%d",&n);
N=n*n;
for(int i=1;i<=n;i++)
{
scanf("%s",carr[i]+1);
}
int id=0,xi,yi;
for(int i=1;i<=n;i++)
for(int j=1;j<=n;j++)
{
id++;
arr[id][(i-1)*n+j]=1;
for(int k=0;k<4;k++)
{
xi=i+dd[k][0],yi=j+dd[k][1];
if(xi>0&&xi<=n&&yi>0&&yi<=n)
{
arr[id][(xi-1)*n+yi]=1;
}
}
arr[id][N+1]=(carr[i][j]!='y');
}
int v=gauss();
if(v==-1)
printf("inf\n");
else
{
reback();
printf("%d\n",minans);
}
}
return 0;
}
浙公网安备 33010602011771号