Speaking as a non-expert, I should think that for n much larger than, say, log(p) (or maybe even 4), the probability should approach 1/2. Are you interested in the answer for n << p, or p << n, or is there some other relationship between n and p which would make answering the question easier? Gerhard "Email Me About System Design" Paseman, 2011.07.05
–
Gerhard PasemanJul 5 '11 at 22:34

Also, I assume p is prime from your title, but it would help to call that out in the body of the question. Gerhard "Email Me About System Design" Paseman, 2011.07.05
–
Gerhard PasemanJul 5 '11 at 22:37

Yes, $p$ prime — and odd — is implicit both in the title and in the use of "quadratic residue".
–
Noam D. ElkiesJul 5 '11 at 22:53

This is my intuition as well. n should be at most polynomial in log(p) and p is meant to be a prime. Is there a standard theorem for this?
–
user16203Jul 5 '11 at 22:54

As you can see, if $p \rightarrow \infty$ then already $n=3$ is enough to get near-equidistribution (and $n=2$ fails only because 0 is over- or under-represented depending on whether $p$ is $+1$ or $−1 \bmod 4$).
–
Noam D. ElkiesJul 6 '11 at 2:36

4 Answers
4

This probability can be calculated exactly, and indeed it approaches $1/2$ rather quickly — more precisely, for each $p$ it approaches the fraction $(p-1)/(2p)$ of quadratic residues $\bmod p$. This can be proved by elementary means, but perhaps the nicest way to think about it is that if you choose $n$ numbers $a_i$ independently and sum $a_i^2 \bmod p$, the resulting distribution is the $n$-th convolution power of the distribution of a random single square — so its discrete Fourier transform is the $n$-th power of the D.F.T., call it $\gamma$, of the distribution of $a^2 \bmod p$. For this purpose $\gamma$ is normalized so $\gamma(0)=1$. Then for $k \neq 0$ we have $\gamma(k) = (k/p) \gamma(1)$ [where $(\cdot/p)$ is the Legendre symbol], and
$$
p \gamma(1) = \sum_{a \bmod p} \exp(2\pi i a^2/p),
$$
which is a Gauss sum and is thus a square root of $\pm p$. It follows that $|\gamma(k)| = p^{-1/2}$, from which we soon see that each value of the convolution approaches $1/p$ at the exponential rate $p^{-n/2}$, and the probability you asked for approaches $(p-1)/(2p)$ at the same rate.

As noted above, this result, and indeed the exact probability, can be obtained by elementary means, yielding a (known but not well-known) alternative proof of Quadratic Reciprocity(!). But that's probably too far afield for the present purpose.

From your answer I infer 0 does not count as a quadratic residue mod p. (I definitely agree it is not as interesting as a quadratic residue, but should it be really ostracized from the set of squares?) Gerhard "Supports The Rights of Zero" Paseman, 2011.07.05
–
Gerhard PasemanJul 5 '11 at 23:12

5

While 0 is indeed a square, it does not count as a "quadratic residue" (nor as a "quadratic nonresidue" — though this term should really have been "nonquadratic residue"). Properly 0 should count as half square and half non-square (think about the number of square roots), and then the limit would really be 1/2.
–
Noam D. ElkiesJul 5 '11 at 23:55

5

For quadratic reciprocity: let $n$ be an odd prime prime $l \neq p$, and $N$ the number of solutions of $\sum_{i=1}^l a_i^2 = 1 \bmod p$. Cyclic permutation of the coordinates has two fixed points or none according as $N \equiv 2$ or $0 \bmod l$. From the formula for $N$ it soon follows that $(l/p) \equiv {p^*}^{(l-1)/2} \bmod l$ where $p^* = \gamma(1)^2 = p$ or $-p$ according as $p \equiv 1$ or $-1 \bmod 4$. By Legendre's formula this means $(l/p) = (p^*/l)$. Exercise: modify this to determine $(2/p)$ from the count of solutions of $x_1^2 + x_2^2 \equiv 1 \bmod p$. [Original source?]
–
Noam D. ElkiesJul 6 '11 at 0:01

1

:-) If no reference turns up here I'll post a reference request as a question... It worked quite well for my one previous M.O. query (on the integral for $\frac{22}{7} - \pi$).
–
Noam D. ElkiesJul 6 '11 at 3:17

5

This idea can also be used to get at a few memorable cases of higher reciprocity. For example: if $p \equiv 1 \bmod 3$, and you already know that the number of solutions of $x^3 + y^3 = 1 \bmod p$ is $p - 2 - a$ where $4p = a^2 + 27 b^2$ for some integers $a,b$, then it follows at once that $2$ is a cubic residue iff $a$ is even, that is, iff $p$ itself can be written as $a^2 + 27b^2$. [NB the count is $p-2-a$, not the familiar $p+1-a$, because for this purpose we must exclude the three points at infinity.]
–
Noam D. ElkiesJul 6 '11 at 14:41

Let
$$S=\sum_{a_1=0}^{p-1}\dots\sum_{a_n=0}^{p-1}\sum_{t=0}^{p-1}\sum_{m=0}^{p-1}e^{2\pi im(a_1^2+\cdots+a_n^2-t^2)/p}
$$
The innermost sum is $p$ if $a_1^2+\cdots+a_n^2-t^2\equiv0\pmod p$ and zero otherwise,
so $S$ counts $2p$ whenever $a_1^2+\cdots+a_n^2$ is a (nonzero) quadratic residue, $p$
whenever it's zero. On the other hand,
$$
S=\sum_{m=0}^{p-1}\sum_{a_1=0}^{p-1}\dots\sum_{a_n=0}^{p-1}\sum_{t=0}^{p-1}e^{2\pi im(a_1^2+\cdots+a_n^2-t^2)/p}
$$
so
$$
S=p^{n+1}+\sum_{m=1}^{p-1}\sum_{a_1=0}^{p-1}\dots\sum_{a_n=0}^{p-1}\sum_{t=0}^{p-1}e^{2\pi im(a_1^2+\cdots+a_n^2-t^2)/p}
$$
so
$$
S=p^{n+1}+\sum_{m=1}^{p-1}\left(\left(\sum_{a_1=0}^{p-1}e^{2\pi ima_1^2/p}\right)\cdots\left(\sum_{a_n=0}^{p-1}e^{2\pi ima_n^2/p}\right)\left(\sum_{t=0}^{p-1}e^{2\pi imt^2/p}\right)\right)
$$
Each of those inner sums is a Gauss sum and known to equal $\sqrt p$ in modulus (more detail: the sum is ${m\overwithdelims()p}\sqrt{{-1\overwithdelims()p}p}$), so $|S-p^{n+1}|\le(p-1)p^{(n+1)/2}$. For $n\gt1$, the main term beats the error term, and you get a good estimate.

Here is a slightly different argument: Let $Q$ be a non degenerate quadratic form over $\mathbb{F}_q$ of rang $n$ and determinant $d$. Let
$A(n,d)=|\{x\in \mathbb{F}_q^n:Q(x)=0\}|$. The claim is that $A(n,d)=q^{n-1}+O(q^{n/2})$.
For $n>2$ we can write $Q(X)=Q_0(X_1,\ldots,X_{n-2}) +X_{n-1}X_n$, where $Q_0$ is a form of rank $n-2$
in the variables $X_1,\ldots,X_{n-2}$. This decomposition shows instantly that
that $A(n,d)=(2 q-1) A(-d,n-2) +(q-1) (q^{n-2}-A(-d,n-2))$. Proceeding by induction we get the estimate $A(n,d)=q^{n-1}+O(q^{n/2})$.
(The error term can be computed exactly using Gauss sums). Applying this to the forms
$X_1^2+\cdots +X_n^2-X_{n+1}^2$ and $X_1^2+\cdots+ X_n^2$ we get that the desired probability is $(A(n+1,-1)-A(n,1))/(2 q^n)=
(q-1)/(2q) +O(q^{-n/2})$.