Lưu ý: Các code mẫu dưới đây chỉ mang tính tham khảo và có thể không AC được bài tập này

Code mẫu của happyboy99x

#include<bits/stdc++.h>
using namespace std;

class SumOfPowers {
    typedef vector<long long> Row;
    typedef vector<Row> Matrix;

    static const int MOD = 1e9 + 7;

    Matrix multiply(const Matrix &a, const Matrix &b) {
        int m = a.size(), n = a[0].size(), p = b[0].size();
        Matrix res (m, Row(p, 0));
        for(int i = 0; i < m; ++i)
            for(int j = 0; j < p; ++j)
                for(int k = 0; k < n; ++k)
                    res[i][j] = (res[i][j] + a[i][k] * b[k][j]) % MOD;
        return res;
    }

    Matrix power(Matrix a, int n) {
        Matrix res (a.size(), Row(a.size(), 0));
        for(unsigned i = 0; i < res.size(); ++i) res[i][i] = 1;
        for(; n > 0; n >>= 1) {
            if(n & 1) res = multiply(res, a);
            a = multiply(a, a);
        }
        return res;
    }

public:
    int value(int n, int k) {
        Matrix a (1, Row(k + 2, 0)); a[0][0] = 1;
        Matrix p (k + 2, Row(k + 2, 0));
        p[0][0] = 1;
        for(int j = 1; j <= k; ++j) {
            p[0][j] = 1;
            for(int i = 1; i <= j; ++i)
                p[i][j] = (p[i-1][j-1] + p[i][j-1]) % MOD;
        }
        for(int i = 0; i <= k; ++i) p[i][k+1] = p[i][k];
        p[k+1][k+1] = 1;
        return multiply(a, power(p, n)).front().back();
    }
};

int main() {
    ios::sync_with_stdio(false);
    for(int n, k; cin >> n >> k; cout << SumOfPowers().value(n, k) << '\n');
    return 0;
}

Code mẫu của ladpro98

#include <cstring>
#include <vector>
#include <list>
#include <map>
#include <set>
#include <queue>
#include <deque>
#include <stack>
#include <bitset>
#include <algorithm>
#include <functional>
#include <numeric>
#include <utility>
#include <sstream>
#include <iostream>
#include <iomanip>
#include <cstdio>
#include <cmath>
#include <climits>
#include <cstdlib>
#include <ctime>
#include <memory.h>
#include <cassert>
#define FOR(i, a, b) for(int i = (a); i < (b); i++)
#define REP(i, a, b) for(int i = (a); i <=(b); i++)
#define FORD(i, a, b) for(int i = (a); i > (b); i--)
#define REPD(i, a, b) for(int i = (a); i >=(b); i--)
#define SZ(a) (int((a).size()))
#define ALL(a) (a).begin(), (a).end()
#define PB push_back
#define MP make_pair
#define LL long long
#define LD long double
#define II pair<int, int>
#define X first
#define Y second
#define VI vector<int>
const int D = 53;
const int MOD = 1000000007;
using namespace std;
int n, k, dim;

struct mat {LL a[D][D];};
mat base;

mat operator * (mat a, mat b) {
    mat c;
    FOR(i, 0, dim) FOR(j, 0, dim) {
        c.a[i][j] = 0;
        FOR(k, 0, dim)
            c.a[i][j] = (c.a[i][j] + a.a[i][k] * b.a[k][j]) % MOD;
    }
    return c;
}

void setDim() {
    memset(base.a, 0, sizeof base.a);
    dim = k + 2;
    base.a[0][0] = base.a[0][k + 1] = 1;
    REP(i, 1, k + 1) {
        base.a[i][1] = base.a[i][i] = 1;
        FOR(j, 2, i)
            base.a[i][j] = (base.a[i - 1][j] + base.a[i - 1][j - 1]) % MOD;
    }
}

int main() {
    ios :: sync_with_stdio(0); cin.tie(0);
    while (cin >> n >> k) {
        setDim();
        mat a = base;
        VI bin;
        while (n) {
            bin.PB(n & 1);
            n >>= 1;
        }
        REPD(i, SZ(bin) - 2, 0) {
            a = a * a;
            if (bin[i]) a = a * base;
        }
        int ans = 0;
        REP(i, 1, k + 1) ans = (ans + a.a[0][i]) % MOD;
        cout << ans << '\n';
    }
    return 0;
}

Code mẫu của RR

#include <iostream>
#include <algorithm>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <string>
#include <vector>
#include <queue>
#include <stack>
#include <map>
#include <set>
#include <iomanip>
#include <bitset>
#include <complex>

#define FOR(i,a,b) for(int i = a; i <= b; ++i)
#define FORD(i,a,b) for(int i = a; i >= b; --i)
#define REP(i,a) for(int i = 0; i < a; ++i)
#define MP make_pair
#define PB push_back
using namespace std;

const long long BASE = 1000000007LL;

struct Poly {
    int deg;
    long long a[55];

    Poly multX() {
        Poly res;
        res.deg = deg + 1;
        memset(res.a, 0, sizeof res.a);
        FORD(i,res.deg,1)
            res.a[i] = a[i-1];
        return res;
    }

    Poly operator + (Poly p) {
        Poly res;
        res.deg = max(p.deg, deg);
        memset(res.a, 0, sizeof res.a);

        FOR(i,0,res.deg) {
            res.a[i] = (p.a[i] + a[i]) % BASE;
        }
        return res;
    }

    Poly operator - (Poly p) {
        Poly res;
        res.deg = max(p.deg, deg);
        memset(res.a, 0, sizeof res.a);

        FOR(i,0,res.deg) {
            res.a[i] = (a[i] + BASE - p.a[i]) % BASE;
        }
        return res;
    }

    Poly operator * (long long x) {
        Poly res;
        res.deg = deg;
        memset(res.a, 0, sizeof res.a);

        FOR(i,0,res.deg) {
            res.a[i] = (a[i] * x) % BASE;
        }
        return res;
    }

    void print() {
        FOR(i,0,deg) cout << a[i] << ' '; cout << endl;
    }
} f[55];

long long power(long long x, int k) {
    if (k == 0) return 1;
    if (k == 1) return x % BASE;
    long long mid = power(x, k >> 1);
    mid = (mid * mid) % BASE;
    if (k & 1) return (mid * x) % BASE;
    else return mid;
}

long long inverse(long long x) {
    return power(x, BASE - 2);
}

void init() {
    f[0].deg = 1;
    f[0].a[0] = 0;
    f[0].a[1] = 1;

    f[1].deg = 2;
    f[1].a[0] = 0;
    f[1].a[1] = inverse(2);
    f[1].a[2] = inverse(2);

    FOR(k,2,50) {
        f[k] = f[k-1].multX() + f[k-1];
        FOR(i,0,k-1) {
            f[k] = f[k] - f[i] * f[k-1].a[i];
        }
        f[k] = f[k] * inverse((f[k-1].a[k] + 1) % BASE);
    }
}

long long n, k;

int main() {
    init();
    while (cin >> n >> k) {
        long long res = 0;
        long long pN = 1;
        FOR(i,0,k+1) {
            res = (res + f[k].a[i] * pN) % BASE;
            pN = (n * pN) % BASE;
        }
        cout << res << endl;
    }
    return 0;
}

Code mẫu của hieult

#include<cstdio>
#include<cmath>
#include<math.h>
#include<cstring>
#include<cstdlib>
#include<cassert>
#include<ctime>
#include<algorithm>
#include<iterator>
#include<iostream>
#include<cctype>
#include<string>
#include<vector>
#include<map>
#include<set>
#include<queue>
#include<list>
#define ep 0.00001
#define maxn 100111
#define oo 1111111111
#define modunlo 35000
#define mod 1000000007
#define TR(c, it) for(typeof((c).begin()) it=(c).begin(); it!=(c).end(); it++)
//#define g 9.81
double const PI=4*atan(1.0);

//#include<conio.h>

using namespace std;

typedef pair<int, int> II;
typedef vector<int> VI;
typedef vector<II> VII;
typedef vector<VI> VVI;
typedef vector<VII> VVII;

long long mu(long long x,long long n){
      if(n==1) return x;
      long long t = mu(x,n/2);
      if(n&1) return (((t*t)%mod)*x)%mod;
      return (t*t)%mod;
};

long long f[35005][51],C[51][51],tong[35005][51],n,k;

int main(){
    // freopen("C11CAL.in","r",stdin);
    // freopen("C11CAL.out","w",stdout);   
     memset(C,0,sizeof(C));
     memset(tong,0,sizeof(tong));

     for(int i = 1;i<=35000;i++){
          long long tinh = 1;   
          for(int j = 0;j<=50;j++){
               f[i][j] = tinh;
               tong[i][j] = (tong[i-1][j]+tinh)%mod;
               tinh = (tinh*i)%mod; 
          }
     }

     for(int i = 0;i<=50;i++) C[0][i] = 1;
     for(int i = 1;i<=50;i++) for(int j = 1;j<=i;j++) C[j][i] = (C[j][i-1]+C[j-1][i-1])%mod;

     while(cin>>n>>k){
          long long u = n/modunlo - 1;
          long long v = n%modunlo;
          long long kq = 0;
          kq = (tong[modunlo][k]*(u+1))%mod;
          if(u>0)
          for(int i = k;i>=1;i--){
               kq = (kq+(((((tong[u][i]*C[i][k])%mod)*f[modunlo][i])%mod)*tong[modunlo][k-i])%mod)%mod;
          }
          kq = (kq+tong[v][k])%mod;
          for(int i = k;i>=1;i--){
               kq = (kq+(((((f[u+1][i]*C[i][k])%mod)*f[modunlo][i])%mod)*tong[v][k-i])%mod)%mod;
          }
          printf("%lld\n",kq);        
     }

}

Code mẫu của skyvn97

#include<cstdio>
#define MAX   55
typedef long long ll;
const ll MOD=1e9+7;
struct matrix {
    int m,n;
    ll d[MAX][MAX];
    matrix(){}
    matrix operator * (const matrix &x) {
        matrix res (*this);
        int a=m;
        int b=n;
        int c=x.n;
        res.m=a;
        res.n=c;
        int i,j,k;
        for (i=0;i<a;i=i+1)
            for (j=0;j<c;j=j+1) {
                res.d[i][j]=0;
                for (k=0;k<b;k=k+1)
                    res.d[i][j]=(res.d[i][j]+d[i][k]*x.d[k][j])%MOD;
            }
        return (res);
    }
    matrix operator ^ (const int &k) {
        if (k==1) return (*this);
        matrix r=(*this)^(k/2);
        r=r*r;
        if (k%2==1) r=r*(*this);
        return (r);
    }
};
matrix fst,mul;
int k,n;
ll c[MAX][MAX];
void finit(void) {
    int i,j;
    c[0][0]=1;
    for (i=1;i<=50;i=i+1) {
        c[0][i]=1;
        c[i][0]=0;
    }
    for (i=1;i<=51;i=i+1)
        for (j=1;j<=51;j=j+1) {
            if (i>j) c[i][j]=0;
            if (i==j) c[i][j]=1;
            if (i<j) c[i][j]=(c[i-1][j-1]+c[i][j-1])%MOD;
        }
}
void construct(void) {
    mul.m=k+2;
    mul.n=k+2;
    int i,j;
    for (i=0;i<k+1;i=i+1)
        for (j=0;j<k+1;j=j+1)
            mul.d[i][j]=c[i][j];
    for (i=0;i<k+1;i=i+1) mul.d[k+1][i]=0;
    for (i=0;i<k;i=i+1) mul.d[i][k+1]=0;
    for (i=0;i<2;i=i+1) mul.d[k+i][k+1]=1;
}
void process(void) {
    scanf("%d",&k);
    construct();
    fst.m=1;
    fst.n=k+2;
    int i;
    for (i=0;i<k+1;i=i+1) fst.d[0][i]=1;
    fst.d[0][k+1]=0;
    if (n>1) fst=fst*(mul^(n-1));   
    printf("%lld\n",(fst.d[0][k]+fst.d[0][k+1])%MOD);
}
int main(void) {
    finit();
    while (scanf("%d",&n)==1) process();
    return 0;

}