He aquí un enfoque elemental. $\#(M\cap M_i)$ es igual al número de soluciones $y-x=i$ con $x,y\in M$ . Todo residuo cuadrático es el cuadrado de dos residuos distintos de cero, por lo que $\#(M\cap M_i)$ es igual a $1/4$ veces el número de soluciones de $a^2-b^2=i$ con $a,b\in\mathbb{F}_p^\times$ . Para $p>2$ podemos hacer el cambio biyectivo de variables $u:=a+b$ y $v:=a-b$ lo que demuestra que estamos contando el número de soluciones $uv=i$ con $u,v\in\mathbb{F}_p$ tal que $u\neq\pm v$ . Sin la restricción $u\neq\pm v$ hay exactamente $p-1$ (de hecho, dejemos que $u\in\mathbb{F}_p^\times$ sea arbitraria y ponga $v:=i/u$ ). Si $i$ es un residuo cuadrático, entonces hay precisamente $4$ soluciones con $u=\pm v$ de lo contrario no existen tales soluciones.
Resumiendo, $\#(M\cap M_i)$ es igual a $\frac{p-5}{4}$ o $\frac{p-1}{4}$ dependiendo de si $i$ es un residuo cuadrático o no.
Para $p\equiv 3\pmod{3}$ hay precisamente $2$ soluciones de $uv=i$ con $u=\pm v$ (porque precisamente uno de $i$ y $-i$ es un residuo cuadrático), por lo que en este caso $\#(M\cap M_i)$ siempre es igual a $\frac{p-3}{4}$ .
Añadido. En general, es fácil contar el número de soluciones de $ a_1x_1^2+\dots +a_kx_k^2=b $ en un campo finito. Para un tratamiento de $k=2$ (que es fácil de extender a $k>2$ ) véase mi respuesta anterior aquí .