In the great answers given so far, I didn't find a simple direct description of the natural comultiplication operation on symmetric functions. So I'll just copy and paste my comment from this question which tries to explain it in as few characters as are allowed in a comment. It shows why infinitely many variables are needed: we need the Hilbert hotel (twice the number of variables must equal the number of variables). OK, since I see copy-and-paste from comments doesn't work without manual repair anyway, I'll add a few words of clarification as well.

Symmetric functions are (degree-bounded power series) in infinitely many variables, and order doesn't matter (due to the "symmetric" part). Now rename the variables $x_0,y_0,x_1,y_1,x_2,\ldots$ and decompose a symmetric function $s=\sum_i u_iv_i$ as a sum of products of a symmetric function $u_i$ in the $x$'s and one $v_i$ in the $y$'s; then $\Delta(s)=\sum_i u_i\otimes v_i$. For instance for elementary symmetric function one has $\Delta(e_k)=\sum_{i+j=k}e_i\otimes e_j$ since the monomials can be arbitrarily spread across the $x$'s and $y$'s (and similarly for complete homogeneous symmetric functions), while for power sums $\Delta(p_k)=p_k\otimes 1+1\otimes p_k$, since the monomials involve a single $x_i$ or a single $y_i$, but the two cannot mix in a power sum.

I may add that multiplication and comultiplication are dual in operations in this case; Zelevinky call is calls it a (positive) self-adjoint Hopf-algebra.

1

In the great answers given so far, I didn't find a simple direct description of the natural comultiplication operation on symmetric functions. So I'll just copy and paste my comment from this question which tries to explain it in as few characters as are allowed in a comment. It shows why infinitely many variables are needed: we need the Hilbert hotel (twice the number of variables must equal the number of variables). OK, since I see copy-and-paste from comments doesn't work without manual repair anyway, I'll add a few words of clarification as well.

Symmetric functions are (degree-bounded power series) in infinitely many variables, and order doesn't matter. Now rename the variables $x_0,y_0,x_1,y_1,x_2,\ldots$ and decompose a symmetric function $s=\sum_i u_iv_i$ as a sum of products of a symmetric function $u_i$ in the $x$'s and one $v_i$ in the $y$'s; then $\Delta(s)=\sum_i u_i\otimes v_i$. For instance for elementary symmetric function one has $\Delta(e_k)=\sum_{i+j=k}e_i\otimes e_j$ since the monomials can be arbitrarily spread across the $x$'s and $y$'s (and similarly for complete homogeneous symmetric functions), while for power sums $\Delta(p_k)=p_k\otimes 1+1\otimes p_k$, since the monomials involve a single $x_i$ or a single $y_i$, but the two cannot mix in a power sum.

I may add that multiplication and comultiplication are dual in operations in this case; Zelevinky call is a (positive) self-adjoint Hopf-algebra.