Explanation of a method to compute $\sum_{k \le n} k^2$ I was searching for methods to compute $\sum_{k\le n} k^2$. I stumbled across this (which is an answer provided by Gareth Rees to this question).

"..Represent $k^2$ in terms of falling powers (easy by inspection in this case, but you can use Stirling subset numbers to convert): $$ k^2 = k^{\underline 2} + k^{\underline 1}$$
Sums of falling powers are easy, just like integration of ordinary powers, except for the treatment of limits: $$ \sum_{k=1}^n k^{\underline 2} + k^{\underline 1} = \bigg({1\over 3}k^{\underline 3} + {1\over 2}k^{\underline 2}\bigg)\ \bigg|^{n+1}_0$$
And then convert back into ordinary powers (by expansion, or using signed Stirling cycle numbers): $$ {1\over 3}((n+1)^3 - 3(n+1)^2 + 2(n+1)) + {1\over 2}((n+1)^2 - (n+1))..."$$

Could someone explain what is going on in this answer and provide some useful and explanatory links?
Thank you.
 A: Partial sums of sequences $a_n$ that can be expressed as polynomials in $n$ are easily found using discrete calculus.
We start with the discrete version of the Fundamental Theorem of Integral Calculus:
\begin{align*}
    \sum\limits_{n=1}^k a_n &= \sum\limits_{n=1}^k (\Delta b)_n \\
         &= (b_2 - b_1)+(b_3 - b_2)+ \ldots + (b_k - b_{k-1}) + (b_{k+1} - b_k) \\
         &= b_{k+1} - b_1
\end{align*}
where $(\Delta b)_n = b_{n+1} - b_n$ is the forward difference. Finding the partial sum has now been reduced to finding a sequence $b_n$ such that $(\Delta b)_n = a_n$.
We will find $b$, the antiderivative of $a$, using falling powers, which are defined by
$$
    n^{\underline{k}} = n(n-1)(n-2)\ldots (n-k+1)
$$
where $k$ is an integer and, by a second definition,
$n^{\underline{0}}=1$. For example
$$
    n^{\underline{2}} = n(n-1).
$$
We now need one more result: the discrete derivative of $n^{\underline{k}}$ is given by
\begin{align*}
\Delta n^{\underline{k}} &= (n+1)^{\underline{k}} - n^{\underline{k}} \\
             &= (n+1)n^{\underline{k-1}} - n^{\underline{k-1}}(n-k+1) \\
             &= kn^{\underline{k-1}}
\end{align*}
Let's now find the partial sum:
\begin{align}
\sum^{k}_{n=1}n(n+1) &= \sum^{k}_{n=1} (n+1)^{\underline{2}} \\
&= \sum^{k}_{n=1} \Delta \left[\frac{1}{3} (n+1)^{\underline{3}}\right] \\
&= \frac{(k+2)(k+1)k}{3} - \require{cancel}\cancelto{0}{\frac{(1+1)(1+0)(1-1)}{3}}
\end{align}
edit: I recommend watching the following lectures by Robert Ghrist: Sequences, Sequences_bonus, Differences, Discrete Calculus. You can also read the corresponding wiki sections here. And, if you are willing to read second-rate material (written by me), then read this which discusses a more difficult series.
A: Instead of falling powers I might use Binomial Coefficients:
$$ k^2 = 2\, \frac{k(k-1)}{2} + k = 2\binom{k}{2}  + \binom{k}{1}$$
Now if we add over all $k$ we can get other binomial coefficients:
$$ \sum_{k \leq n} \binom{k}{2} = \binom{n+1}{3} \text{ and } \sum_{k \leq n} \binom{k}{1} = \binom{n+1}{2}$$
This feels just like the forumla $\int x^k = \tfrac{1}{k+1} x^{k+1} + C$ and finally 
$$ \sum_{k \leq n} k^2 
= 2\sum_{k \leq n} \binom{n}{2}  + \sum_{k \leq n} \binom{n}{1}
= 2\binom{n+1}{3} + \binom{n+1}{2}
= 2\cdot  \frac{(n+1)n(n-1)}{6} + \frac{(n+1)n}{2}$$
This recovers the correct formula since:
$$ 2\cdot  \frac{(n+1)n(n-1)}{6} + \frac{(n+1)n}{2} = \frac{n(n+1)(2n+1)}{6}$$
A: The method uses the "falling powers" concept: instead of the usual powers like $k^3:=k\cdot k\cdot k$, consider the "partial factorials" $k^{\underline 3}=k\cdot(k-1)\cdot(k-2)$.
The falling powers fulfill an interesting summation formula that reminds of the antiderivative of the powers ($\int x^3\,dx=x^4/4$):
$$\sum_{j=1}^kj^{\underline3}=\frac{(k+1)^{\underline 4}}4.$$
As you can check,
$$3\cdot2\cdot1=\frac{4\cdot3\cdot2\cdot1}4\\
3\cdot2\cdot1+4\cdot3\cdot2=\frac{5\cdot4\cdot3\cdot2}4\\
3\cdot2\cdot1+4\cdot3\cdot2+5\cdot4\cdot3=\frac{6\cdot5\cdot4\cdot3}4\\
\cdots$$
The property holds for other integer exponents.
So the spirit of the method is to convert ordinary powers to a sum of falling powers, apply the summation formula and convert back to ordinary powers.
$$j^2=j\cdot(j-1)+j=j^{\underline 2}+j^{\underline 1},$$
yielding after summation
$$\frac{(k+1)^{\underline3}}3+\frac{(k+1)^{\underline2}}2=\frac{(k+1)\cdot k\cdot(k-1)}3+\frac{(k+1)\cdot k}2=(k+1)\cdot k\cdot\frac{2k+1}6.$$
A: Here's an alternative proof:
$$\sum_{r=1}^n r^2=\frac 14\sum_{r=1}^n(2r)^2=\frac 14\sum_{r=1}^{n}\binom{2r}2+\binom{2r+1}2=\frac14\sum_{r=1}^{2n+1}\binom r2=\frac 14\binom{2n+2}3$$
This can also be written with falling factorials using the fact that $k^2=\frac12\left[ k^\underline{2}+(k+1)^\underline{2}\right]$, as shown below:
$$\sum_{r=1}^n r^2
=\frac 14\sum_{r=1}^n\color{blue}{(2r)^2}
=\frac 14\cdot\color{blue}{\frac12}\sum_{r=1}^{n}\color{blue}{\bigg[(2r)^\underline{2}+(2r+1)^\underline{2}\bigg]}
=\frac18\sum_{r=1}^{2n+1} r^\underline{2}
=\frac 18\cdot \frac{(2n+2)^\underline{3}}{3}
=\frac{(2n+2)^\underline{3}}{24}$$
Both expressions above expand directly to
$$\color{lightgray}{\sum_{r=1}^nr^2=}\frac{(2n+2)(2n+1)(2n)}{24}=\frac16 n(n+2)(2n+1)$$
