Codeforces 717A Festival Organization(组合数学：斯特林数+Fibonacci数列+推公式)

Codeforces 717A Festival Organization(组合数学：斯特林数+Fibonacci数列+推公式)

牛逼题。。。。。推公式非常的爽。。。虽然我是看了别人的博客才推出来的。。。

0.1 斯特林数

下面要用到的是带符号的第一类斯特林数。

\(x^{n\downarrow}=\prod_{i=0}^{n-1}(x-i)=\sum_{k=0}^ns(n,k)x^k\)

有递推公式\(s(n,m)=s(n-1,m-1)-(n-1)*s(n-1,m)\)

0.2 斐波那契数列的通项公式

\(f_i=\frac{1}{\sqrt5}[(\frac{1+\sqrt5}{2})^i-(\frac{1-\sqrt5}{2})^i]\)

1 求出长度为\(n\)的合法0/1序列数

先考虑dp做法，\(a_{i,0}\)表长度为\(i\),结尾为\(1\)的方案数，\(a_{i,1}\)同理。

很明显有\(a_{0,0}=a_{0,1}=1,a_{i,0}=a_{i-1,1},a_{i,1}=a_{i-1,0}+a_{i-1,1}\)

数学归纳法很轻松就能证出\(a_{i,0}=F_{i},a_{i,1}=F_{i+1}\),则长度为\(n\)的方案数为\(F_{i+2}\),其中\(F_i\)表斐波那契数列的第\(i\)项\((F_1=F_2=1)\)

2 \(Ans,S_n\)

\(Ans=\sum_{i=l}^r C_{F_{i+2}}^k\)

记\(S_n=\sum_{i=1}^nC_{F_{i+2}}^k=\sum_{i=3}^{n+2}C_{F_i}^k\)

则有\(Ans=S_r-S_{l-1}\)

3 化简

\(S_n\)

\(=\frac{1}{k!}\sum_{i=3}^{n+2}\prod_{j=0}^{k-1}(F_i-j)\)

\(=\frac{1}{k!}\sum_{i=3}^{n+2}\sum_{j=0}^ks(k,j)F_i^j\)

\(=\frac{1}{k!}\sum_{j=0}^k\sum_{i=3}^{n+2}s(k,j)F_i^j\)

\(=\frac{1}{k!}\sum_{j=0}^ks(k,j)\sum_{i=3}^{n+2}F_i^j\)

\(=\frac{1}{k!}\sum_{j=0}^ks(k,j)\sum_{i=3}^{n+2}(\frac{1}{\sqrt5}[a^i-b^i])^j\) 其中\(a=\frac{1+\sqrt5}{2},b=\frac{1-\sqrt5}{2}\)

\(=\frac{1}{k!}\sum_{j=0}^ks(k,j)\sum_{i=3}^{n+2}\sum_{p=0}^j(\frac{1}{5}\sqrt5)^jC_j^pa^{ip}(-1)^{j-p}b^{i(j-p)}\)

\(=\frac{1}{k!}\sum_{j=0}^ks(k,j)\sum_{p=0}^j(\frac{1}{5}\sqrt5)^jC_j^p(-1)^{j-p}\sum_{i=3}^{n+2}a^{ip}b^{i(j-p)}\)

\(=\frac{1}{k!}\sum_{j=0}^ks(k,j)\sum_{p=0}^j(\frac{1}{5}\sqrt5)^jC_j^p(-1)^{j-p}\sum_{i=3}^{n+2}(a^pb^{j-p})^i\)

\(=\frac{1}{k!}\sum_{j=0}^ks(k,j)\sum_{p=0}^j(\frac{1}{5}\sqrt5)^jC_j^p(-1^{j-p})*T(j,p)\)

其中,令\(q(j,p)=a^pb^{j-p}\)

\[T(j,p)=\begin{cases} n & q(j,p)=1 \\ \frac{q^3(j,p)(1-q^n(j,p))}{1-q(j,p)} & q(j,p)\not=1 \end{cases}\]

枚举\(j,p\)就能得到答案了。

4 代码

#include <bits/stdc++.h>
#define pii pair<long long,long long>
#define LL long long
#define x first
#define y second
#define MAXN 7000
using namespace std;
const LL mod = 1000000007;
LL k, L, R;
LL qp(LL a, LL b)
{
    LL ret = 1;
    while (b)
    {
        if (b & 1)
            ret = ret * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return ret;
}


LL Inv[208];
LL inv(LL a)
{
    if (a > 200)
        return qp(a, mod - 2);
    else if (!Inv[a])
        return Inv[a] = qp(a, mod - 2);
    else
        return Inv[a];
}

LL CC[208][208];
LL S[208][208];
void init()
{
    CC[0][0] = 1;
    for (LL i = 1; i <= 200; ++i)
    {
        CC[i][0] = 1;
        for (LL j = 1; j <= i; ++j)
        {
            CC[i][j] = CC[i][j - 1] * inv(j) % mod * (i - j + 1) % mod;
        }
    }
    S[0][0] = 1;
    for (LL i = 1; i <= 200; ++i)
    {
        S[i][0] = 0;
        for (LL j = 1; j < i; ++j)
        {
            S[i][j] = ((S[i - 1][j - 1] - (i - 1) * S[i - 1][j] % mod) % mod + mod) % mod;
        }
        S[i][i] = 1;
    }
}



struct tls
{
    LL a, b;
    tls(){}
    tls(LL _a, LL _b){a = (_a % mod + mod) % mod; b = (_b % mod + mod) % mod;}
    tls operator+(const tls& t) const
    {
        return tls((a + t.a) % mod, (b + t.b) % mod);
    }
    tls operator-(const tls& t) const
    {
        return tls((a - t.a + mod) % mod, (b - t.b + mod) % mod);
    }
    tls operator*(const tls& t) const
    {
        return tls((a*t.a%mod + 5*b%mod*t.b%mod) % mod, (a*t.b%mod + b*t.a%mod) % mod);
    }
    tls operator/(const tls& t) const
    {
        LL r = inv((t.a*t.a%mod - 5*t.b%mod*t.b%mod + mod) % mod);
        if ((t.a*t.a%mod - 5*t.b%mod*t.b%mod + mod) % mod == 0)
            throw;
        return tls(r * (((a*t.a%mod) - (b*t.b%mod*5%mod) + mod) % mod) % mod, (r * (((t.a*b%mod) - (a*t.b%mod) + mod) % mod) % mod));
    }
};

tls qp(tls a, LL b)
{
    tls ret(1, 0);
    while (b)
    {
        if (b & 1)
            ret = ret * a;
        a = a * a;
        b >>= 1;
    }
    return ret;
}

ostream& operator<<(ostream& out, tls p)
{
    out << p.a << "+" << p.b << "sqrt(5)";
    return out;
}

LL gao(LL n, LL k)
{
    tls a(inv(2)%mod, inv(2)%mod), b(inv(2)%mod,(mod-inv(2)%mod)%mod);
    tls ret(0, 0);
    for (LL j = 0; j <= k; ++j)
    {
        tls tot(0, 0);
        for (LL p = 0; p <= j; ++p)
        {
            tls temp(1, 0);
            tls q = qp(a/b, p) * qp(b, j);
            if (q.a == 1 && q.b == 0)
            {
                temp = temp * tls(n, 0);
            }
            else
            {
                temp = temp * qp(q, 3) * (tls(1, 0) - qp(q, n));
                temp = temp / (tls(1, 0) - q);
            }
            temp = temp * tls(qp(mod - 1, j - p), 0) * tls(CC[j][p], 0);

            tot = tot + temp;
        }
        tot = tot * tls(S[k][j], 0) * qp(tls(0, inv(5)), j);
        ret = ret + tot;
    }
    for (LL i = 1; i <= k; ++i)
    {
        ret = ret * tls(inv(i), 0);
    }
    return ret.a;
}


int main()
{
    /*cout << inv(25) << endl;
    cout << qp(tls(0, inv(5)), 4) << endl;*/
    init();
    cin >> k >> L >> R;
    cout << (gao(R, k) - gao(L - 1, k) + mod) % mod << endl;
    return 0;
}