组合数求法讲解

  首先对于c(n,m),如果n,m比较大的话,这个值会很大,为了简化变成复杂度,也是为了更好的求解,都要求求c(n,m) mod p的值,我们由最简单的问题慢慢提高难度。

  要求求解c(n,m) mod p,p为质数,且n,m<=10^5,这时我们可以存储fac[i]=i! mod p.这样就可以直接用公式求了。

  当n,m的数据范围为10^9,且p为质数的时候且<=10^5,显然没法用数组存储10^9的信息,但是我们可以记录10^5的阶乘信息。然后我们引入lucas定理。

c(n,m)=lucas(n,m)=lucas(n div p,m div p)*c(n mod p,m mod p)。这样我们可以递归的将10^9的数据范围缩小到每次的10^5。

lucas:

function lucas(x,y,p:int64):int64;
var
    a, b                    :int64;
begin
    if y=0 then exit(1);
    a:=x mod p;
    b:=y mod p;
    if a<b then exit(0) else lucas:=lucas(x div p,y div p,p)*combine(a,b,p);
end;   

  当p不是质数,但是p为互不相同的质数的乘积,且每个质数小于10^5候,我们可以将p分解质因数。分解成k个质数相乘。那么ans=c(n,m) mod p,我们可以得到c(n,m) mod pi=ansi这样的k个方程组,每个方程组解的方法与上一问题相同。利用crt(中国剩余定理)可以将这k个方程组合并,得出原问题的答案。

  当p为若干质数的乘积,且pi^c<10^5时,我们可以由于各个pi^c互质,所以我们假设能求出c(n,m) mod pi^c这k个等式的解,就可以类似上一问题用crt合并出当前问题的解。那么问题就转化成了求解c(n,m) mod pi^c 因为mod的数是非质数,所以无法类似第一个方法求解(因为没有办法求逆元),观察答案形式,ans=c(n,m)=n!/(m!*(n-m)!),我们可以提出来分子和分母中含有pi因子的项的pi因子,上下可以消掉,剩下的就是和pi^ci互质的数,就可以求逆元了。

  最后一个问题bzoj 2142 http://61.187.179.132/JudgeOnline/problem.php?id=2142

//By BLADEVIL
var
    m, n                                 :longint;
    pj, c                               :array[0..110] of longint;
    pi, s, a                            :array[0..110] of int64;
    p                                   :int64;
    tot                                 :longint;
 
procedure divide(p:int64);
var    
    i, j                                :longint;
     
begin
    tot:=0;
    for i:=2 to trunc(sqrt(p)) do
    if p mod i=0 then
    begin
        inc(tot);
        pj[tot]:=i;
        while p mod i=0 do
        begin
                inc(c[tot]);
                p:=p div i;
        end;
    end;
    if p>1 then
    begin
        inc(tot);
        pj[tot]:=p;
        c[tot]:=1;
    end;
    for i:=1 to tot do
    begin
        pi[i]:=1;
        for j:=1 to c[i] do pi[i]:=pi[i]*pj[i];
    end;
end;
 
function ex_gcd(a,b:int64;var x,y:int64):int64;
var    
    t                                   :int64;
begin
    if (b=0) then
    begin
        x:=1;y:=0;
        exit(a);
    end;
    ex_gcd:=ex_gcd(b,a mod b,x,y);
    t:=x;
    x:=y;
    y:=t-(a div b)*y;
end;
 
function gcd(a,p:int64):int64;
var    
    x, y                                :int64;
     
begin
    x:=0;y:=0;
    ex_gcd(a,p,x,y);
    x:=(x mod p+p)mod p;
    exit(x);
end;
 
function mi(x,y,q:int64):int64;
var
    rec                                 :int64;
     
begin
    rec:=1;
    while (y>0) do
    begin
        if y and 1=1 then rec:=rec*x mod q;
        x:=x*x mod q;
        y:=y shr 1;
    end;
    exit(rec);
end;
 
function fac(n,p,q:int64):int64;
var
    cnt                                 :int64;
    i                                    :longint;
begin
    cnt:=1;
    for i:=1 to n do
        if (i mod p>0) then
            cnt:=cnt*i mod q;
    exit(cnt);
end;
 
function fact(n:int64;var sum:int64;p,q:int64):int64;
var
    cnt, rec                            :int64;
     
begin
    rec:=1;
    cnt:=fac(q,p,q);
    while n>=p do
    begin
        sum:=sum+n div p;
        if (n div q>0) then rec:=rec*(mi(cnt,n div q,q) mod q)mod q;
        if (n mod q>0) then rec:=rec*(fac(n mod q,p,q)mod q) mod q;
        n:=n div p;
    end;
    if n>1 then rec:=rec*fac(n,p,q) mod q;
    exit(rec);
end;
 
function combine(n,m,p,q:int64):int64;
var    
    ans1, ans2, ans3, ans               :int64;
    a, b, c                             :int64;
     
begin
    a:=0;b:=0;c:=0;
    ans1:=fact(n,a,p,q);
    ans2:=fact(m,b,p,q);
    ans3:=fact(n-m,c,p,q);
    a:=a-(b+c);
    ans:=mi(p,a,q);
    ans:=ans*ans1 mod q;
    ans:=ans*gcd(ans2,q) mod q;
    ans:=ans*gcd(ans3,q) mod q;
    exit(ans);
end;
 
function doit(n,m:longint):int64;
var
    i                                   :longint;
    x, y, sum                           :int64;
         
begin
    sum:=0;
    for i:=1 to tot do
            a[i]:=combine(n,m,pj[i],pi[i]);
    for i:=1 to tot do
    begin
        x:=0;y:=0;
        ex_gcd(s[i],pi[i],x,y);
        x:=(x mod pi[i]+pi[i])mod pi[i];
        sum:=(sum+((x*s[i] mod p)*a[i])mod p)mod p;
    end;
    exit(sum mod p);
end;
 
procedure main;
var    
    i                                   :longint;
    w                                   :array[0..100] of longint;
    ans                                 :int64;
    sum                                 :int64;
     
begin
    readln(p);
    divide(p);
    for i:=1 to tot do s[i]:=p div pi[i];
    readln(n,m);
    sum:=0;
    for i:=1 to m do
    begin
        readln(w[i]);
        inc(sum,w[i]);
    end;
    if sum>n then
    begin
        writeln('Impossible');
        exit;
    end;
    ans:=1;
    for i:=1 to m do
    begin
        ans:=ans*doit(n,w[i]) mod p;
        n:=n-w[i];
    end;
    writeln(ans);
end;
 
begin
    main;
end.

 

转载于:https://www.cnblogs.com/BLADEVIL/p/3511190.html

### C语言中的回溯法算法讲解 #### 定义基本概念 回溯法是一种通过构建解空间树并尝试每一种可能路径的方法来解决问题的技术。这种方法适用于那些需要探索所有可能性才能找到解决方案的情况,比如排列、组合等问题。 #### 参数设计原则 在编写基于回溯法的程序时,函数通常被设置成`void`类型[^2]。这是因为回溯过程中并不直接返回某个特定的结果给调用者;相反,它会修改传入的数据结构或者全局变量作为输出的一部分。因此,在定义这样的递归函数时,重点在于如何传递必要的状态信息以及控制流程继续否的标准。 #### 终止条件设定 当采用回溯策略解决具体问题时,必须明确规定何时停止进一步深入搜索。这往往依赖于当前处理的状态是否满足预设的目标或者是达到了某种边界情况。例如,在计算全排列的过程中,一旦形成了完整的序列就可以结束本轮迭代,并且回退至上一层级重新选择其他分支进行试探。 #### 实际应用案例——素数环问题 考虑这样一个例子:构造一个由n个连续整数组成的圆圈(即所谓的“素数环”),使得任意相邻两数相加均为质数。可以利用一维数组`home[]`表示这个圆形布局,并指定首个位置上的数以简化后续操作[^4]。 以下是针对上述描述的具体实现: ```c #include <stdio.h> #define MAX_SIZE 100 int isPrime(int num){ if (num <= 1) return 0; for (int i = 2; i * i <= num; ++i) if (num % i == 0) return 0; return 1; } // 判断新加入节点后的链表是否合法 bool isValid(int* path, int size){ if (!isPrime(path[size - 1] + path[(size - 2 + size) % size])) { return false; } return true; } void findPrimerCircle(int n, int pos, bool used[], int circle[]){ if(pos >= n && isPrime(circle[n-1]+circle[0])){ printf("Solution found:\n"); for(int i=0;i<n;++i){ printf("%d ",circle[i]); } putchar(&#39;\n&#39;); return ; } for(int i=1;i<=n;++i){ if(!used[i]){ circle[pos]=i; used[i]=true; if(isValid(circle,pos+1)){ findPrimerCircle(n,pos+1,used,circle); } // 回滚上一步的选择 used[i]=false; } } } ``` 这段代码展示了怎样运用回溯的思想去寻找符合条件的素数环配置方案。其中包含了几个重要的组成部分: - `isValid()`用于验证新增添成员之后的整体合法性; - 主体部分则负责按照既定规则逐步填充候选列表直至完成整个循环或是发现矛盾而提前中断。
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值