您的位置:首页 > 其它

哈尔滨理工大学-CPC23 2014-3-K-喵喵的神·数

2014-10-04 15:59 323 查看
K.喵喵的神·数
Time Limit: 1000 MSMemory Limit: 2560000 K
Total Submit: 193 (49 users)Total Accepted: 27 (26 users)Special Judge: No
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;
}






内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: