CF717A Festival Organization

https://www.luogu.com.cn/problem/CF717A
先把斐波那契数列的通项公式写出来，然后推式子
先把组合数展开，然后把下降幂再用第一类斯特林数展开

发现5在模意义下没有二次剩余，所以扩个域就好了

code:

#include<bits/stdc++.h>
#define mod 1000000007
#define ll long long
#define N 505
using namespace std;
int add(int x, int y) { x += y;
    if(x >= mod) x -= mod;
    return x;
}
int sub(int x, int y) { x -= y;
    if(x < 0) x += mod;
    return x;
}
int mul(int x, int y) {
    return 1ll * x * y % mod;
}
int qpow(int x, ll y) {
    int ret = 1;
    for(; y; y >>= 1, x = mul(x, x)) if(y & 1) ret = mul(ret, x);
    return ret;
}
const int inv2 = 500000004;
const int inv5 = 400000003;
struct cp {
    int x, y;
    cp(const int x = 0, const int y = 0) : x(x), y(y) {};
    cp operator + (const cp &o) const { return cp(add(x, o.x), add(y, o.y));}
    cp operator - (const cp &o) const { return cp(sub(x, o.x), sub(y, o.y));}
    cp operator * (const cp &o) const { return cp(add(mul(x, o.x), mul(5, mul(y, o.y))), add(mul(x, o.y), mul(y, o.x)));}
    cp inv() const {
        cp ret = cp(x, (mod - y) % mod);
        int d = qpow(sub(mul(x, x), mul(5, mul(y, y))), mod - 2);
        return ret * d;
    }
    cp operator / (const cp &o) const { return *this * o.inv();}
} a(0, inv5), b(inv2, inv2), c (0, mod - inv5), d(inv2, mod - inv2);
cp qpoww(cp x, ll y) {
    cp ret(1, 0);
    for(; y; y >>= 1, x = x * x) if(y & 1) ret = ret * x;
    return ret;
}
cp sum(cp x, ll l, ll r) {
    if(x.x == 1 && !x.y) return (r - l + 1) % mod;
    return (qpoww(x, r + 1) - qpoww(x, l)) / (x - 1);
}
int n, C[N][N], s[N][N];
ll l, r;
int main() {
    scanf("%d%lld%lld", &n, &l, &r); l += 2, r += 2;
    for(int i = 0; i <= n; i ++) C[i][0] = 1;
    C[1][1] = s[1][1] = 1;
    for(int i = 2; i <= n; i ++)
        for(int j = 1; j <= i; j ++)
            C[i][j] = add(C[i - 1][j], C[i - 1][j - 1]),
            s[i][j] = add(s[i - 1][j - 1], mul(i - 1, s[i - 1][j]));
    
    
    int ans = 0;
    for(int i = 0; i <= n; i ++) {
        int ret = 0;
        for(int j = 0; j <= i; j ++) {
            cp fn = qpoww(a, j) * qpoww(c, i - j), ha = qpoww(b, j) * qpoww(d, i - j);
            //printf("(%d %d)   (%d %d) \n", fn.x, fn.y, ha.x, ha.y);
            fn = fn * sum(ha, l, r);
            ret = add(ret, mul(C[i][j], fn.x)); //printf("%d\n", mul(1, fn.x));
        }
    //    printf("%d ", ret);
        ret = mul(ret, s[n][i]);
        
        if((n - i) & 1) ret = (mod - ret) % mod;
        ans = add(ans, ret);
    }
    int ninv = 1;
    for(int i = 1; i <= n; i ++) ninv = mul(ninv, i);
    ninv = qpow(ninv, mod - 2);
    printf("%d", mul(ans, ninv));
    return 0;
}