hdu5526

小明和他的矩阵

题目描述：

Tom放学回家的路上，看到天空中出现一个矩阵。Tom发现，如果矩阵的行、列从0开始标号，第i行第j列的数记为ai,j，那么ai,j=Cji
如果i < j，那么ai,j=0
Tom突发奇想，想求一个矩形范围内所有数的和。Tom急着回家，当然不会自己算，所以就把任务交给你了。
因为数可能很大，答案对一个质数p取模
输入：输入包含多组数据（大约8组）。每组数据只有一行五个非负整数，x1、y1、x2、y2、p，你要求的是∑x2i=x1∑y2j=y1ai,j模p后的值。
x1≤x2≤1e5,y1≤y2≤1e5,2≤p≤1e9

题解：

发现一列一列的看最好，利用【l，r】区间的方法，算出【0，r】的就好，发现每一列是上标相同的组合数，那么之前加一项，之后再减去就能变成一个组合数。然后每一列算一次累加起来就好了。剩下关键是怎么求组合数并且MODp。本来是用欧拉定理算逆元推，但是要求x和p互质，p是质数，但是如果p太小，比如p是2，x是4，就跪了。于是隆重推出专门算组合数的lucas定理：C【n】[m] = C[n/p][m/p]*C[n%p][m%p]%p;
这样一来，对于C[n%p][m%p]来说，都小于p，预处理出来阶层，可以用逆元直接搞，C[n/p][m/p]是一直递归下去，logn就搞定。

重点：

首先，发现每一列的规律，这是这道题的关键思路。其次，认识到当MOD太小的时候，组合数要用预处理阶层+Lucas定理

代码：

#include <iostream>
#include <cstdio>
#include <cstring>
#include <string>
#include <cmath>
#include <ctype.h>
#include <limits.h>
#include <cstdlib>
#include <algorithm>
#include <vector>
#include <queue>
#include <map>
#include <stack>
#include <set>
#include <bitset>
#define CLR(a) memset(a, 0, sizeof(a))
#define REP(i, a, b) for(ll i = a;i < b;i++)
#define REP_D(i, a, b) for(ll i = a;i <= b;i++)

typedef long long ll;

using namespace std;

const ll maxn = 1e5 + 100;
ll jiecheng[maxn];
ll x_1, x_2, y_2, y_1;
ll p;
void getJiecheng()//预处理阶层
{
    jiecheng[0] = 1;
    REP_D(i, 1, 100000)
    {
        jiecheng[i] = jiecheng[i-1]*i%p;
    }
}
ll pow_mod(ll x, ll n, ll M)//求逆元
{
    if(n==0)
        return 1;
    ll xx = x*x%(M);
    ll nn = n/2;
    ll ans = pow_mod(xx, nn, M);
    if(n%2==1)
    {
        ans = (ans*x)%(M);
    }
    return ans;
}


ll getC(ll n, ll m)//递归版lucas定理求C
{
    if(m > n)//特判
        return 0;
    if(m == 0||n == m)//小边界
        return 1;
    ll a = n%p, b = m%p;
    if(a < b)//特判
        return 0;
    else
        return getC(n/p, m/p)*jiecheng[a]%p*pow_mod(jiecheng[b], p-2, p)%p*pow_mod(jiecheng[a-b], p-2, p);
}
ll gao(ll a, ll b, ll x)//算从列a到列b然后x已经加过一了。具体公式自己推一下就好
{
    if(a + 1 > x)
        return 0;
    ll ans = 0;
    //ll tmp = 1;
//    REP_D(i, 1, a)
//    {
//        tmp = (tmp*(x-i+1)%p*pow_mod(i, p-2, p)%p);
//    }
    REP_D(i, a + 1, b + 1)
    {
        if(i > x)
            break;
        //tmp = (tmp*(x-i+1)%p*pow_mod(i, p-2, p)%p);
        ans = (ans + getC(x, i))%p;
    }
    return ans;
}

void solve()
{
    ll ans = 0;
    ans = gao(y_1, y_2, x_2 + 1);//【l，r】方法
    if(x_1 != 0)
        ans = ((ans - gao(y_1, y_2, x_1))%p + p)%p;
    printf("%I64d\n", ans);
}


int main()
{
    freopen("3Cin.txt", "r", stdin);
    //freopen("3Cout.txt", "w", stdout);
    while(scanf("%I64d%I64d%I64d%I64d%I64d", &x_1, &y_1, &x_2, &y_2, &p) != EOF)
    {
        getJiecheng();
        solve();
    }
    return 0;
}