杜教筛BM(找规律)

简介: 杜教筛BM(找规律)

代码来自学长

#include<bits/stdc++.h>
using namespace std;
#define rep(i,a,n) for (int i=a;i<n;i++)
#define pb push_back
typedef long long ll;
#define SZ(x) ((ll)(x).size())
typedef vector<ll> VI;
typedef pair<ll,ll> PII;
const ll mod=1000000007;
ll powmod(ll a,ll b) {ll res=1;a%=mod; assert(b>=0); for(;b;b>>=1){if(b&1)res=res*a%mod;a=a*a%mod;}return res;}
// head
ll _,n;
namespace linear_seq {
    const ll N=10010;
    ll res[N],base[N],_c[N],_md[N];
    vector<ll> Md;
    void mul(ll *a,ll *b,ll k) {
        rep(i,0,k+k) _c[i]=0;
        rep(i,0,k) if (a[i]) rep(j,0,k) _c[i+j]=(_c[i+j]+a[i]*b[j])%mod;
        for (ll i=k+k-1;i>=k;i--) if (_c[i])
            rep(j,0,SZ(Md)) _c[i-k+Md[j]]=(_c[i-k+Md[j]]-_c[i]*_md[Md[j]])%mod;
        rep(i,0,k) a[i]=_c[i];
    }
    ll solve(ll n,VI a,VI b) { 
        ll ans=0,pnt=0;
        ll k=SZ(a);
        assert(SZ(a)==SZ(b));
        rep(i,0,k) _md[k-1-i]=-a[i];_md[k]=1;
        Md.clear();
        rep(i,0,k) if (_md[i]!=0) Md.push_back(i);
        rep(i,0,k) res[i]=base[i]=0;
        res[0]=1;
        while ((1ll<<pnt)<=n) pnt++;
        for (ll p=pnt;p>=0;p--) {
            mul(res,res,k);
            if ((n>>p)&1) {
                for (ll i=k-1;i>=0;i--) res[i+1]=res[i];res[0]=0;
                rep(j,0,SZ(Md)) res[Md[j]]=(res[Md[j]]-res[k]*_md[Md[j]])%mod;
            }
        }
        rep(i,0,k) ans=(ans+res[i]*b[i])%mod;
        if (ans<0) ans+=mod;
        return ans;
    }
    VI BM(VI s) {
        VI C(1,1),B(1,1);
        ll L=0,m=1,b=1;
        rep(n,0,SZ(s)) {
            ll d=0;
            rep(i,0,L+1) d=(d+(ll)C[i]*s[n-i])%mod;
            if (d==0) ++m;
            else if (2*L<=n) {
                VI T=C;
                ll c=mod-d*powmod(b,mod-2)%mod;
                while (SZ(C)<SZ(B)+m) C.pb(0);
                rep(i,0,SZ(B)) C[i+m]=(C[i+m]+c*B[i])%mod;
                L=n+1-L; B=T; b=d; m=1;
            } else {
                ll c=mod-d*powmod(b,mod-2)%mod;
                while (SZ(C)<SZ(B)+m) C.pb(0);
                rep(i,0,SZ(B)) C[i+m]=(C[i+m]+c*B[i])%mod;
                ++m;
            }
        }
        return C;
    }
    ll gao(VI a,ll n) {
        VI c=BM(a);
        c.erase(c.begin());
        rep(i,0,SZ(c)) c[i]=(mod-c[i])%mod;
        return solve(n,c,VI(a.begin(),a.begin()+SZ(c)));
    }
};
int main() {
    int T;
    scanf("%d",&T);
    while (T--)
    {
        scanf("%lld",&n);
        vector<ll>v;
        v.push_back(3);
        v.push_back(9);
        v.push_back(20);
        v.push_back(46);
        v.push_back(106);
        v.push_back(244);
        v.push_back(560);
        v.push_back(1286);
        v.push_back(2956);
        v.push_back(6794);
        printf("%lld\n",linear_seq::gao(v,n-1));
    }
}
目录
相关文章
|
8月前
|
算法 安全 测试技术
[组合数学]LeetCode:2954:统计感冒序列的数目
[组合数学]LeetCode:2954:统计感冒序列的数目
|
机器学习/深度学习 算法
【算法基础】筛质数
【算法基础】筛质数
74 0
|
算法
初阶OI素数算法——埃拉托尼斯筛
时间复杂度比较优秀且易于理解的素数筛选法
92 0
|
7月前
|
机器学习/深度学习 资源调度 算法
杜教筛
【6月更文挑战第9天】
32 2
|
存储
欧拉筛&&埃氏筛
欧拉筛&&埃氏筛
107 0
|
算法
质数筛法:朴素素数筛,埃氏筛,欧式筛
质数筛法:朴素素数筛,埃氏筛,欧式筛
173 0
|
算法 C++
容斥原理算法的实现
容斥原理算法的实现
容斥原理算法的实现
|
存储
【每日一题Day95】LC1815得到新鲜甜甜圈的最多组数 | 状态压缩dp 记忆化搜索
子问题、哪些操作会影响数据:余下的甜甜圈数量left,以及剩余可以选的元素个数 cnt[x]【dfs函数的两个参数->使用状态压缩至一个int类型变量中】
124 0
【每日一题Day95】LC1815得到新鲜甜甜圈的最多组数 | 状态压缩dp 记忆化搜索
【每日一题Day63】LC1753移除石子的最大得分 | 贪心 + 数学
当对a,b,c进行升序排序后,如果a+b≥c,应先优先让a 和b去匹配,匹配的次数为a+b−c,得分为count=⌈(a+b−c)/2⌉
104 0
【CCCC】L3-018 森森美图 (30分),计算几何+判断三点共线+bfs最短路
【CCCC】L3-018 森森美图 (30分),计算几何+判断三点共线+bfs最短路
154 0