Sum of sum of $k$th power of first $n$ natural numbers. I was working on a problem which involves computation of $k$-th power of first $n$ natural numbers.
Say $f(n) = 1^k+2^k+3^k+\cdots+n^k$ we can compute $f(n)$ by using Faulhaber's Triangle also by fast computation of Bernoulli numbers.
Here we have to compute $G(n)= f(1)+f(2)+f(3)+\cdots+f(n)$
I did some work on this to simplify but unable to come up with any easy equation.
Can someone help me how to compute $G(n)$ efficiently , I can efficiently compute $f(n)$.
Here $n<123456789$
and $k<321$.
Thanks.
 A: $$G_k(n)= \sum_{i=1}^n f_k(i)= \sum_{i=1}^n \sum_{j=1}^i j^k = \sum_{j=1}^n (n+1-j) j^k=(n+1)\sum_{j=1}^n  j^k - \left(\sum_{j=1}^n j^{k+1} \right)\\
=(n+1)f_k(n)-f_{k+1}(n)$$
Since you said you can efficiently calculate 
$$f_{k}(n) =\sum_{j=1}^n j^k$$
you can also calculate efficiently $$(n+1)f_k(n)-f_{k+1}(n)=g_k(n)$$
A: You might find it useful to visualize the sum:
$$
\begin{matrix}
1^k \\
1^k & 2^k \\
1^k & 2^k & 3^k \\
1^k & 2^k & 3^k & 4^k \\
\vdots&\vdots&\vdots&\vdots&\ddots \\
1^k & 2^k & 3^k & 4^k & \dots & n^k \\
\end{matrix}
$$
$$
G_k(n) = n\cdot1^k + (n-1)\cdot2^k + (n-2)\cdot3^k + \dots + 2\cdot(n-1)^k + 1\cdot n^k
     = \sum_{s=1}^{n} (n-s+1)s^k
$$
This will take $O(n\lg k)$ multiplications to compute directly. I'm trying to find a better way.
EDIT: I was on the right track. For glorious ending see answer posted by N. S.
A: You can compute the sum of first $n$ $k$-th powers in $O(k^2)$
We define $$F_k(n) = \sum_{i=1}^{n} i^k$$
Forse example, $F_0(n) = n$ and $F_1(n) = \frac{n(n+1)}{2} $.
Now, we can prove that
$$ F_k(n) = (n+1)^{k+1} -1 - \sum_{i=2}^{k+1} \binom{k+1}{i} \cdot F_{k-1+i}(n)$$
You simply compute the first $k$ rows of Pascal Triangle, and then for every $j \le k$ you compute $F_j(n)$
