convolution of all arrangements of a multiset You are given the set of numbers 2,7,8,9,9,23,33,41,41,57 and want to find the sum of all convolutions of the possible arrangement of these ten terms.  Each arrangement is convolved with itself.  One can find the answer by examining case by case each arrangement, but this becomes impracticable if the number of terms is larger than ten.  Is there any faster method to find the sum of all these convolutions?
 A: Suppose the multiset $M$ has unique terms $x_1,x_2,x_3,\ldots,x_k$ with these terms appearing $n_1,n_2,n_3,\ldots,n_k$ times each.
Then the total number of unique elements is $\displaystyle T(M)=\sum_{i=1}^k n_i$ and the total number of arrangements is $\displaystyle A(M)=\frac{T(M)!}{\prod_i n_i!}$ with the caveat that $A(M)=0$ if any of the $n_i$ are negative.  
Also consider the multisets: $M_{-i}$ where $n_i$ is reduced by $1$ and $M_{-i-j}$ where $n_i$ and $n_j$ are each reduced by $1$ ($i$ and $j$ can be the same or different: if the same then this means $n_i$ is reduced by $2$). 
Then I suspect that your result is something like 
$$\displaystyle 2\big\lfloor \tfrac{T(M)}{2} \big\rfloor\left( \sum_{i=1}^k \sum_{j=1}^k x_i x_j A(M_{-i-j})    \right) + \left(T(M)-2\big\lfloor\tfrac{T(M)}{2}\big\rfloor \right)\left( \sum_{i=1}^k x_i^2 A(M_{-i})   \right) $$ and you can factor out $A(M)$ using  $A(M_{-i}) = A(M)\frac{n_i}{T(M)}$, $A(M_{-i-i}) = A(M)\frac{n_i(n_i-1)}{T(M)(T(M)-1)}$ and $A(M_{-i-j}) = A(M)\frac{n_i n_j}{T(M)(T(M)-1)}$ when $i\not = j$. [Added] The $2\big\lfloor \tfrac{T(M)}{2} \big\rfloor$ term is $T(M)$ when this is even and $T(M)-1$ when $T(M)$ is odd, while the  $\left(T(M)-2\big\lfloor\tfrac{T(M)}{2}\big\rfloor \right)$ term is $0$ when $T(M)$ is even and $1$ when $T(M)$ is odd. 
If you do that and then let  $\displaystyle S(M)=\sum_{i=1}^k n_i x_i$, i.e. the sum over the multiset, and let $\displaystyle Q(M)=\sum_{i=1}^k n_i x_i^2$, i.e. the sum of squares over the multiset, then I think you get the fairly simple
$$\dfrac{A(M)\left(S(M)^2 - Q(M)  \right)  }{(T(M)-1)} \text{ when } T(M) \text{ is even}$$  
$$\dfrac{A(M)\,S(M)^2 }{T(M)} \text{ when } T(M) \text{ is odd}$$  
For the toy example of $M=\{8,9,9\}$, this has $x_1=8,x_2=9,n_1=1,n_2=2$ and so $T(M)=1+2=3$ and $A(M)=\frac{3!}{1! \times 2!}=3$, $S(M)=8+9+9=26$, $Q(M) = 8^2+9^2+9^2 = 226$,  so the result becomes $\frac{3\times 26^2}{3} = 676$ as expected 
A: [Note: I think my solution is slightly different from Henry's and I had much of it written already, so I'll post anyway.]
Assume first that the multiset does not contain any repeated elements, so that its elements $\{a_1, \ldots, a_n\}$ are distinct.
Suppose first $n$ is even. Let $a_i$ and $a_j$ be any pair of elements in the multiset, and count the number of arrangements in which the product $a_ia_j$ is used in the convolution. There are $n$ locations for $a_i$, and then only one choice for $a_j$; the remaining elements can be distributed arbitrarily in $(n-2)!$ ways. For each of these $n(n-2)!$ arrangements, this pair contributes a summand of $2a_ia_j$ to the convolution; hence to the sum of all convolutions, this pair contributes a total of $2a_ia_jn(n-2)!$. Thus the sum of all convolutions is:
$$\displaystyle 2n(n-2)!\sum_{1 \le i < j \le n} a_ia_j$$
If $n$ is odd, we can again let $a_i$ and $a_j$ be any pair of elements in the multiset, and count the number of arrangements where the product $a_ia_j$ is used in the convolution. Since $n$ is odd, there are only $n-1$ locations for $a_i$, since it cannot go in the middle position. But the rest of the calculation is the same, meaning this pair contributes a total of $2a_ia_j(n-1)!$ to the sum of all convolutions. Unlike the even case, a single element $a_i$ can contribute to the convolution without a partner by being in the middle position, where it contributes $a_i^2$ to the convolution. This term appears in $(n-1)!$ arrangements, adding a total of $a_i^2(n-1)!$ to the sum of all convolutions. Thus the total convolution sum in this case is:
$$\displaystyle 2(n-1)!\sum_{1 \le i < j \le n} a_ia_j + (n-1)!\sum_{i=1}^n a_i^2$$
Now suppose our multiset consists of distinct elements $b_1 \ldots b_m$, where $b_i$ has multiplicity $n_i$, and let $n = \sum n_i$. For the moment, treat the $n$ elements of the multiset as being distinct objects, though some might have the same numerical value. Then the analysis above would still give the total convolution sum over all $n!$ arrangements.
By wishing to treat objects of the same numerical value as being equivalent, we are essentially creating equivalence classes over the set of $n!$ arrangements of our distinct objects such that two arrangements are equivalent if both arrangements have the same ordering of numerical values. Clearly, each equivalence class has $\prod_{i=1}^m n_i!$ elements, and moreover, the convolution of each element in the equivalence class is the same. Since each contributing summand to the total convolution sum is counted exactly this many times, we can get the desired answer simply by dividing the above quantities by $\prod_{i=1}^m n_i!$, to get the formulas:
For $n$ even:
$$\frac{2n(n-2)!}{\prod_{i=1}^m n_i!} \left( \displaystyle \sum_{1 \le i < j \le m}b_in_ib_jn_j + \sum_{i=1}^m\frac{n_i(n_i-1)}{2}b_i^2 \right)$$
For $n$ odd:
$$\frac{(n-1)!}{\prod_{i=1}^m n_i!} \left( \displaystyle \sum_{1 \le i < j \le m}2b_in_ib_jn_j + \sum_{i=1}^m n_i^2b_i^2 \right)$$
A: Suppose for the moment that the multiset $S$ contains $n$ unique integers and let $x$ be a vector of these integers.
Let $R$ denote the $n\times n$ reverse order matrix and let $P$ denote an $n\times n$ permutation matrix corresponding to
any of the $n!$ permutations $\sigma$ of $\{1,\ldots,n\}$. Then the convolution of any arrangement of $S$ may be expressed
as $x^TP^TRPx$. Let $Q$ denote the sum of all such matrices of the form $P^TRP$. Then it follows that
the sum of all possible convolutions is given by $x^TQx$. All that remains is to determine the elements of $Q$. 
Consider a matrix of the form $P^TRP$ corresponding to the permutation $\sigma$. By the nature of $P$ and $R$ it follows that
$$
(P^TRP)_{ij} = p_i^TRp_j = 
\left\{
\begin{array}{ll}
1 & \text{if}~\sigma^{-1}(i) = n - \sigma^{-1}(j) + 1\\
0 & \text{otherwise}
\end{array}\right.
$$
for all $i,j = 1,\ldots,n$.
Hence, by a simple counting argument
$$
q_{ij} = \left\{
\begin{array}{ll}
0 & n~\text{even},~i = j\\
n(n-2)! & n~\text{even},~i \neq j\\
(n-1)! & n~\text{odd}\\
\end{array}\right.
$$
Now consider the case in which $S$ consists of $k$ distinct elements repeated $n_1,\ldots,n_k$ times.
Let $n = n_1 + \cdots + n_k$. Then there are
$$
\frac{n!}{n_1!n_2!\cdots n_k!}
$$
unique permutations of $S$ and so we have overcounted by a factor of $N = (n_1!n_2!\cdots n_k!)$. Thus, the general solution is
$$
\frac{1}{N}x^TQx.
$$
Lastly we note that $Q$ is symmetric, which further reduces the work to compute the sum. Clearly, a computational routine would not form $Q$ explicitly but would rather hard code the product $x^TQx$ to reduce the memory overhead.
