  In this paper we study the sum of powers in the Gaussian integers
$\mathbf{G}_k(n):=\sum_{a,b \in [1,n]} (a+b i)^k$. We give an explicit formula
for $\mathbf{G}_k(n) \pmod n $ in terms of the prime numbers $p \equiv 3 \pmod
4$ with $p \mid \mid n$ and $p-1 \mid k$, similar to the well known one due to
von Staudt for $\sum_{i=1}^n i^k \pmod n$. We apply this formula to study the
set of integers $n$ which divide $\mathbf{G}_n(n)$ and compute its asymptotic
density with six exact digits: $0.971000\ldots$.
