I am trying to find the poles an residues of an expression of the form
$$\left(\sum_i \zeta(m_i) \right)^2$$
where $m_i\in \mathcal{M}$ some set of integers including 1. Now I know that on the we only have a single simple pole at 1, so I would imagine that the square I have above has both simple poles arising from
$$\zeta(1)\sum_{m_i \neq1} \zeta(m_i) $$
with residue $\sum_{m_i \neq1} \zeta(m_i)$ as well as a second order pole from $\zeta(1)\zeta(1)$.
On the other hand
$$\left(\sum_i \zeta(m_i) \right)^2 = \left(\sum_i\sum_n n^{-m_i} \right)^2 = \sum_{i,j}\sum_n n^{-m_i-m_j} = \sum_{i,j}\zeta(m_i+m_j) = \sum_i p^2_{\mathcal{M}} (m_i)\zeta(m_i)$$
where $p^2_{\mathcal{M}} (m_i)$ is the number of solutions to $m+n = m_i$ with $m,n \in \mathcal{M}$. In any case this only has a simple pole at $m_i=1$ with residue $p^2_{\mathcal{M}} (1)$. .
I dont see why these two should be equivalent (in fact, barring some mathemagic, I cant see how these are), but I also dont see which of these is correct.
The error in the calculation is in the second equality. You correctly turned the index $i$ inside the square into two indices $i,j$ when expanding the square out, but you forgot to do the same with the single index $n$—it should become two indices $k,n$ after squaring out. (Note that if this calculation were correct, it would show in particular the incorrect identity $\zeta(s)^2=\zeta(2s)$.)