Find the number of subsets $\{a_1,\cdots, a_k\}$ of $\{1,\cdots, p\}$ so that $a_1+a_2+\cdots + a_k\equiv 0\mod p.$ 
Let $p$ be an integer exceeding 1 and suppose $1\leq k\leq p-1.$ Find the number of $k$-element subsets $\{a_1,\cdots, a_k\}$ of $\{1,\cdots, p\}$ so that $a_1+a_2+\cdots + a_k\equiv 0\mod p.$

I've only managed to solve the above problem when p is prime.
If $p=2,$ then the number of 1-element subsets of $\{1,2\}$ with sum congruent to 0 mod 2 is $1.$ Henceforth assume $p>2$. Then note that for every $1\leq k\leq p-1$, we can find a k element subset whose sum is congruent to 0 mod p. Indeed,

*

*if k is odd and $k=2y+1$ for some $0\leq y\leq p/2 - 1$, we can take the sum $p+(p-1)+1 + \cdots (p-y) + y$ as $p-i, i < p, p-i\neq i$ for $i < p$, $p-y > y$ since $y < p/2,$ and $y$ is the highest second element among the pairs $(p-1,1),(p-2,2),\cdots, (p-y,y)$ while the first terms are strictly decreasing.

*If $k=2y$ is even, then we have the sum $(p-1 + 1) + \cdots + (p-y+y),$ which is a sum of distinct elements by similar reasoning to above as $y \leq (p-1)/2.$
Now given such a subset, say $B := \{a_1,\cdots, a_k\},$ define for $0\leq i\leq p-1$ the subset $B_i$ by
$$B_i := \{f(a_1+i), \cdots, f(a_k+i)\}$$
where $f(x)$ is the unique element $y$ of $\{1,\cdots, p\}$ so that $x\equiv y\mod p.$
Then note that all the elements $a_j + i$ are distinct modulo $p$ and since $f(a_j+i)\equiv a_j+i\bmod p$ for all $a_j+i$, the terms $f(b)$ are distinct where $b$ ranges over $B.$
Also, $\sum_{a\in B_i} a \equiv a_1+\cdots + a_k + ik \equiv ik\mod p.$ Each residue $ik\mod p$ is distinct because if $ik \equiv jk$ for some $0\leq i < j < p,$ then $p$ would divide $(i-j)k$ and hence $p$ would divide $(i-j)$. Since $p$ is coprime to $k$, this implies $i-j$ is divisible by $p$, which is a contradiction since $0<|i-j|<p.$ Hence each $B_i$ is distinct.
Now note that each $k$-subset of $\{1,\cdots, p\}$ belongs to exactly one $B_i, 0\leq i\leq p-1$ for some $B$. Indeed, let for each $B,$ the set of $B_i$'s defined above be the cycle of B. Note that if two cycles intersected, then they would be the same cycle. Note that here it's crucial that $p$ is prime (the claims above could still work if $p$ was replaced by $\dfrac{p}{\gcd(p,k)}$). Suppose the cycles of $B=\{b_1,\cdots, b_k\}$ and $A=\{a_1,\cdots, a_k\}$ intersected. Then there exist $s,t$ so that  $A_t = B_s$. So $a_i + t \equiv b_{\sigma(i)} + s$ for some permutation $\sigma$ of $\{1,\cdots, k\}.$ Then $a_i\equiv b_{\sigma(i)} + s-t$ for all i. But then A equals $B_{s-t}$ ($B_{s-t}$ may not be defined if we replaced $p$ with $\dfrac{p}{\gcd(p,k)}$ in our construction of the $B_i$'s). Hence for all $0\leq i\leq p-1, A_i$ equals $B_{s-t + i\mod p}$ and since the cycle elements $B_i$ are distinct for all i, this shows that the map $A_i\mapsto B_{(s-t+i)\mod p}$ defines a bijection between the cycle of A and the cycle of B. Hence distinct cycles are disjoint. But in each cycle, each of the $B_i$'s has a distinct sum modulo p as proved above, so exactly one $B_i$ has sum $0$ mod p.
Thus, the required answer is $\dfrac{{p\choose k}}{p}.$

I'm not sure if generating functions could be useful when p is not prime.
Let $\xi$ be a primitive p-th root of unity (note that these are precisely the numbers of the form $e^{2\pi i k/p}$ where k is coprime to p).

Consider the product $\prod_{j=1}^p (1+\xi^j).$ Note that when this product is written in the form $b_1 \xi + b_2\xi^2+\cdots + b_p \xi^p$ (which is possible since $\xi^p=1\Rightarrow \xi^{p+k} = \xi^k$ for all k), $b_p$ equals the number of subsets of $x$ with sum divisible by p.
 A: Consider $$f(x,y)=(1+xy)(1+x^2y)\cdots(1+x^py)=\sum_{n,k}a_{n,k}x^ny^k.$$
Then $a_{n,k}$ counts the number of subsets of $\{1,2,\dots,p\}$ of size $k$ and sum $n.$
Then if $\zeta$ is a primitive $p$th root of unity:
$$g(y)=\dfrac{\sum_{j=0}^{p-1}f(\zeta^j,y)}{p}=\sum_k y^k\sum_{p\mid n} a_{n,k}$$
Now, $$f(\zeta^i,y)=(-y)^p\prod_{j=1}^{p}\left(-\frac1y-\zeta^{ij}\right).$$
When $\gcd(i,p)=d,$ then this becomes:
$$f(\zeta^i,y)=(-y)^p\left(\left(-\frac1y\right)^{p/d}-1\right)^d=\left(1-(-y)^{p/d}\right)^d.$$
there are, for each $d\mid p,$ exactly $\phi(p/d)$ values $i$ with $\gcd(i,p)=d.$ So we get: $$\begin{align}g(y)&=\frac1p\sum_{d\mid p}\phi(p/d)\left(1-(-y)^{p/d}\right)^{d}\\&=\frac1p\sum_{d\mid p}\phi(d)\left(1-(-y)^{d}\right)^{p/d}\end{align}$$
So the coefficient of $y^k$ in $g(y)$ is: $$\frac1p\sum_{d\mid \gcd(p,k)}(-1)^{k+k/d}\phi(d)\binom{p/d}{k/d}$$
This agrees with your result for $p$ prime, and more generally, if $\gcd(k,p)=1,$ you get $\frac1p\binom pk.$
But I'd check this further before treating this as definitive - I might have made an algebra error.

Examples for checking
If $p=4,k=2,$ this gives: $$\frac14\left(\binom{4}{2}-\binom21\right)=1,$$ and there is only one subset of size $2$ with sum a multiple of $4.$

When $p=k,$ we know the value should be $1$ when $p$ is odd and $0$ when $p$ is even. The above formula becomes:
$$\frac1p\sum_{d\mid p}(-1)^{p+p/d}\phi(d)\tag1$$
When $p$ is odd, $p+p/d$ is always even, and we get $1$ for $(1)$ because $\sum_{d\mid p}\phi(d)=p$
When $p=2^kp',$ for $p'$ odd and $k>0$ then you get:
$$\sum_{d'\mid p'}\phi(d')\left(-\phi(2^k)+\sum_{i=0}^{k-1}\phi(2^i)\right)\tag2$$
But  $$\begin{align}\phi(1)+\phi(2)+\cdots \phi(2^{k-1})&=1+(2-1)+(2^2-2)+\cdots +(2^{k-1}-2^{k-2})\\&=2^{k-1}=\phi(2^k), \end{align}$$ so each term in the sum $(2)$ is zero.
