zoukankan      html  css  js  c++  java
  • ACM学习历程—SNNUOJ 1116 A Simple Problem(递推 && 逆元 && 组合数学 && 快速幂)(2015陕西省大学生程序设计竞赛K题)

    Description

    Assuming a finite – radius “ball” which is on an N dimension is cut with a “knife” of N-1 dimension. How many pieces will the “ball” be cut into most?
    However, it’s impossible to understand the following statement without any explanation.
    Let me illustrate in detail.
    When N = 1, the “ball” will degenerate into a line of finite length and the “knife” will degenerate into a point. And obviously, the line will be cut into d+1 pieces with d cutting times.
    When N = 2, the “ball” will degenerate into a circle of finite radius and the “knife” will degenerate into a line. Likewise, the circle will be cut into (d^2+d+2)/2 pieces with d cutting times.
    When N = 3, the “ball” will degenerate into a ball on a 3-dimension space and the “knife” will degenerate into a plane.
    When N = 4, the “ball” will degenerate into a ball on a 4-dimension space and the “knife” will degenerate into a cube on a 3 dimension.
    And so on.
    So the problem is asked once again: Assuming a finite-radius “ball” which is on an N dimension is cut with a “knife” of N-1 dimension. How many “pieces” will the “ball” be cut into most?

    Input

    The first line of the input gives the number of test cases T. T test cases follow. T is about 300.
    For each test case, there will be one line, which contains two integers N, d(1 <= N <= 10^5, 1 <= d <= 10^6).

    Output

    For each test case, output one line containing “Case #x: y”, where x is the test case number (standing from 1) and y is the maximum number of “pieces” modulo 10^9+7.

    Sample Input

    3

    3 3

    3 5

    4 5

    Sample Output

    Case #1: 8

    Case #2: 26

    Case #3: 31

    HINT

    Please use %lld when using long long

     

    题目大意就是n维空间切d刀,最多能分成几个部分。

    基本上通过推倒三位的大概就能很快推出整个的递推式。

    设p(n, d)表示n维空间切d刀。

    假设已经切了d-1刀,最后一刀自然切得越多越好。于是最后一刀如果和所有d-1到切的话自然是最好。但是可以逆过来看,相当于d-1到切最后一刀这个n-1维空间。

    于是p(n, d) = p(n, d-1) + p(n-1, d-1)

    然而这个式子虽然出来了,但是根据n和d的范围打表是不可能的。也不能直接暴力递推求解,自然考虑到可能要直接求表达式。

    然而,表达式我求了好久没求出来,不过看了最后表达式后,大概能有以下思路来求通项:

    首先有以下事实:

    1:手写打表的话:

    d->

    0

    1

    2

    3

    4

    5

    n

    1

    1

    2

    3

    4

    5

    6

    2

    1

    2

    4

    7

    11

    16

    3

    1

    2

    4

    8

    15

    26

    4

    1

    2

    4

    8

    16

    31

    5

    1

    2

    4

    8

    16

    32

    6

    1

    2

    4

    8

    16

    32

    会发现当n >= d时,通项是2^d,其实稍微考虑一下确实如此。因为第一列都是1,自然第二列从第二项开始都是2,同理往后从对角线往后都是乘2,自然是2^d。

    2:设p(n, d)的差数列为a(n, d)的话,

    自然a(n, d) = p(n, d) – p(n-1, d)

    由原式得p(n-1, d) = p(n-1, d-1) + p(n-2, d-1)

    三式式消去p得

    a(n, d) = a(n, d-1) + a(n-1, d-1)

    说明p的差数列也是满足这个递推式,同理p的任意k阶差数列都满足这个式子。

    然而让这些差数列最后通项不同的因素自然应该是前几项导致的

    有了上面两个结论,于是只用求n < d的情况,可以从下面两个角度考虑

    1:利用组合数式子:C(n, m) = C(n-1, m) + C(n-1, m-1),其中C(n, m)表示从n个中取m个。

    由于这个式子和题目递推式非常形似。 于是猜测C(n, m)为p的某一阶差数列。根据前几列和前几行的计算,C(n, m)为p的第一阶差数列。于是p(n, d) = sum(C(d, i)) (0 <= i <= n)

    2:根据第一个结论:列出第一阶的差数列

    d->

    0

    1

    2

    3

    4

    5

    n

    1

    1

    2

    3

    4

    5

    2

    0

    0

    1

    3

    6

    10

    3

    0

    0

    0

    1

    4

    10

    4

    0

    0

    0

    0

    1

    5

    5

    0

    0

    0

    0

    0

    1

    6

    0

    0

    0

    0

    0

    0

    基本上可以找规律,发现第一阶差数列是C(n, m)。

    然后就是求C(d, i)的和了,由于d很大,考虑C(d, i) = A(d, i) / i!,然后就是求分子和分母在模10^9+7的情况下的商了。自然需要考虑到逆元。

    这里对于逆元的处理可以预处理打表,经测试直接在线求exgcd逆元会T掉。

    这里预处理用了网上的一个神奇的递推式,还有一种是我大连海事一个同学的做法。

    代码(神奇式子):

    #include <iostream>
    #include <cstdio>
    #include <cstdlib>
    #include <cmath>
    #include <cstring>
    #include <algorithm>
    #define LL long long
    #define N 1000000007
    
    using namespace std;
    
    //快速幂
    //O(logn)
    LL quickPower(LL x, int n)
    {
        x %= N;
        LL a = 1;
        while (n)
        {
            a *= n&1 ? x : 1;
            a %= N;
            n >>= 1;
            x = (x*x) % N;
        }
        return a;
    }
    
    LL c[100005], a[100005], inv[100005];
    int n, d;
    
    void init()
    {
        //***预处理所有i在质数MOD下的逆元
        inv[1] = 1;
        for (int i = 2; i <= 100000; i++)
            inv[i] = inv[N%i]*(N-N/i) % N;
    
        a[0] = 1;
        for (int i = 1; i <= 100000; ++i)
            a[i] = (inv[i]*a[i-1]) % N;
    }
    
    void work()
    {
        if (n >= d)
        {
            printf("%lld
    ", quickPower(2, d));
            return;
        }
        LL now = d, ans = 0;
        c[0] = 1;
        for (int i = 1; i <= n; ++i)
        {
            c[i] = (now*c[i-1]) % N;
            now--;
        }
        for (int i = 0; i <= n; ++i)
        {
            ans += c[i]*a[i];
            ans %= N;
        }
        printf("%lld
    ", ans);
    }
    
    int main()
    {
        //freopen("test.in", "r", stdin);
        //freopen("test.out", "w", stdout);
        init();
        int T;
        scanf("%d", &T);
        for (int times = 1; times <= T; ++times)
        {
            printf("Case #%d: ", times);
            scanf("%d%d", &n, &d);
            work();
        }
        return 0;
    }
    View Code

    代码二(exgcd):

    #include <iostream>
    #include <cstdio>
    #include <cstring>
    #include <cmath>
    #include <algorithm>
    #include <queue>
    #define FOR(i,x,y)  for(int i = x;i < y;i ++)
    #define IFOR(i,x,y) for(int i = x;i > y;i --)
    #define ll  long long
    #define N   111111
    #define D   1111111
    #define MOD 1000000007
    
    using namespace std;
    
    ll c[N],mu[N];
    ll n,d;
    
    ll quickpow(ll a,ll n,ll m){
        ll ans=1;
        while(n){
            if(n&1) ans = (ans*a)%m;
            a = (a*a)%m;
            n>>=1;
        }
        return ans;
    }
    
    void ex_gcd(ll a,ll b,ll& d,ll& x,ll& y){
        if(!b)  {d = a;x = 1;y = 0;return;}
        ex_gcd(b,a%b,d,y,x);
        y -= x*(a/b);
    }
    
    ll inv(ll a,ll n){
        ll d,x,y;
        ex_gcd(a,n,d,x,y);
        return d == 1 ? (x+n)%n : -1;
    }
    
    void init(){
        FOR(i,1,N){
            mu[i] = inv(i,MOD);
        }
    }
    
    void C(){
        c[0] = 1;
        FOR(i,1,n+1){
            ll tem = (d+1-i)*mu[i]%MOD;
            c[i] = (tem*c[i-1]) % MOD;
        }
    }
    
    ll solve(){
        ll res = 0;
        FOR(i,0,n+1){
            res += c[i];
            res %= MOD;
        }
        return res;
    }
    
    int main()
    {
        //freopen("test.in","r",stdin);
        int t,tCase = 0;
        scanf("%d",&t);
        init();
        while(t--){
            printf("Case #%d: ",++tCase);
            scanf("%lld%lld",&n,&d);
            ll ans = 0;
            if(n >= d){
                ans = quickpow(2,d,MOD);
            }
            else{
                C();
                ans = solve();
            }
            printf("%lld
    ",ans);
        }
        return 0;
    }
    View Code
  • 相关阅读:
    2020.10.25【NOIP提高A组】模拟 总结
    6831. 2020.10.24【NOIP提高A组】T1.lover
    枚举一个数$n$的所有质因子
    gmoj 6832. 2020.10.24【NOIP提高A组】T2.world
    2020.10.24【NOIP提高A组】模拟 总结
    2020.10.17【NOIP提高A组】模拟 总结
    jQuery EasyUI Portal 保存拖动位置,仿谷歌DashBoard效果的
    SQLMAP注入教程-11种常见SQLMAP使用方法详解
    Windows下sqlmap的安装图解
    swap file "*.swp" already exists!的解决方法
  • 原文地址:https://www.cnblogs.com/andyqsmart/p/4571352.html
Copyright © 2011-2022 走看看