[数论第三节]高斯消元法/求组合数/卡特兰数

• 某行乘以一个非零数
• 交换两行
• 某行加上另一行的若干倍

• 完美阶梯型：唯一解
• 非完美阶梯型：
  • 0 == 非0：无解
  • 0 == 0：无穷解

• 枚举每一列
  • 找到这一列系数的绝对值最大的一行
  • 将这一行与第一行交换
  • 将改行的第一个数变成一(方程两边同乘某数)
  • 把下面所有行的当前列的系数消成0(某行加上第一行的若干倍)

const int N = 110;
const double esp = 1e-6; //x < esp,则x = 0,否则 x != 0,由于浮点数精度问题,0可能是0.000001
double a[N][N];
int n;

int gauss(){
	int r, c;
	for(r = 1, c = 1; c <= n; ++ c){//遍历每一列
		int t = r;
		for(int i = r; i <= n; ++ i) //找当前列的系数最大的行
			if(fabs(a[t][c]) < fabs(a[i][c]))
				t = i; //记录最大行
		
		if(fabs(a[t][c]) < esp) continue; //最大行系数为0说明该列系数均为0
		
		for(int i = c; i <= n + 1; ++ i) swap(a[r][i], a[t][i]);//否则交换第一行与系数最大行
		for(int i = n + 1; i >= c; -- i) a[r][i] /= a[r][c]; //将第一行当前列系数变为1
		
		for(int i = r + 1; i <= n; ++ i) //将下面所有行的当前列的系数变为0
			if(fabs(a[i][c]) > esp) //当前列系数非0
				for(int j = n + 1; j >= c; -- j) 
					a[i][j] -= a[r][j] * a[i][c];

		r ++ ;
	}
	
	if(r <= n){
		for(int i = r; i <= n; ++ i)
			if(fabs(a[i][n + 1]) > esp) return 0; //无解
		return 1; //无数解
	}
	
	for(int i = n; i >= 1; -- i) //从后往前求解
		for(int j = i + 1; j <= n; ++ j)
			a[i][n + 1] -= a[i][j] * a[j][n + 1];
			
	return 2; //唯一解
}
int main(){
	cin >> n;
	
	for(int i = 1; i <= n; ++ i)
		for(int j = 1; j <= n + 1; ++ j)
			cin >> a[i][j];
	
	int t = gauss();
	
	if(t == 1) cout << "Infinite group solutions";
	else if(t == 0) cout << "No solution";
	else if(t == 2){
		for(int i = 1; i <= n; ++ i) 
			printf("%.2fn", a[i][n + 1]);
	}
	
	return 0;
}

//询问10000次
//a,b <= 2000
const int N = 2e3 + 10;
const int mod = 1e9 + 7;
int c[N][N];
int 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 - 1] + c[i - 1][j]) % mod;
}

• 开二维数组会爆内存，于是利用公式，预处理出阶乘以及阶乘的逆(O(Nlog(N)))
• 代码：

  typedef long long LL;
  const int N = 1e5 + 10;
  const int mod = 1e9 + 7;
  int fact[N], infact[N];
  int n;
  
  //快速幂
  int qmi(int a, int b, int p){
      int res = 1;
      while(b){
          if(b & 1) res = (LL)res * a % p;
          a = (LL)a * a % p;
          b >>= 1;
      }
      return res;
  }
  
  //预处理各阶乘以及阶乘的逆
  void init(){
      fact[0] = infact[0] = 1;
      for(int i = 1; i < N; ++ i){
          fact[i] = (LL)fact[i - 1] * i % mod; //递推式n! = (n - 1)! * n;
          infact[i] = (LL)infact[i - 1] * qmi(i, mod - 2, mod) % mod; //递推式(n!)^-1 = ((n - 1)!)^-1 * n^-1，取模意义下的逆就是逆元
      }
  }
  
  int main(){
      cin >> n;
      init();
      while(n -- ){
          int a, b;
          cin >> a >> b;
          cout << (LL)fact[a] * infact[b] % mod * infact[a - b] % mod << endl;//两个int最大值相乘不会爆longlong，但是三个相乘会爆
      }
      return 0;
  }
  

• 两数相乘会爆longlong，利用卢卡斯定理预处理组合数 (O(plog(N)log(p)))
• 定理：(C_a^bequiv C_{a mod p}^{b mod p}·C_{a/p}^{b/p} pmod p)
• 代码：

  typedef long long LL;
  int p;
  int n;
  
  //快速幂求逆元
  int qmi(int a, int b){
      int res = 1;
      while(b){
          if(b & 1) res = (LL)res * a % p;
          a = (LL)a * a % p;
          b >>= 1;
      }
      return res;
  }
  
  int C(int a, int b){ //求组合数C_a^b
      if(a < b) return 0; //可以不要
      int x = 1, y = 1; //分子分母
      for(int i = 0; i < b; ++ i){
          x = (LL)x * (a - i) % p; //求分子
          y = (LL)y * (i + 1) % p; //求分母
      }
      return (LL)x * qmi(y, p - 2) % p;
  }
  
  //卢卡斯定理
  int lucas(LL a, LL b){
      if(a < p && b < p) return C(a, b); //不必用卢卡斯，可以直接利用公式求出
      return (LL)C(a % p, b % p) * lucas(a / p, b / p) % p; //用卢卡斯定理
  }
  
  int main(){
      cin >> n;
      while(n -- ){
          LL a, b;
          cin >> a >> b >> p;
          cout << lucas(a, b) << endl;
      }
      return 0;
  }
  

• 此时组合数较大，需要用高精度存储结果
• 阶乘中素数的个数：
  • (v_p(n!) = sum_{i=1}^{infty} leftlfloor frac{n}{p^i} rightrfloor)
  • 证明：将 (n!) 记为 (1times 2times cdots times ptimes cdots times 2ptimes cdots times lfloor n/prfloor ptimes cdots times n) 那么其中 (p) 的倍数有 (ptimes 2ptimes cdots times lfloor n/prfloor p=p^{lfloor n/prfloor }lfloor n/prfloor !) 然后在 (lfloor n/prfloor !) 中继续寻找 (p) 的倍数即可，这是一个递归的过程
• 代码：

  const int N = 5e3 + 10;
  int primes[N], cnt;
  int st[N];
  int sum[N];
  
  //线性筛，1-n中质数就是n!中的质因子
  void get_primes(int a){
      for(int i = 2; i <= a; ++ i){
          if(!st[i]) primes[ ++ cnt] = i;
          for(int j = 1; primes[j] <= a / i; ++ j){
              st[primes[j] * i] = 1;
              if(i % primes[j] == 0) break;
          }
      }
  }
  
  //求a!中质数p的个数
  int get(int a, int p){
      int res = 0;
      while(a){
          res += a / p;
          a /= p;
      }
      return res;
  }
  
  //高精度乘法
  vector<int> mul(vector<int> &A, int b){
      vector<int> c;
      int t = 0;
      for(int i = 0; i < A.size() || t; ++ i){
          if(i < A.size()) t += A[i] * b;
          c.push_back(t % 10);
          t /= 10;
      }
      while(c.size() > 1 && c.back() == 0) c.pop_back();
      return c;
  }
  
  int main(){
  	int a, b;
      cin >> a >> b;
      get_primes(a);
      for(int i = 1; i <= cnt; ++ i){
          int p = primes[i];
          sum[i] = get(a, p) - get(b, p) - get(a - b, p); //a!中质因子个数-b!中质因子个数-(a-b)!中质因子个数
      }
      
      vector<int> res;
      res.push_back(1);
      
      for(int i = 1; i <= cnt; ++ i) //累乘
          for(int j = 1; j <= sum[i]; ++ j)
              res = mul(res, primes[i]);
      
      for(int i = res.size() - 1; i >= 0; -- i)//输出
          cout << res[i];
      
      return 0;
  }
  

#include <iostream>

using namespace std;

typedef long long LL;
const int N = 1e5 + 10;
const int p = 1e9 + 7;
int n;

int qmi(int a, int k){
    int res = 1;
    while(k){
        if(k & 1) res = (LL)res * a % p;
        a = (LL)a * a % p;
        k >>= 1;
    }
    return res;
}

int main(){
    cin >> n;
    
    int res = 1;
    
    int a = 2 * n, b = n;
    int x = 1, y = 1;
    for(int i = 0; i < b; ++ i){
        x = (LL)x * (a - i) % p;
        y = (LL)y * (i + 1) % p;
    }
    
    res = (LL)x * qmi(y, p - 2) % p;
    res = (LL)res * qmi(n + 1, p - 2) % p;
    cout << res;
    
    return 0;
}