How to prove $\sum\limits_{i=1}^n \sum\limits_{j=1}^n \min\{s_i,s_j\}=\sum\limits_{k=1}^M c_k^2$ Let us assume that $s_1,\dots,s_n$ are some given non-negative integers. Let $M=\max\{s_1,\dots,s_n\}$ and for any $k$ let use denote $c_k=|\{i=1,\dots,n; s_i\ge k\}|$. (In the other words, $c_k$ is the number of elements of the given set which have size at least $k$.) Then
$$\sum\limits_{i=1}^n \sum\limits_{j=1}^n \min\{s_i,s_j\}=\sum_{k=1}^M c_k^2.\tag{*}$$
(Notice that $c_0=n$, however the term $c_0^2$ is not included on the RHS.)
How can we prove this?
This equality seems to be rather elementary (and I have included a proof by induction below). Still, I suppose that there are many interesting ways to prove this, so I thought that it might be reasonable to ask here to see various approaches to this sum and various ways to formulate the arguments why this equality holds.
Motivation. This expression appears naturally in some considerations about commuting matrices. See, for example, Robert Israel's answer to Dimension of a subspace of $M_n(\mathbb C)$. (The sum here is the formula he mentions in the last paragraph, in the case when there is only a single eigenvalue.)
Example. Just to make sure that we are not trying to prove something which does not hold, let's have a look at least at one example. (After all, having some examples in front of us might be also useful when thinking about proofs.)
Let us try the numbers $1$, $1$, $2$, $3$, $3$, i.e., $s_1=s_2=1$, $s_3=2$, $s_4=s_5=3$.
We can write down the following table, where we have $\min\{s_i,s_j\}$ in the position $(i,j)$ and the last row/column contain the sums.
$$
\begin{array}{ccccc|c}
  1 & 1 & 1 & 1 & 1 & 5 \\
  1 & 1 & 1 & 1 & 1 & 5 \\
  1 & 1 & 2 & 2 & 2 & 8 \\
  1 & 1 & 2 & 3 & 3 &10 \\
  1 & 1 & 2 & 3 & 3 &10 \\\hline
  5 & 5 & 8 &10 &10 &38
\end{array}
$$
At the same time we get $c_1^2+c_2^2+c_3^2=5^2+3^2+2^2=25+9+4=38$.
A special case. In particular, if we take $s_i=i$, then we get the sum
$$\sum\limits_{i=1}^n \sum\limits_{j=1}^n \min\{i,j\}=\sum_{k=1}^n (n-k+1)^2 = \sum_{k=1}^n k^2 = \frac{n(n+1)(2n+1)}6.\tag{⋄}$$
This sum seems to me as a rather natural exercise to give students for manipulation with sums and for trying out combinatorial proofs.
I have tried searching whether there is some question on this site about $(\diamond)$. I did not find such question - the closest one I've seen was this one: How prove this identity $\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{n}\min{\{i,j\}}(a_{i}-a_{i+1})(b_{j}-b_{j+1})=\sum\limits_{i=1}^{n}a_{i}b_{i}$?
EDIT: After searching a bit more I found Combinatorial proof of $\sum_{1\le i\le n,\ 1\le j\le n}\min(i,j)=\sum_{i=1}^ni^2=\frac{n(n+1)(2n+1)}6$ and A Combinatorial proof for the identity $\sum_i \sum_j \min(i,j) = \sum_k k^2$
Proof by induction. There are certainly various more clever proofs, but as a starting point we can try induction on $M$.
$1^\circ$ If $M=0$, then both sides of $(*)$ are equal to zero.
$2^\circ$  Let us assume that the claim is true whenever we have number with maximum less than $M$.
W.l.o.g. we can assume that $s_i\ge1$ for each $i$. (If suffices to notice that by including also $s_i$'s for which $s_i=0$ we add zero to both sides. The RHS is not influenced at all, since such $s_i$'s do not count in $c_k$ for $k\ge1$. On the LHS we are only adding/omiting the terms of the form $\min\{0,s_i\}=0$.)
With this assumption we can use $s'_i=s_i-1$. We can then apply the induction hypothesis on the numbers $s'_i$. By noticing that $\min\{s_i,s_j\}=\min\{s'_i,s'_j\}+1$ and $c_k=c'_{k=1}$ we get
\begin{align*}
S
&=\sum\limits_{i=1}^n \sum\limits_{j=1}^n \min\{s_i,s_j\} \\
&=\sum\limits_{i=1}^n \sum\limits_{j=1}^n (1+\min\{s'_i,s'_j\}) \\
&=n^2+\sum\limits_{i=1}^n \sum\limits_{j=1}^n \min\{s'_i,s'_j\} \\
&=n^2+(c'_1)^2+\dots+(c'_{M-1})^2= \\
&=m^2+c_2^2+\dots+c_M^2= \\
&=c_1^2+c_2^2+\dots+c_M^2
\end{align*}
 A: Let $X$ be a random variable such that $\Pr(X=s_i)=\frac 1n k_i$, where $k_i$ is the number of times where $s_i$ appears in the $n$-uple $\left(s_1,\dots,s_n\right)$. Let $Y$ be an independent copy of $X$. Then 
$$
\frac 1{n^2}\sum\limits_{i=1}^n \sum\limits_{j=1}^n \min\{s_i,s_j\} = \mathbb E\left[\min\{X,Y\}\right].$$
One the other hand, the latter expectation can be computed by the formula 
$$
\mathbb E\left[Z\right]=\int_0^{+\infty}\Pr(Z\geqslant t)dt, \quad Z\geqslant 0.
$$
With $Z:=\min\{X,Y\}$, we get by independence and the fact that $X$ and $Y$ have the same distribution that 
$$
\Pr\left(\min\{X,Y\}\geqslant t\right)=\Pr\left(X\geqslant t,Y\geqslant t\right)=\Pr\left(X\geqslant t\right)\Pr\left(Y\geqslant t\right)=\Pr\left(X\geqslant t\right)^2.
$$
Then cutting the integral at integers gives the wanted formula.
A: One way to show the equality is to show that both sides count the number of the triples
$$\{(i,j,k)\in\mathbb N^3; 1\le i \le n, 1\le j\le n, k\le s_i, 1\le k\le s_j\}.$$
Basically, the difference between the LHS is just the order of sums. We have
\begin{align*}
\sum_{i=1}^n \sum_{j=1}^n \sum_{k\le s_i,s_j} 1
&= \sum_{k=1}^M \sum_{\substack{1\le i\le n \\ s_i\ge k}}
\sum_{\substack{1\le j\le n \\ s_j\ge k}} 1 \\
&= \sum_{k=1}^M \sum_{\substack{1\le i\le n \\ s_i\ge k}} c_k \\
&= \sum_{k=1}^M c_k^2
\end{align*}
Exchanging the order of sums is similar to Martin R's answer.
Looking at both sides as counting some triples is perhaps similar to the combinatorial proofs (or picture proofs) posted in the linked questions A combinatorial proof for the identity $\sum_i \sum_j \min(i,j) = \sum_k k^2$ and Combinatorial proof of $\sum_{1\le i\le n,\ 1\le j\le n}\min(i,j)=\sum_{i=1}^ni^2=\frac{n(n+1)(2n+1)}6$.
A: Using Iverson brackets, i.e. $[\cdots]$, define
$$
c_k=\sum_{i=1}^n[s_i\ge k]\tag1
$$
then
$$
\begin{align}
\sum_{k=1}^Mc_k
&=\sum_{k=1}^M\sum_{i=1}^n[s_i\ge k]\tag{2a}\\
&=\sum_{i=1}^n\sum_{k=1}^M[s_i\ge k]\tag{2b}\\
&=\sum_{i=1}^ns_i\tag2
\end{align}
$$
Explanation:
$\text{(2a)}$: apply $(1)$
$\text{(2b)}$: change order of summation
$\phantom{a}(2)$: there are $s_i$ values of $k$ for which $s_i\ge k$
Similarly, define
$$
\begin{align}
b_k
&=\sum_{i=1}^n\sum_{j=1}^n[\min(s_i,s_j)\ge k]\tag{3a}\\
&=\sum_{i=1}^n\sum_{j=1}^n[s_i\ge k][s_j\ge k]\tag{3b}\\
&=\sum_{i=1}^n[s_i\ge k]\sum_{j=1}^n[s_j\ge k]\tag{3c}\\[6pt]
&=c_k^2\tag3
\end{align}
$$
Explanation:
$\text{(3a)}$: definition
$\text{(3b)}$: $\min(s_i,s_j)\ge k\iff s_i\ge k\text{ and }s_j\ge k$
$\text{(3c)}$: rearrange the product
$\phantom{a}(3)$: apply $(1)$
Then
$$
\begin{align}
\sum_{k=1}^Mb_k
&=\sum_{k=1}^M\sum_{i=1}^n\sum_{j=1}^n[\min(s_i,s_j)\ge k]\tag{4a}\\
&=\sum_{i=1}^n\sum_{j=1}^n\sum_{k=1}^M[\min(s_i,s_j)\ge k]\tag{4b}\\
&=\sum_{i=1}^n\sum_{j=1}^n\min(s_i,s_j)\tag4
\end{align}
$$
Explanation:
$\text{(4a)}$: apply definition $\text{(3a)}$
$\text{(4b)}$: change order of summation
$\phantom{a}(4)$: there are $\min(s_i,s_j)$ values of $k$ for which $\min(s_i,s_j)\ge k$
Equations $(3)$ and $(4)$ say
$$
\sum_{k=1}^Mc_k^2=\sum_{i=1}^n\sum_{j=1}^n\min(s_i,s_j)\tag5
$$
A: We can compute the sum of non-negative integers $a_1, \ldots,a_n$ in terms of their “distribution”:
$$
 \sum_{i=1}^n a_i = \sum_{i=1}^n \sum_{k=1}^\infty \Bbb 1_{k \le a_i}
= \sum_{k=1}^\infty \sum_{i=1}^n\Bbb 1_{k \le a_i} =\sum_{k=1}^\infty | \{ i=1,\ldots, n : a_i \ge k \} |
$$
where only finitely many terms are non-zero. In your case
$$
c_k=|\{i=1,\dots,n :  s_i\ge k\}|
$$
implies that
$$
c_k^2=|\{i,j=1,\dots,n :  \min\{s_i,s_j\}\ge k\}|
$$
and therefore
$$
\sum_{i,j=1}^n  \min\{s_i,s_j\}=\sum_{k=1}^\infty c_k^2 \, .
$$
