Special equations
Time Limit: 2000/1000 MS (Java/Others) Memory Limit: 32768/32768 K (Java/Others)Total Submission(s): 20 Accepted Submission(s): 12
Special Judge
Then comes T lines, each line starts with an integer deg (1<=deg<=4), meaning that f(x)'s degree is deg. Then follows deg integers, representing a n to a 0 (0 < abs(a n) <= 100; abs(a i) <= 10000 when deg >= 3, otherwise abs(a i) <= 100000000, i<n). The last integer is prime pri (pri<=10000).
Remember, your task is to solve f(x) 0 (mod pri*pri)
4 2 1 1 -5 7 1 5 -2995 9929 2 1 -96255532 8930 9811 4 14 5458 7754 4946 -2210 9601
Case #1: No solution! Case #2: 599 Case #3: 96255626 Case #4: No solution!
容易看出,若x0是同余方程
f(x) º 0 (mod pa) (1)
的解,则它必是方程
f(x) º 0 (mod pa- 1) (2)
的解,因此,必有与x0相应的方程(2)的某个解x1,使
x0º x1 (mod pa- 1),x0 = x1+pa- 1t0,
此处,t0是某个适当的整数。
这提示我们:可以从方程(2)的解中去求方程(1)的解。于是,现在的问题是,对于方程(2)的每个解x1,是否必有方程(1)的解x0与之对应?若有,如何去确定它?
定理 设p是素数,a³ 2是整数,f(x) = anxn+L+ a1x + a0是整系数多项式,又设x1是同余方程(2)的一个解。以f¢(x)表示f(x)的导函数。
(ⅰ) 若f¢(x1)0 (mod p),则存在整数t,使得
x = x1+pa- 1t (3)
是同余方程(1)的解。
(ⅱ) 若f¢(x1) º 0 (mod p),并且f(x1)º 0 (mod pa),则对于t = 0,1, 2, L, p - 1,式(3)中的x都是方程(1)的解。
证明 我们来说明,如何确定式(3)中的t,使x1+ pa- 1t满足同余方程(1),即要使
an(x1+pa- 1t)n+an- 1(x1+pa- 1t)n- 1+L+a1(x1+pa- 1t)+a0º 0 (mod pa), (4)
即
f(x1) +pa- 1tf¢(x1)º 0 (mod pa),
tf¢(x1)º-(mod p)。 (5)
下面考虑两种情形。
(ⅰ) 若f¢(x)0 (mod p),则关于t的同余方程(5)有唯一解t º t0 (mod p),即t = t0+ pk(kÎZ),于是
x º x1+ pa- 1t0 (mod pa)
是同余方程(1)的解。
(ⅱ) 若f¢(x1) º 0 (mod p),并且f(x1)º 0 (mod pa),则式(5)对于任意的整数t成立,即同余方程(5)有p个解
t º i(mod p),0 £ i £ p - 1。
于是x º x1+ pa- 1i (mod pa),0 £ i £ p - 1,都是同余方程(1)的解。证毕。
在定理中,没有对f¢(x1)º 0 (mod p)并且 f(x1)0 (mod pa)的情形进行讨论。事实上,此时,同余方程(5)无解。即,我们无法从同余方程(2)的解x1出发去求得同余方程(1)的解。
由定理,可以把解同余方程(1),转化为解同余方程
f(x) º 0 (mod p)。 (6)
事实上,由方程(6)的解,利用定理,可以求出方程f(x) º 0 (mod p2)的解,再利用定理,又可以求出方程f(x)º 0 (mod p3)的解,LL,直到求出方程(1)的解。
#include<stdio.h>
#include<stdlib.h>
#include<vector>
using namespace std;
typedef __int64 II;
int gcd(int a,int b){
return b==0? a:gcd(b,a%b);
}
typedef struct {
int deg,a[15];
}function;
function f,diff;
II op(function f, II i){
II ret=0;
for(int loop=f.deg;loop>=0;loop--){
ret=ret*i+(II)f.a[loop];
}
return ret;
}
II op(function f, II i, II mod){
II ret=0;
for(int loop=f.deg;loop>=0;loop--){
ret=ret*i+f.a[loop];
ret=ret%mod;
}
return ret;
}
II mm(II a,II b,II c){
if(b==0) return 1%c;
II tmp=mm(a,b/2,c);
return tmp*tmp%c*(b%2? a:1)%c;
}
II ext_gcd(II a,II &x,II b,II &y){//a*x+b*y=(x,y)=d
if(b==0){
x=1;y=0;
return a;
}
II xx,yy,d=ext_gcd(b,xx,a%b,yy);
x=yy;
y=xx-a/b*yy;
return d;
}
vector<II> linear(II a,II b,II m){//ax=b(mod m)
a%=m,b%=m;
vector<II> ret;
ret.clear();
II x,y,d=ext_gcd(a,x,m,y);
//if(d<0) d=-d;
if(b%d==0){
II e=(x*(b/d)%m+m)%m;//基础解
//for(int i=0;i<(d<0? -d:d);i++) ret.push_back((e+i*(m/d))%m);
ret.push_back(e);
}
return ret;
}
int main(){
//freopen("in_std.txt","r",stdin);
//freopen("out_std.txt","w",stdout);
int T;
II pri,p2;
scanf("%d",&T);
for(int cas=1;cas<=T;cas++){
scanf("%d",&f.deg);
for(int i=f.deg;i>=0;i--) scanf("%d",&f.a[i]);
scanf("%I64d",&pri);
p2=pri*pri;
II i;
for(i=0;i<pri;i++)
if(op(f,i,pri)==0) break;
if(i>=pri) {printf("Case #%d: No solution!\n",cas);continue;}
diff.deg=f.deg-1;
for(int loop=f.deg-1;loop>=0;loop--){
diff.a[loop]=f.a[loop+1]*(loop+1);//?
}
II fx1=op(diff,i),t;
if(fx1%pri!=0){//fx1与pri互质
II fx=op(f,i);
vector<II> fuck=linear(fx1,-fx/pri,pri);
t=fuck[0];
printf("Case #%d: %I64d\n",cas,((i+t*pri)%p2+p2)%p2);
}else{
II fx=op(f,i);
if(fx%p2==0) printf("Case #%d: %I64d\n",cas,i);
else printf("Case #%d: No solution!\n",cas);
}
}
//while(1);
return 0;
}