哈尔滨理工大学-CPC23 2014-3-K-喵喵的神·数
2014-10-04 15:59
323 查看
K.喵喵的神·数 | |||||
| |||||
Description | |||||
喵喵对组合数比较感兴趣,并且对计算组合数非常在行。同时为了追求有后宫的素质的生活,喵喵每天都要研究质数。 我们先来复习一下什么叫做组合数。对于正整数P、T 然后我们再来复习一下什么叫质数。质数就是素数,如果说正整数N的约数只有1和它本身,N就是质数;另外,1不是质数。 今天,喵喵想要知道 | |||||
Input | |||||
输入第一行是一个整数N(N<=1000)。 接下来N行,每行包括一个正整数T和一个质数P(1<=P<=T<231)。 | |||||
Output | |||||
包括N行,根据输入的顺序,每一行为一个整数: | |||||
Sample Input | |||||
2 3 2 10 3 | |||||
Sample Output | |||||
1 0 正常的算法: 一》C(m, n) % p (p为素数, m, n <= 10^8) 算法伪代码: 1) 预处理 f[x] = (x!) mod p. (0 <= x < p) 2) 求m!中把p因子去掉,余下的数乘起来 mod p. 函数 G(m)求解 2)问题: 伪代码: G(m) 1) if (0 == m) return 1; 2) m = k*p + t; 3) G[m] = f[p-1] ^ k * f[t] * G[m / k] mod p; (其中 f[p-1]^k mod p 二分做) 3) C(m, n) = m! / (n! * (m-n)!); A = m! mod p; B = n!*(m-n)! mod p; C = A / B mod p; 在这里有俩种算法: 1) 求B的逆元; ( O(n^2) 枚举, 或者.. ) 2) 扩展欧几里德求 BC = A(mod p) (= 是同余的意思) 最小的 C 就是答案; 这里的 A, B 转换成 对 2)问题 求解. #include <iostream> #include <map> #include <set> #include <string> #include <algorithm> #include <math.h> #include <queue> using namespace std; template <class T> void out(T x, int n){ for(int i = 0; i < n; ++i) cout << x[i] << ' '; cout << endl; } template <class T> void out(T x, int n, int m){ for(int i = 0; i < n; ++i) out(x[i], m); cout << endl; } #define OUT(x) (cout << #x << " = " << x << endl) #define FOR(i, a, b) for(int i = (int)a; i < (int)b; ++i) #define REP(i, b) FOR(i, 0, b) #define FORD(i, a, b) for(int i = (int)a; i >= (int)b; --i) #define N 10007 #define MOD(x, y) ((x) % (y) + y) % (y) long long f[N+1]; long long ext_gcd(long long a, long long b, long long &x, long long &y) { if (0 == b) { x = 1; y = 0; return a; } else { long long tx, ty, d; d = ext_gcd(b, MOD(a, b), tx, ty); x = ty; y = tx - a / b * ty; return d; } } long long modular_linear(long long a, long long b, long long n) { long long tx, ty, d, x0, i; d = ext_gcd(a, n, tx, ty); if (0 == b % d) { x0 = MOD(tx * b / d, n); return x0; } else { return 1; } } void init(){ f[0] = 1; FOR (i, 1, N) { f[i] = (f[i-1] * i) % N; } } // m! 有多少p因子 long long calP(long long m, long long p){ long long ret = 0; while (m) { ret += m / p; m /= p; } return ret; } long long square(long long x){ return x * x; } long long bigmod(long long b, long long p, long long m) { if (0 == p) return 1; else if (0 == p % 2) return square(bigmod(b, p/2, m)) % m; else return ((b % m) * bigmod(b, p-1, m)) % m; } // 计算m!-mod p long long calM(long long m, long long p){ if (0 == m) return 1; long long ret = 0, k, t; k = m / p; t = m % p; ret = (bigmod(f[p-1], k, N) * f[t]) % p; ret *= calM(m / p, p); return ret % p; } long long calCmnModP(long long m, long long n, long long p){ long long A, B, u, v; u = calP(m, p); v = calP(n, p) + calP(m-n, p); if (u > v) return 0; else { A = calM(m, p); B = (calM(n, p) * calM(m-n, p)) % p; return modular_linear(B, A, p); } } int main(){ int t, n, m; scanf("%d", &t); init(); while (t--) { scanf("%d %d", &m, &n); printf("%d\n", calCmnModP(m, n, N)); } return 0; } 二》<h2 class="entry_title">大整数求组合数取余(Lucas定理)</h2><div class="archive_info"><span class="date">2013年02月25日</span><span class="category"> ⁄ 综合</span>⁄ 共 915字 ⁄ 字号 <span class="font"><a>小</a> <a>中</a> <a>大</a></span><span class="comment"> ⁄ <span>评论关闭</span></span><span class="edit"></span></div><div class="ad_r"></div><p>【卢卡斯(Lucas)定理】</p><p>Lucas定理用来求C(a,b)mod p的值,其中p为素数。</p><p>数学表达式为:</p><p>Lucas(a,b,q)=C(a%q,b%q)*Lucas(a/p,b/p,p);</p><p>Lucas(a,0,q)=0;</p><p>通过这个定理就可以很方便的把大数的组合转化成小数。但其中还是要求C(a%q,b%q)%p,所以这里引入逆元来求。</p><p>【定义】若整数a,b,p, 满足a·b≡1(mod p).则称a 为b 模p 的乘法逆元, 即a=b- 1mod p.其中, p 是模数。</p><p>应用到组合数中来就是:</p><p> a!/[b!*(a-b)!] % p == a! * [b!*(a-b)!]-1 %p</p><p>【逆元求法】:</p><p>应用费马小定理,ap-1=1 mod p ,即 a*ap-2=1 mod p</p><p>也就是说 ap-2就是a的逆元。</p><p>当然这里求出来的逆元是在取模p的逆元,对我们最终目标没有影响。这也是比较方便而且比较好的方法。</p><p>弄个模板来。</p><p> </p><p> </p><p><pre name="code" class="cpp">#include <cstdio> #include <iostream> #include <cmath> #include <cstring> #include <algorithm> using namespace std; #define maxn 100010 typedef long long LL; LL m,n,p; LL Pow(LL a,LL b,LL mod) { LL ans=1; while(b) { if(b&1) { b--; ans=(ans*a)%mod; } else { b/=2; a=(a*a)%mod; } } return ans; } LL C(LL n,LL m) { if(n<m) return 0; LL ans=1; for(int i=1;i<=m;i++) { ans=ans*(((n-m+i)%p)*Pow(i,p-2,p)%p)%p; } return ans; } LL Lucas(LL n,LL m) { if(m==0) return 1; return (Lucas(n/p,m/p)*C(n%p,m%p))%p; } int main() { int t; scanf("%d",&t); while(t--) { scanf("%lld%lld%lld",&n,&m,&p); printf("%lld\n",Lucas(n,m)); } return 0; } 三》<h2 id="t_618194800100l5rv" class="titName SG_txta">大数取余</h2> <div id="sina_keyword_ad_area2" class="articalContent "><p>有一类题目会因为求出的结果太大而只要求输出对某个数m取余后的结果,而且这个m是比较小的数,比如不超过32位整数… 而这类大数都是可以由较小的数经过某些运算得到的… 于是我整理了一下对付几种运算的方法…包括四则运算,指数,组合数,塔函数的应对方法…</p><p>那么就开始吧…慢慢来… 首先是最常识的加减法:</p><pre>add_mod(a,b,m){ return ((a%m)+(b%m))%m; } 别小看加法哦…很多用dp解的题目中靠着加法可是能达到很大的数呢… minus_mod(a,b,m){ return (a-b+m)%m; } 减法…会遇到吗? 接着是依然很简单的乘法: multiply_mod(a,b,m){ return ((a%m)*(b%m))%m; } 这是当m*m不会溢出时可以用的,同时也是通常的情况… 但是如果m*m连long long都会溢出的话…就需要把一个数一位位拆开来做了… multiply_mod(a,b,m){ if(b==0)return 0; return (((b&1)?a:0)+(multiply_mod(a,b>>1,m)<<1)%m)%m; } 然后是除法,但有点限制: (a/b)%m 特殊条件:m和b互质 前提:a能被b整除 这个有点特殊,意为虽然不知道a是多少,但是已知c,而且c=a%m,用c和b来求(a/b)%m的方法 虽然需要m和b互质,但是不互质的话也是可以做的,因为a也一定是gcd(b,m)的倍数,具体可以看看这里 需要用到扩展欧几里德来求… 至于扩展欧几里德是什么…去Google一下吧… extend_euclid(a,b,&x0,&y0){ if(b==0){ x0=1; y0=0; return a; } r=extend_euclid(b,a%b); t=x0; x0=y0; y0=t-a/b*y0; return r; } divide_mod(a,b,m){ extend_euclid(b,m,x0,y0); return (a*(((x0%m)+m)%m))%m; } 这个在求组合数的时候可能用到… 不过似乎很少遇到需要用除法取余的情况呢… 然后是更大却很简单的幂运算: (a^b)%m 这是初学递归或者二分时就会遇上的一个很简单的方法,和之前的乘法差不多 power_mod(a,b,m){ if(b==0)return 1; if(b&1)return (a*power_mod((a*a)%m,b>>1,m))%m; else return power_mod((a*a)%m,b>>1); } 其中乘法取余会运算中溢出的话可以改成之前那个multiply_mod() 恩…开始有趣了…下面是组合数: C(a,b)%m 特殊条件:m是质数 如果b或者a-b比较小的话,可以用之前计算(a/b)%m的方式来把组合数公式展开来计算 不过当b和a-b与m相比很大时,有更好的方法: a,b在m进制下表示为 a=(ak,…,a0),b=(bk,…,b0)0=<ai,bi<m 于是会有这样的性质: C(a,b)=C(ak,bk)*...*C(a0,b0) (mod m) 最后是难以想象的大的数…塔函数: (a↑↑b)%m 这里可以看到其一些性质 比如ProjectEuler 282 Ackermann 函数就是非常恶心的大数 第一层是很小的常数,第二层是n个a相加,也就是乘法,第三层是n个a相乘,也就是幂,第四层是n个a叠着做幂即塔函数,第k+1层是n个a做第k层运算… 用小数字居然也能表示出如此之大的数…佩服啊… 在上一篇中也有提到: n>=phi(m)时,a^n=a^(n%phi(m)+phi(m)) (mod m) 其中phi()是欧拉函数 由于phi(x)<x,所以就算是对a↑↑b这样大的数也总会在足够短的时间内收敛,然后计算出(modm)的值 这个不止在塔函数中,也可以用在各种指数异常大的情况下 特别的,在b>m的情况下,(a↑↑b)%m的值将是定值 恩…就到这里了… PS:逆元素-http://baike.baidu.com/view/3140397.htm?from_id=11054145&type=syn&fromtitle=%E9%80%86%E5%85%83&fr=aladdin 非正常的思路: (找到了余数的规律 = = 。。) 代码还过了。。。 卢卡思、费马、欧几里德怕不是要死的活过来!??? #include<stdio.h> int main() { int i,n,p,t; while (scanf("%d", &n) != EOF) { for (i = 0; i < n; i++) { scanf("%d%d", &t, &p); printf("%d\n", (t / p) % p); } } return 0; } |
相关文章推荐
- CPC23-4 K.喵喵的神·数
- 哈尔滨理工大学-CPC23 2014-2-D-Alice and Bob
- [POJ 2014] Flow Layout · 模拟
- CPC23-4-K. 喵喵的神数 (数论 Lucas定理)
- CPC23 2014-1 B.Crossings(求逆序数)
- CPC23 2014-1 B.Crossings(求逆序数)
- SQL Server 2014 使用(持续更新中···)【目录】
- CPC23-4-K. 喵喵的神数 (数论 Lucas定理)
- 2014游戏开发者大会 · 中国(GDC China)
- 2014 Unity亚洲开发者大会·中国站:技术课程全面解析
- 2014中学教师资格证视频下载 教育学 教育心理学 全真模拟
- android一些可用实例2014收藏
- 漫谈构建运维平台化2014运维线下沙龙(北京站) | 峰云就她了
- BZOJ P3504[Cqoi2014]危桥
- 华为2014机试样题及代码(三)
- MZ Training 2014 #6 C题
- BZOJ 3505: [Cqoi2014]数三角形( 组合数 )
- bzoj3527 [Zjoi2014]力
- ··· Socks Online - 在内部网也能使用QQ ···
- 2014第11周四Eclipse开发问题记