In his book Shift Register Sequences, Solomon Golomb refers to "Gauss's product formula":
The basic device here is “Gauss's product formula” [39] which expresses the “product” of any two cosets as a sum of cosets. (This "product" is in terms of evaluation with 1's and 0's. It has nothing to do with multiplication in the factor group.)
I'm completely lost. What is he talking about? The cosets here in Golomb's example are cyclotomic cosets $\bmod 2^N - 1$ for $N=5$, with
$$\begin{align} C_0 &= \{1,2,4,8,16\} \cr C_1 &= \{3,6,12,24,17\} \cr C_2 &= \{9,18,5,10,20\} \cr C_3 &= \{27,23,15,30,29\} \cr C_4 &= \{19,7,14,28,25\} \cr C_5 &= \{26,21,11,22,13\} \end{align}$$
I understand what these cosets are (elements in the same coset are produced by multiplying by powers of 2, $\bmod 2^N-1$), but how can you possibly "multiply" or "add" them, and write equations regarding $C_iC_j$?

My guess is that Golomb wanted to avoid use of structures from abstract algebra. I would describe this product as follows. If you reach the end you may understand why Golomb wanted to use his simpler description! I do this anyway because this gives us many nice properties of the product (proving e.g. associativity may not be obvious otherwise).
All these arithmetic operations take place in the quotient ring of polynomials with coefficients in $\Bbb{F}_2$. When dealing with cyclotomic cosets modulo $m$ (an odd natural number) we use the ring $$ R_m:=\Bbb{F}_2[x]/(x^m-1). $$ You may have seen this ring used when studying cyclic codes of length $m$ also. Anyway, we view a cyclotomic coset modulo $m$ as an element of $R_m$. If $C$ is such a cyclotomic coset, we can think of it as the polynomial $C(x)$ defined simply as $$ C(x)=\sum_{i\in C}x^i. $$ So for example in $R_{31}$ we have $$ C_0(x)=x+x^2+x^4+x^8+x^{16} $$ and $$ C_1(x)=x^3+x^6+x^{12}+x^{24}+x^{17}. $$ When you multiply those polynomials in the ring $R_{31}$ you multiply them almost like usual polynomials. 1) You need to keep in mind that the arithmetic of the coefficients is done modulo two. 2) You also need to keep in mind that the arithmetic of exponents is done modulo $m$. So here for example $$ x^{16}\cdot x^{17}=x^{16+17}=x^{33}=x^{31+2}=x^2. $$ Taking the quotient by $x^{31}-1$ does exactly that: equating $x^{31}$ with $1$, $x^{32}$ with $x$ et cetera.
On with the example. When we calculate the product $C_0(x)\cdot C_1(x)$ we get first $$ \begin{aligned} &\left(x^{16}+x^8+x^4+x^2+x\right) \left(x^{24}+x^{17}+x^{12}+x^6+x^3\right)\\ =&x^{40}+x^{33}+x^{32}+2 x^{28}+x^{26}+2 x^{25}+x^{22}+x^{21}+x^{20}+2 x^{19}+x^{18}+x^{16}+2 x^{14}+x^{13}+x^{11}+x^{10}+x^8+2 x^7+x^5+x^4. \end{aligned} $$ Reducing the exponents modulo $31$ gives $$ 2 x^{28}+x^{26}+2 x^{25}+x^{22}+x^{21}+x^{20}+2 x^{19}+x^{18}+x^{16}+2 x^{14}+x^{13}+x^{11}+x^{10}+x^9+x^8+2 x^7+x^5+x^4+x^2+x. $$ As the finishing touch we then reduce the coefficients modulo two, and get the result: $$ x^{26}+x^{22}+x^{21}+x^{20}+x^{18}+x^{16}+x^{13}+x^{11}+x^{10}+x^9+x^8+x^5+x^4+x^2+x. $$ You see that this is the sum $$C_0(x)+C_2(x)+C_5(x)$$ as promised.
I introduced a lot extra baggage. This does help if you do research with these things. I trust Golomb to have known his customers, and done what he viewed best pedagogically.