自然数幂和

首先我来介绍一下什么是自然数幂和：

$1+2+3+...+i+...+n=？$

$1^2 +2^2+3^2+...+i^2+...+n^2=？$

$1^k +2^k+3^k+...+i^k+...+n^k=？$

类似上述式子的就是自然数幂和了，那么具体怎么求呢，这就是今天的重点了：

$1+2+3+...+i+...+n=？$ 这个式子大家肯定都会求，这就是一个等差数列求和，然后套公式：

首先求 $(n+1)^{k+1}-n^{k+1}------------------------(1)$

根据牛顿二项式展开定理：

$(n+1)^{k+1}=C(k+1,0)*n^{k+1}+C(k+1,1)*n^k+...+C(k+1,k+1)*n^0$

所以：

$(1)式=C(k+1,1)*n^k+...+C(k+1,k+1)*n^0$

$(n+1)^{k+1}-n^{k+1}=C(k+1,1)*n^k+...+1--------------(2)$

那么我们现在消去了最高次幂的项 $x^{k+1}$ 那么我们现在对 $(2)式从1-n求和$ 得到：

$(n+1)^{k+1}-n^{k+1}=C(k+1,1)*n^k+...+1$

$n^{k+1}-(n-1)^{k+1}=C(k+1,1)*(n-1)^k+...+1$

…

$2^{k+1}-1^{k+1}=C(k+1,1)*1^k+...+1$

左边加左边，右边加右边，得到：

$(n+1)^{k+1}-1=C(k+1,1)*\sum_{i=1}^{n}i^k+...+n$

根据上式我们得到幂指数是 $k$ 的：

$\sum_{i=1}^ni^k=\frac1{k+1}*[(n+1)^{k+1}-...-C(k+1,i)*\sum_{i=0}^ni^{k+1-i}-...-n-1]$

也可以写为：

$\sum_{i=1}^ni^k=\frac1{k+1}*((n+1)^{k+1}-(C(k+1,i)*\sum_{i=0}^ni^{k+1-i}+...+n+1))$

现在 我们将 $\sum_{i=1}^ni^k$ 记作 $S(n,k)$

那么我们现在得到一个递推式：

$S(n,k)=\frac{1}{k+1}*((n+1)^{k+1}-(C(k+1,2)*S(n,k-1)+C(k+1,i)*S(n,k+1-i)+C(k+1,k)*S(n,1)+n+1))$

当 $k == 1$ 的时候是递归出口：

$1^2+2^2+3^2+...+n^2$

我们先对 $(n+1)^3-n^3$ 求和得到：

$(n+1)^3-n^3=3*n^2+3*n+1$

$n^3-(n-1)^3=3*(n-1)^2+3*(n-1)+1$

$...$

$(1+1)^3-b^3=3*1^2+3*1+1$

求和得：

$(n+1)^3-1^3=3*(1^2+2^2+...+n^2)+3*(1+2+...+n)+n$

那么$(1^2+2^2+...+n^2)=\frac{(n+1)^3-(3*(1+2+...n)+n+1)}3$

又因为：

$1+2+...+n=\frac{n*(n+1)}2$

所以：
$(1^2+2^2+...+n^2)=\frac{(n+1)^3-(3*\frac{n*(n+1)}2+n+1)}3$ $=\frac{n*(2*n+1)*(n+1)}{6}$
模板：

/**
2016 - 08 - 06 上午
Author: ITAK

Motto:

今日的我要超越昨日的我，明日的我要胜过今日的我，
以创作出更好的代码为目标，不断地超越自己。
**/

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <cmath>
#include <vector>
#include <queue>
#include <algorithm>
#include <set>
using namespace std;
typedef long long LL;
typedef unsigned long long ULL;
const int INF = 1e9+5;
const int MAXN = 2e3+5;
const LL MOD = 1e9+7;
const double eps = 1e-7;
const double PI = acos(-1);
using namespace std;
LL c[MAXN][MAXN];
void Get_Fac()
{
    for(int i=0; i<MAXN; i++)
        c[i][0] = 1;
    for(int i=1; i<MAXN; i++)
        for(int j=1; j<=i; j++)
            c[i][j] = (c[i-1][j]+c[i-1][j-1])%MOD;
}
LL ans[MAXN];
LL Solve(LL n, LL k)///1^k+2^k+...+n^k
{
    Get_Fac();
    if(ans[k] != -1)
        return ans[k];
    if(k == 1)
        return (n+1)*n/2;
    LL tmp = 1;
    for(int i=0; i<=k; i++)
        tmp = tmp*(n+1);
    tmp = tmp - (n+1);
    LL sum = 0;
    for(int i=1; i<k; i++)
        sum += c[k+1][i+1]*Solve(n,k-i);
    ans[k] = (tmp-sum)/(k+1);
    return ans[k];
}
int main()
{
    LL n, k;
    while(cin>>n>>k)
    {
        memset(ans, -1, sizeof(ans));
        cout<<Solve(n,k)<<endl;
    }
    return 0;
}