求 A249665 前 \(n\) 项的 \(m\) 次方和，对 \(10^9 + 7\) 取模。
\(1 \leq n \leq 10^{18}\)，\(1 \leq m \leq 3\)。
我们可以暴力枚举 \(1\) 到 \(n\) 的排列，由此求出 \(a_n\)。
经过打表，可以得到 \(A\) 的前几项为 \(1,1,1,2,6,14,28,56,118,254,541,1140,2401,5074,10738\)。
由此我们可以得到规律：\(a_n = -a_{n-8}-a_{n-7}+a_{n-5}+a_{n-4}+2 \times a_{n-3}-a_{n-2}+2 \times a_{n-1}\)。
于是可以用矩阵快速幂求出 \(1\) 次方和，对于 \(2\) 次方和 \(3\) 次方可以打更多的表，最终会得到一百多项的递推式。
#include<iostream>
#include<cstdio>
#include<algorithm>
using namespace std;
const long long mod=1000000007;
const int K1=8,K2=36,K3=120;
long long k_1[K1+1]={0,2,1000000006,2,1,1,0,1000000006,1000000006};
long long k_2[K2+1]={0,3,2,11,36,57,999999978,999999697,999999268,999999713,1022,1055,999999922,428,708,999999575,641,1883,1536,577,999999772,999999093,999998982,999999574,999999841,114,31,0,999999993,75,17,13,1000000005,1000000005,0,1,1000000006};
long long k_3[K3+1]={0,6,9,119,776,2928,999999866,999947231,999685929,999532233,1006034,3403103,3091675,18854028,41939739,994839435,123907888,686578570,538510624,463875551,593907298,775112621,862799554,826799974,323568509,272946549,865736932,849802794,311419209,805100035,339107278,456958980,305935140,337495512,230577161,989563382,512084924,610163742,175186535,153394145,542690041,50533570,340316201,88175161,908013867,593526711,370645361,270223951,599117789,984653980,300238059,940672603,619708604,241513668,555065433,447221261,517442267,369809014,727409057,654158750,399582554,676918142,975097682,435015557,940699679,119850330,454035523,899881054,501312719,213643333,937627741,348520370,422818051,912039075,988291205,627079501,197735216,763988949,402495671,193045517,116881412,386288279,506736426,560043784,134890183,826753538,554299933,726711648,124505105,979617265,769179718,964951092,201335067,391976488,964735448,16230894,520533206,281250193,800363870,100229687,972206749,12790510,994869596,1367679,999105328,1188015,999184856,110016,999999173,3761,162,1064,999998446,2081,999999729,7,6,1000000006,1000000005,0,1000000006};
long long f1[K3+1]={0,1,1,1,2,6,14,28,56},f2[K3+1],f3[K3+1],sf1[K3+1],sf2[K3+1],sf3[K3+1];
void init(){for(int i=9;i<=K3;i++){for(int j=1;j<=K1;j++){f1[i]+=f1[i-j]*k_1[j]%mod;f1[i]%=mod;}}for(int i=1;i<=K3;i++){f2[i]=f1[i]*f1[i]%mod;f3[i]=f1[i]*f1[i]%mod*f1[i]%mod;sf1[i]=(sf1[i-1]+f1[i])%mod;sf2[i]=(sf2[i-1]+f2[i])%mod;sf3[i]=(sf3[i-1]+f3[i])%mod;}
}
long long A[K3+1][K3+1],I[K3+1][K3+1],T[K3+1][K3+1];
long long solve1(long long n){if(n<=K1){return sf1[n];}for(int i=1;i<K1;i++){A[i][i-1]=1;}for(int i=1;i<=K1;i++){A[K1-i][K1-1]=k_1[i];A[K1-i][K1]=k_1[i];}A[K1][K1]=1;for(int i=0;i<=K1;i++){I[i][i]=1;}long long pow_num=n-K1;while(pow_num){if(pow_num&1){for(int i=0;i<=K1;i++){for(int j=0;j<=K1;j++){T[i][j]=0;for(int k=0;k<=K1;k++){T[i][j]+=I[i][k]*A[k][j]%mod;}T[i][j]%=mod;}}for(int i=0;i<=K1;i++){for(int j=0;j<=K1;j++){I[i][j]=T[i][j];}}}for(int i=0;i<=K1;i++){for(int j=0;j<=K1;j++){T[i][j]=0;for(int k=0;k<=K1;k++){T[i][j]+=A[i][k]*A[k][j]%mod;}T[i][j]%=mod;}}for(int i=0;i<=K1;i++){for(int j=0;j<=K1;j++){A[i][j]=T[i][j];}}pow_num>>=1;}long long ans=0;for(int i=0;i<K1;i++){ans+=f1[i+1]*I[i][K1]%mod;ans%=mod;}ans+=sf1[K1]*I[K1][K1]%mod;ans%=mod;return ans;
}
long long solve2(long long n){if(n<=K2){return sf2[n];}for(int i=1;i<K2;i++){A[i][i-1]=1;}for(int i=1;i<=K2;i++){A[K2-i][K2-1]=k_2[i];A[K2-i][K2]=k_2[i];}A[K2][K2]=1;for(int i=0;i<=K2;i++){I[i][i]=1;}long long pow_num=n-K2;while(pow_num){if(pow_num&1){for(int i=0;i<=K2;i++){for(int j=0;j<=K2;j++){T[i][j]=0;for(int k=0;k<=K2;k++){T[i][j]+=I[i][k]*A[k][j]%mod;}T[i][j]%=mod;}}for(int i=0;i<=K2;i++){for(int j=0;j<=K2;j++){I[i][j]=T[i][j];}}}for(int i=0;i<=K2;i++){for(int j=0;j<=K2;j++){T[i][j]=0;for(int k=0;k<=K2;k++){T[i][j]+=A[i][k]*A[k][j]%mod;}T[i][j]%=mod;}}for(int i=0;i<=K2;i++){for(int j=0;j<=K2;j++){A[i][j]=T[i][j];}}pow_num>>=1;}long long ans=0;for(int i=0;i<K2;i++){ans+=f2[i+1]*I[i][K2]%mod;ans%=mod;}ans+=sf2[K2]*I[K2][K2]%mod;ans%=mod;return ans;
}
long long solve3(long long n){if(n<=K3){return sf3[n];}for(int i=1;i<K3;i++){A[i][i-1]=1;}for(int i=1;i<=K3;i++){A[K3-i][K3-1]=k_3[i];A[K3-i][K3]=k_3[i];}A[K3][K3]=1;for(int i=0;i<=K3;i++){I[i][i]=1;}long long pow_num=n-K3;while(pow_num){if(pow_num&1){for(int i=0;i<=K3;i++){for(int j=0;j<=K3;j++){T[i][j]=0;for(int k=0;k<=K3;k++){T[i][j]+=I[i][k]*A[k][j]%mod;}T[i][j]%=mod;}}for(int i=0;i<=K3;i++){for(int j=0;j<=K3;j++){I[i][j]=T[i][j];}}}for(int i=0;i<=K3;i++){for(int j=0;j<=K3;j++){T[i][j]=0;for(int k=0;k<=K3;k++){T[i][j]+=A[i][k]*A[k][j]%mod;}T[i][j]%=mod;}}for(int i=0;i<=K3;i++){for(int j=0;j<=K3;j++){A[i][j]=T[i][j];}}pow_num>>=1;}long long ans=0;for(int i=0;i<K3;i++){ans+=f3[i+1]*I[i][K3]%mod;ans%=mod;}ans+=sf3[K3]*I[K3][K3]%mod;ans%=mod;return ans;
}
int main(){init();long long n;int m;scanf("%lld %d",&n,&m);if(m==1){printf("%lld",solve1(n));}else if(m==2){printf("%lld",solve2(n));}else{printf("%lld",solve3(n));}return 0;
}
