poj 3696 &&hdoj 2462 The Luckiest number（同余&欧拉） 好题

The Luckiest number

Time Limit: 1000MS		Memory Limit: 65536K
Total Submissions: 4839		Accepted: 1293

Description

Chinese people think of '8' as the lucky digit. Bob also likes digit '8'. Moreover, Bob has his own lucky numberL. Now he wants to construct his luckiest number which is the minimum among all positive integers that are a multiple ofL and
consist of only digit '8'.

Input

The input consists of multiple test cases. Each test case contains exactly one line containingL(1 ≤
L ≤ 2,000,000,000).

The last test case is followed by a line containing a zero.

Output

For each test case, print a line containing the test case number( beginning with 1) followed by a integer which is the length of Bob's luckiest number. If Bob can't construct his luckiest number, print a zero.

Sample Input

8
11
16
0

Sample Output

Case 1: 1
Case 2: 2
Case 3: 0

我就不多说了，直接上大神的思路。

0. 回顾
欧拉函数： φ(n)= n(1 – 1/p1)(1 – 1/p2)…(1 – 1/pm) 。
欧拉定理：若 (a,p) =1 ，则 a^φ(p)≡ 1 (mod p) 。
不大于 n(n>=2) 且与 n 互质的所有正整数的和为 n*φ(n)/2 。
1. 阶
满足 a^x≡ 1 (mod p) 的最小的 x 称为 a 关于 p 的 阶 ，记为 δ p (a) ，简写为 δ (a) 。
性质 1 ： a^x≡ 1 (mod p)  ó   δ p (a) | x 。
性质 2 ： δ (a^k)=  δ (a) / ( δ (a), k)
性质 3 ： ( δ (a),  δ (b))= 1  ó   δ (ab) =  δ (a) δ (b) 。
【阶的求法】
由于 a^φ(p) ≡ 1 (mod p) ，所以阶一定是 φ(p) 的约数。
求出 φ(p) 以及它的所有约数，放到一个数组中排序，从小到大用快速幂依次验证，第一个满足的就是答案。
时间复杂度上界为 O(sqrt(p)log(p)) ，一般情况下远远达不到。
By ——lyd尽管我很想pingback，可惜没这个选项(TNT)言归正传：
进入列方程Mode：
8/9*(10^m-1)= k*L
8*(10^m-1) [b]=9kL[/b]
[b][b]8*(10^m-1)/gcd(8,9L)=9kL/gcd(8,9L)[/b][/b]
8/gcd(8,9L)*(10^m-1) [b]=9L/gcd(8,9L) *k[/b]
[b]由定义(a/gcd(a,b),b/gcd(a,b))=1：( 8/gcd(8,9L),9L/gcd(8,9L))=1[/b]
10^m-1 | 9L/gcd(8,9L)
10^m  ≡1 (mod 9L/gcd(8,9L))
(8,9)=1 ->gcd(8,9L)=gcd(8,L)
->10^m  ≡1 (mod 9L/gcd(8,L))
值得称道的是MulMod 写法：
拆成2进制， a*b=a*2^c1+a*2^c2+a*2^c3+...+a82^cn
 

解题思路：因为M全部由8组成，即M=（10^x -1）*8/9=k*N;
则 (10^x-1)*8/gcd(8,N)=9*k*N/gcd(8,N);
令p=8/gcd(8,N); q=9*N/gcd(8,N); 即 (10^x-1)*p=k*q;
由于p和q互质，则(10^x-1)%q==0;
根据同余定理可知，10^x ≡1（mod q）
根据欧拉定理可知当gcd(a,b)==1时，a^φ(b)≡1(mod b);
即可得出：当gcd(10,q)==1时 10^φ(q)≡1(mod q) 即通过枚举φ(q)的因子（最小因子）就能得出结果

比较靠谱的推法。首先，888...=111...*8=(10^0+10^1+...+10^m-1)*8=(10^m - 1)/9*8，PS：m代表888...的长度。
好吧，终于化成指数了，现在有8*(10^m-1)/9=K*L，最小的m就是我们要求的答案啦。
   方式1：
   => 8 * (10^m-1) = 9 * k * L
   => 8/d*(10^m-1)=9*k*L/d，d=gcd(8,9L)
   => 10^m-1 = 0 % 9 * L / gcd(8, 9L) = 0 % 9*L/gcd(8,L)，(由于gcd(8/d,9L/d)=1，那么10^m-1必然是9*L/d的倍数了)。
   => 10^m = 1 % 9 * L / gcd(8,L)
   方式2：
   => 8*(10^m-1)/9 = 0 % L
   => 8*(10^m-1) = 0 % 9*L(这步的推出，比如x/9 = k*n，那么x=9*k*n了，显然成立)
   => 10^m-1 = 0 % 9*L/gcd(9*L,8)，假如，d = gcd(9*L,8)，那么有8/d*(10^m-1)=k*9*L/d，因为8/d不可能是9  *L / d
的倍数，所以10^m-1必定是9*L/d的倍数，所以10^m-1 = 0 % 9*L/gcd(9*L,8))，=>，10^m - 1 = 0 % 9 * L / gcd(L, 8),
(因为gcd(9,8)=1)。
   => 10^m = 1 % 9*L/gcd(8,L) 
   至此，2种方式都推出了，10^m = 1 % 9*L/gcd(8,L) 。
   那么怎么解答这个问题了，这个就用到了欧拉定理了。令p = 9 * L / gcd(8,L)，那么有10^m = 1 % p。由欧拉定理知,Z*p中所有的
数字a均满足a^euler(p) = 1 % p。那么，10只要是p的乘法群中就肯定有解了。如果，10不在Z*p中了，那么有10^m= 2^m * 5^m。
而且10和p有公告因子2或者5，所以p = 2 * k或者 p = 5 * k，2^m=0%p或者5^m=0%p，那么10^m就永远不可能是1%p了。
   综上所述，要满足式子a^m=1%p，必须gcd(p,a)=1，即a必须是p的乘法群中的数字。
   现在的问题是求最小的m，由欧拉定理知道a^euler(p)=1%p，m再大就开始循环了。但是m可能会更小。比如，我们现在知道最小的m
是min，那么有a^min=1%p，因为要满足a^euler(p)=1%p，那么a^euler(p)肯定能变换成(a^min)^k,至于k是多少就不知道了，当然
也可以求出来。那么min就是euler(p)的一个因子，而且是最小的一个满足a^min=1%p的因子了。
   现在就可以通过枚举euler(p)的因子，找到最小的因子min满足式子a^min = 1 % p就能解决本问题了。

#include<stdio.h>
#include<string.h>
#include<algorithm>
#define ll __int64
using namespace std;
ll a[10010];
ll gcd(ll a,ll b)
{
return b?gcd(b,a%b):a;
}
ll eular(ll n)
{
ll res=1,i;
for(i=2;i*i<=n;i++)
{
if(n%i==0)
{
res*=(i-1);
n/=i;
while(n%i==0)
{
res*=i;
n/=i;
}
}
}
if(n!=1)
res*=(n-1);
return res;
}
ll mumod(ll a,ll b,ll m)
{
ll ans=0;
while(b)
{
if(b&1)
ans=(ans+a)%m;
a=(2*a)%m;
b>>=1;
}
return ans;
}
ll exmod(ll p,ll n,ll m)
{
ll ans=1;
p%=m;
while(n)
{
if(n&1)
ans=mumod(ans,p,m);
p=mumod(p,p,m);
n>>=1;
}
return ans%m;
}
int main()
{
ll T=1;
ll n,m,g,mm,i,j,k;
ll kk;
while(scanf("%I64d",&m),m)
{
printf("Case %I64d: ",T++);
g=gcd(m,8);
mm=9*m/g;
kk=0;
if(gcd(mm,10)!=1)
{
printf("0\n");
continue;
}
k=eular(mm);
ll cnt=k;
for(i=1;i*i<=k;i++)
{
if(k%i==0)
{
if(exmod(10,i,mm)==1)
{
cnt=i;
break;
}
if(exmod(10,k/i,mm)==1)
cnt=min(cnt,k/i);
}
}
printf("%I64d\n",cnt);
}
return 0;
}