题目
点这里看题目。
分析
首先可以娴熟地推倒一发式子：
[egin{split}&f_k(n)&overset{mathrm{def}}{=}f(k,n)\Rightarrow &f_k(n)&=sum_{i=1}^{n-1}f_k(i)+n^k\end{split}
]
定义(S_k(n)=sum_{i=1}^{n}f_k(i))，做差计算可以得到：
[egin{split}Rightarrow &&S_k(n)&&=2S_k(n-1)+n^k\&&&&=dots\&&&&=sum_{i=1}^n2^{n-i}i^k\&&&&=2^nsum_{i=1}^n2^{-i}i^k\&&&a&overset{mathrm {def}}{=}2^{-1}\Rightarrow &&S_k(n)&&=2^nsum_{i=1}^na^ii^kend{split}
]
自然我们可以得到：
[f_k(n)=n^k+2^{n-1}sum_{i=1}^{n-1}a^ii^k
]
为了方便，我们就可以定义出：
[F_k(n)=sum_{i=1}^{n-1}a^ii^k
]
问题于是转化为了：如何快速求出(F_k(n))？
60 pts ~ 80 pts
可以发现(F_k(n))就是一个前缀和，我们尝试错位相减：
[egin{split}(1-a)F_k(n)&=sum_{i=1}^na^ii^k-sum_{i=1}^{n+1}a^i(i-1)^k\&=sum_{i=1}^na^i[i^k-(i-1)^k]-a^{n+1}n^k\&=sum_{i=1}^na^isum_{j=1}^kinom k j(-1)^{j+1}i^{k-j}-a^{n+1}n^k\&=sum_{j=1}^kinom k j(-1)^{j+1}sum_{i=1}^na^ii^{k-j}-a^{n+1}n^k\&=sum_{j=1}^kinom k j(-1)^{j+1}F_{k-j}(n)-a^{n+1}n^kend{split}
]
并且当(k=0)时，(F_0(n)=sum_{i=1}^na^i=frac{1-a^{n}}{1-a})，可以快速计算。通过这个递推关系我们可以(O(k^2))计算出(F_k(n))。用一个任意模数 NTT 说不定可以优化到(O(klog_2k))。
说不定就是说我也没试过
100 pts
考虑一个(f_k(n))的不存在求和符的递推式：
[egin{split}f_k(n)-f_k(n-1)&=n^k+F_k(n-1)\f_k(n)&=2f_k(n-1)+n^k-(n-1)^kend{split}
]
尝试构造出一个(g_k(n))满足：
[f_k(n)+g_k(n)=2(f_k(n-1)+g_k(n-1))
]
显然一个可行的解为：
[g_k(n)=2g_k(n-1)-n^k+(n-1)^k
]
即可以解出：
[f_k(n)=2^{n-1}(1+g_k(1))-g_k(n)
]
这也就意味着：
[F_k(n)=-a^{n-1}(g_k(n)+n^k)+(g_k(1)+1)
]
定义(G_k(n)=-2(g_k(n)+n^k))，我们就有：
[F_k(n)=a^nG_k(n)-aG_k(1)
]
又可以通过(g_k(n))的通项公式推导得出：
[G_k(1)=2G_k(0)
]
于是就可以得到
[F_k(n)=a^nG_k(n)-G_k(0)
]
根据某些奇奇怪怪的东西我们可以知道(G_k(n))是关于(n)的不超过(k)次的多项式。
i.e.作者自己也不会证明
这就意味着，我们只需要知道了(G_k(n))的(k+1)个点值，就可以为所欲为使用拉格朗日插值法计算出(G_k)在任何一个位置的点值。
另外，我们不难写出如下式子：
[G_k(n)=2^nG_k(0)+2^nF_k(n)
]
这意味着，我们只需要计算出(F_k(n))，我们就可以得到(G_k(n))。鉴于 (F_k(n))有其简洁的递推形式，我们可以快速而方便地得到(G_k(1)dots G_k(k+1))与(G_k(0))的关系。
考虑到如下的高阶差分公式：
[Delta^kf(k)=sum_{i=0}^kinom{k}{i}(-1)^{k-i}f(k+i)
]
推导难度不大，自行定义算子(mathrm Rf(n)=f(n+1))并重定义差分，展开式子即可。
不过有一个很有用的性质——任何(k)次的多项式在经过(k+1)次的差分后，结果必然是(0)。
这就是很像是求导了（差分本身就可以理解为离散中的导）。于是应有：
[Delta^{k+1}G_k(0)=sum_{i=0}^{k+1}inom{k+1}{i}(-1)^{k+1-i}G_k(i)=0
]
根据这个式子我们可以列出(G_k(1)dots G_k(k+1))关于(G_k(0))的关系。而(G_k(1)dots G_k(k+1))本身可以使用(G_k(0))表示，因而我们就相当于得到了一个关于(G_k(0))的一次方程，直接解出(G_k(0))即可。
另外，还有一种列方程的方法。我们直接利用拉格朗日插值，可以得出：
[G_k(0)=sum_{i=1}^{k+1}G_k(i)prod_{j
ot=i}frac{-j}{(i-j)}
]
这同样是(G_k(1)dots G_k(k+1))和(G_k(0))的关系，解方程就好。
进而我们可以得到(G_k(1)dots G_k(k+1))，然后就可以插值得到(G_k(n))，最后计算出(F_k(n))和(f_k(n))。
时间复杂度(O(k)sim O(klog_2m))，其中(m)为模数。
感觉......本题的核心在其数学构造。 60 pts ~ 80 pts 中出现的错位相减是一个常见的推前缀和通项的方法。但是考场上好像把减的对象搞错了。 100 pts 中的构造很巧妙。对于(G_k(n))的多项式性质似乎有点 “ 猜想 ” 的意味，直接证明好像有难度，如果可以证明可以评论，谢谢。
代码
#include <cstdio>
const int mod = 1e9 + 7, phi = mod - 1;
const int MAXK = 1e6 + 5;
template<typename _T>
void read( _T &x )
{
x = 0;char s = getchar();int f = 1;
while( s > '9' || s < '0' ){if( s == '-' ) f = -1; s = getchar();}
while( s >= '0' && s <= '9' ){x = ( x << 3 ) + ( x << 1 ) + ( s - '0' ), s = getchar();}
x *= f;
}
template<typename _T>
void write( _T x )
{
if( x < 0 ){ putchar( '-' ); x = ( ~ x ) + 1; }
if( 9 < x ){ write( x / 10 ); }
putchar( x % 10 + '0' );
}
int inv( const int );
int qkpow( int, int );
struct func
{
int a, b;
func() { a = b = 0; }
func( const int C ) { a = 0, b = C; }
func( const int A, const int B ) { a = A, b = B; }
func operator + ( const func &g ) const { return func( ( a + g.a ) % mod, ( b + g.b ) % mod ); }
func operator - ( const func &g ) const { return func( ( a - g.a + mod ) % mod, ( b - g.b + mod ) % mod ); }
func operator * ( const func &g ) const { return func( ( 1ll * a * g.b % mod + 1ll * b * g.a % mod ) % mod, 1ll * b * g.b % mod ); }
func operator / ( const int &g ) const { int t = inv( g ); return func( 1ll * a * t % mod, 1ll * b * t % mod ); }
void operator += ( const func &g ) { *this = *this + g; }
void operator -= ( const func &g ) { *this = *this - g; }
void operator *= ( const func &g ) { *this = *this * g; }
void operator /= ( const int &g ) { *this = *this / g; }
int getX( const int y ) const { return 1ll * ( y - b + mod ) * inv( a ) % mod; }
int getY( const int x ) const { return ( 1ll * a * x % mod + b ) % mod; }
};
func coe[MAXK];
int G[MAXK], neg[MAXK];
int fac[MAXK], finv[MAXK];
int prime[MAXK], pn, pwk[MAXK];
char S[MAXK];
int N, K;
bool isPrime[MAXK];
int qkpow( int base, int indx )
{
int ret = 1;
while( indx )
{
if( indx & 1 ) ret = 1ll * ret * base % mod;
base = 1ll * base * base % mod, indx >>= 1;
}
return ret;
}
int C( const int n, const int m ) { return 1ll * fac[n] * finv[m] % mod * finv[n - m] % mod; }
int inv( const int n ) { return n <= K + 10 ? 1ll * finv[n] * fac[n - 1] % mod : qkpow( n, mod - 2 ); }
void EulerSieve( const int siz )
{
pwk[1] = 1;
for( int i = 2 ; i <= siz ; i ++ )
{
if( ! isPrime[i] ) prime[++ pn] = i, pwk[i] = qkpow( i, K );
for( int j = 1 ; j <= pn && 1ll * i * prime[j] <= siz ; j ++ )
{
isPrime[i * prime[j]] = true, pwk[i * prime[j]] = 1ll * pwk[i] * pwk[prime[j]] % mod;
if( ! ( i % prime[j] ) ) break;
}
}
}
void init( const int siz )
{
fac[0] = 1;
EulerSieve( siz );
for( int i = 1 ; i <= siz ; i ++ ) fac[i] = 1ll * fac[i - 1] * i % mod;
finv[siz] = qkpow( fac[siz], mod - 2 );
for( int i = siz - 1 ; ~ i ; i -- ) finv[i] = 1ll * finv[i + 1] * ( i + 1 ) % mod;
}
int Lagrange( const int n )
{
#define M K + 1
if( n <= M ) return G[n];
neg[0] = 1; int up = 1, ans = 0;
for( int i = 2 ; i <= M ; i ++ ) up = 1ll * up * ( n - i ) % mod;
for( int i = 1 ; i <= M ; i ++ ) neg[i] = 1ll * neg[i - 1] * inv( mod - i ) % mod;
for( int i = 1 ; i <= M ; i ++ )
{
ans = ( ans + 1ll * up * finv[i - 1] % mod * neg[M - i] % mod * G[i] % mod ) % mod;
if( i < M ) up = 1ll * up * ( n - i ) % mod * inv( n - i - 1 ) % mod;
}
return ans;
#undef M
}
int main()
{
int ind = 0;
scanf( "%s %d", S, &K );
for( int i = 0 ; S[i] ; i ++ )
N = ( 10ll * N + S[i] - '0' ) % mod,
ind = ( 10ll * ind + S[i] - '0' ) % phi;
init( K + 10 ), coe[0] = func( 1, 0 );
int F = 0, alph = inv( 2 ), pwa = 1, pw2 = 2;
for( int i = 1 ; i <= K + 1 ; i ++ )
{
F = ( F + 1ll * pwa * pwk[i - 1] % mod ) % mod;
coe[i] = func( pw2, 1ll * pw2 * F % mod );
pwa = 1ll * pwa * alph % mod, pw2 = 2ll * pw2 % mod;
}
func tmp;
for( int i = 0 ; i <= K + 1 ; i ++ )
{
if( ( K + 1 - i ) & 1 ) tmp -= coe[i] * C( K + 1, i );
else tmp += coe[i] * C( K + 1, i );
}
G[0] = tmp.getX( 0 );
for( int i = 1 ; i <= K + 1 ; i ++ )
G[i] = coe[i].getY( G[0] );
int GN = Lagrange( N );
int FN = ( 1ll * qkpow( 2, phi - ind ) * GN % mod - G[0] + mod ) % mod;
write( ( qkpow( N, K ) + 1ll * FN * qkpow( 2, ind - 1 + phi ) % mod ) % mod ), putchar( '
' );
return 0;
}