第四章 数学知识（三）（高斯消元、求组合数、卡特兰数）

一、高斯消元

1、O（n^3）时间复杂度求解n元线性方程组

2、三种初等行变换不影响方程组的解（通过初等行变换得到一个上三角形式）

变换之后的三种情况 

变换过程：枚举每一列找到绝对值最大的一行，交换到第一行，变换为1，将其他行这一列的数化为0。每次进行完操作后在剩余的行中再去找某一列绝对值最大的行进行相同的变换。 

对于n=3的数，进行两次变换就可以得到上三角。（如下得到完美的阶梯型）

 从下往上消去未知数： 

高斯消元法求解线性方程组 

#include<bits/stdc++.h>
//883 高斯消元
using namespace std;
typedef long long LL;
const int N=110;
const double eps=1e-8;
double a[N][N];
int n;
//t用来存绝对值最大的一行
int gauss()
{
    int r,c;
    //逐列逐行化简
    for(c=0,r=0;c<n;c++)
    {
        int t=r;//找绝对值最大的行
        for(int i=r;i<n;i++)
            if(fabs(a[i][c])>fabs(a[t][c]))t=i;
        //注意有一种情况整个c列都是0
        if(fabs(a[t][c])<eps)continue;

        //那一列绝对值最大的一行移到第一行
        //而且只需要用到当前列后面的那些数
        for(int i=c;i<=n;i++)swap(a[t][i],a[r][i]);
        //将绝对值最大的那个数变成1
        for(int i=n;i>=c;i--)a[r][i]/=a[r][c];
        //将那一列的其余行变化为0

        for(int i=r+1;i<n;i++)
        {
            //因为这一行的每一个数都要用到c列对应的那个数
            //要最后改那个数
            //先判断有没有必要去遍历
            if(fabs(a[i][c])>eps)
                for(int j=n;j>=c;j--)
                   a[i][j]-=a[i][c]*a[r][j];
            //想一下:首先看c列当前行比t行是多少倍，然后将r行的这些倍减到i行
        }
        r++;

    }

    //当所有的列都遍历了一遍之后，由于可能存在中间某一列全为0，导致后面r没有加加
    //导致最后c到n了，而r还没有到n
    //由于前面已经得到了所有的变量
    //剩余式子的参数全部为0了
    //r到n中间的式子如果两边都是0有无数解。
    //有一个不是0，就无解
    if(r<n)
    {
        for(int i=r;i<n;i++)
            if(fabs(a[i][n])>eps)return 2;
        return 1;
    }


    //从后往前逐步更新。依次得到xn ……
    //只想利用后面的式子得到本行需要求出的未知数
    for(int i=n-1;i>=0;i--)
    {
        for(int j=i+1;j<n;j++)
        {
            //从i+1到n依次为当前行后面的未知数的下标，
            //同时也是当前行下面的所有行的行号
            //一个未知数一个未知数地处理
            a[i][n]-=a[i][j]*a[j][n];
            //因为每一行的要求的未知数的参数已经变成1了。
            //每次更新只更新a[i][n]的位置。
        }
    }
    return 0;


}
int main()
{
    cin>>n;
    for(int i=0;i<n;i++)
        for(int j=0;j<n+1;j++)
             cin>>a[i][j];
    int t=gauss();
    if(t==1)cout<<"Infinite group solutions"<<endl;
    else if(t==2)cout<<"No solution"<<endl;
    else
    {
      for(int i=0;i<n;i++)
        printf("%.2lf\n",a[i][n]);
    }
    return 0;

}

二、 求组合数

1、组合数的定义

2、组合数可以分解为两部分

a个苹果中有b个是红苹果，第一次拿了一个苹果，第二次拿的时候就有两种情况。将两种情况的次数相加得到。

、

（一）递推法求组合数

ab范围：

#include<bits/stdc++.h>
//885 求组合数
using namespace std;
const int N=2010,mod=1e9+7;
int c[N][N],n;
void init()
{
   for(int i=0;i<N;i++)
    for(int j=0;j<=i;j++)
       if(!j)c[i][j]=1;
       else c[i][j]=(c[i-1][j]+c[i-1][j-1])%mod;

}

int main()
{
    init();
    cin>>n;
    while(n--)
    {
        int a,b;
        cin>>a>>b;
        cout<<c[a][b]<<endl;
    }


}

（二）逆元法求组合数

n和ab的范围：

首先知道：

组合数的逆元组成：

#include<bits/stdc++.h>
//885 求组合数
using namespace std;
typedef long long LL;
const int N=100010,mod=1e9+7;
int fact[N],infact[N];

int qmi(int a,int k,int p)
//快速幂
{
    int res=1;
    while(k)
    {
        if(k&1)res=(LL)res*a%p;
        a=(LL)a*a%p;//k的相邻位置平方的关系
        k>>=1;
    }
    return res;
}
int main()
{
   fact[0]=1,infact[0]=1;
   for(int i=1;i<N;i++)
   {
       fact[i]=(LL)fact[i-1]*i%mod;
       //阶乘的逆
       infact[i]=(LL)infact[i-1]*qmi(i,mod-2,mod)%mod;
   }
    int n;
    cin>>n;
    while(n--)
    {
     int a,b;
     cin>>a>>b;
     cout<<(LL)fact[a]*infact[b]%mod*infact[a-b]%mod<<endl;
    }
    return 0;

}

（三）卢卡斯定理求组合数

（1）询问次数n和ab的范围：

（2）采用ab分别取模的方式化简后求组合数:

（3）

#include<bits/stdc++.h>
//887 求组合数
using namespace std;
typedef long long LL;
const int N=100010,mod=1e9+7;
int fact[N],infact[N];

//使用快速幂的方法求逆，当a和p互逆的时候a^(p-2)模p为1
int qmi(int a,int k,int p)
{
  int res=1;
  while(k)
  {
      if(k&1)res=(LL)res*a%p;
      a=(LL)a*a%p;
      k>>=1;
  }
  return res;
}

//定义法求组合数
int C(int a,int b,int p)
{
    if(b>a)return 0;
    int res=1;
    for(int i=1,j=a;i<=b;i++,j--)
    {
       //根据定义分子分母的乘积的个数是一样的，
       res=(LL)res*j%p;//先乘以j
       res=(LL)res*qmi(i,p-2,p)%p;
    }
    return res;
}

//卢卡斯定理
int lucas(LL a,LL b,int p)
{
    if(a<p&&b<p)return C(a,b,p);
    return (LL)C(a%p,b%p,p)*lucas(a/p,b/p,p)%p;
}


int main()
{
    int n;
    cin>>n;
    while(n--)
    {
     LL a,b;
     int p;
     cin>>a>>b>>p ;
     cout<<lucas(a,b,p)<<endl;
    }
    return 0;

}

（四）高精度组合数求解

（运算过程中不取模）

使用质因数分解的思想。以及一个定理：

将一个数的阶乘用多个质数指数次幂的乘积表示，而某个质因子在阶乘中的指数次幂的计算方式如下所示。

#include<bits/stdc++.h>
//888 求组合数
using namespace std;
typedef long long LL;
const int N=100010,mod=1e9+7;
int primes[N],st[N],sum[N];
//高精度乘法
vector<int> mul(vector<int> a,int b)
{
    int t=0;
    vector<int>c;
    for(int i=0;i<a.size();i++)
    {
        t+=a[i]*b;
        c.push_back(t%10);
        t/=10;
    }
    while(t)
    {
      c.push_back(t%10);
      t/=10;
    }
    return c;

}

//n里面有几个p
int get(int n ,int p)
{
    int res=0;
    while(n)
    {
        res+=n/p;
        n/=p;
    }
    return res;
}

//求素数
int get_Prime(int n)
{
    int res=0;
    for(int i=2;i<=n;i++)
    {
        if(!st[i]) primes[res++]=i;
        //注意这个地方的边界是等于号
        for(int j=0;primes[j]<=n/i;j++)
        {
           st[i*primes[j]]=1;
           if(i%primes[j]==0)break;
        }
    }
    return res;
}

int main()
{
    int a,b;
    cin>>a>>b;
    //求a以内的所有prime
    int res=get_Prime(a);
    //求每个因子在最后的结果里面的次数
    for(int i=0;i<res;i++)
    {
        int p=primes[i];
        sum[i]=get(a,p)-get(b,p)-get(a-b,p);
    }
    vector<int>ans;
    //ans用来乘高精度乘法结果从低位存到高位
    ans.push_back(1);//记得先初始化为1
    for(int i=0;i<res;i++)
    {
        for(int j=1;j<=sum[i];j++)
            ans=mul(ans,primes[i]);
    }
    for(int i=ans.size()-1; i>=0;i--)
    {
        cout<<ans[i];
    }
    puts("");
    return 0;
}

（五）卡特兰数

889满足条件的01序列，要求n个1和n个0组成一个序列，并且要求所有前缀都是0的个数大于1的个数。

当n为6的时候，向右为0，向上为1，结果就是所有红线下面的到（6,6）这个点的路径。

只要经过红线的，在第一次经过红线的点关于红线做轴对称，最后一定会经过(5,7)。所有到（6,6）的路径的条数-所有到（5,7)的条数即为所有合法条数

化为n

最后化简

#include<bits/stdc++.h>
//889 01序列
using namespace std;
typedef long long LL;
const int mod = 1e9+7;
//p为质数的时候逆元的算法
int qmi(int a,int b ,int p)
{
    int res=1;
    while(b)
    {
        if(b&1)res=(LL)res*a%mod;
        a=(LL)a*a%mod;
        b>>=1;
    }
    return res;
}
int main()
{
    int n;
    cin>>n;
    int a=2*n,b=n;
    int res=1;
    for(int i=a,j=1;j<=b;j++,i--)
    {
        res=(LL)res*i%mod;
        res=(LL)res*qmi(j,mod-2,mod)%mod;
    }
    res=(LL)res*qmi(n+1,mod-2,mod)%mod;
    cout<<res<<endl;
}