Let $p \equiv 3 \pmod{4}$ be prime, and define $$ f(p)=\sum_{m=1}^{\frac{p-1}{2}}\sum_{n=1}^{\frac{p-1}{2}}\left(\frac{m+n-\frac{p+1}{4}}{p}\right), $$ where $\left(\frac{a}{p}\right)$ denotes the Legendre symbol (i.e. +1,-1,0, depending on whether $a$ is a quadratic residue/non-residue/zero).
The construction of $f(p)$ is admittedly artificial and was essentially done by trial-and-error to create a "large" f(p). What I'm interested in showing (which I'm not even sure is true) is that $$ f(p) = \Theta(p \sqrt{p}). $$
Now, I've programmatically computed $f(p)$ for all $p < 50000$, and it is from this that I'm conjecturing $$ \frac{p \sqrt{p}}{6} \le f(p) \le \frac{p \sqrt{p}}{4}. $$ In particular, it is the lower bound that I'm actually interested in (though the upper bound might also be interesting in its own way), since I would have expected $f(p)$ to be close to zero.
More generally, are there any papers that research something similar to this? The closest things I can find are related to clique numbers of Paley graphs, or the smallest quadratic non-residue, but neither of these seem to relate directly to this question.
First step, we introduce a function to roughly describe the number of the term $\big(\frac{k}{p}\big)$. For convenience purpose, We define it on $[0,2\pi]$ as $$f(x)=\cases{x, & \text{$0\leqslant x\leqslant \pi/2$}\\ \pi-x, & \text{$\pi/2\leqslant x\leqslant3\pi/2$}\\ x-2\pi, & \text{$3\pi/2\leqslant x\leqslant 2\pi$}}.\\$$
It is straightforward to show that the Fourier series of $f(x)$ is $$\frac{4}{\pi}\sum_{n=1}^\infty\frac{(-1)^{n-1}}{(2n-1)^2}\sin((2n-1)x),$$ which convergent to $f(x)$ absolutely and uniformly.
Second step, we use the result of Gauss sum for prime $p\equiv 3\pmod{4}$ $$\sum_{k=0}^{p-1}\biggl(\frac{k}{p}\biggr)e^{2\pi ik/p}=i\sqrt{p},$$ which implies $$\sum_{k=0}^{p-1}\biggl(\frac{k}{p}\biggr)\sin(2\pi km/p)=\biggl(\frac{m}{p}\biggr)\sqrt{p}.$$
Thus we can calculate the sum \begin{align} S(p)&:=\frac{p}{2\pi}\sum_{k=0}^{p-1}\biggl(\frac{k}{p}\biggr)f(2\pi k/p)\\ &=\frac{2p}{\pi^2}\sum_{k=0}^{p-1}\sum_{n=1}^{\infty}\frac{(-1)^{n-1}}{(2n-1)^2}\biggl(\frac{k}{p}\biggr)\sin(2\pi k(2n-1)/p)\\ &=\frac{2p}{\pi^2}\sum_{n=1}^{\infty}\frac{(-1)^{n-1}}{(2n-1)^2}\sum_{k=0}^{p-1}\biggl(\frac{k}{p}\biggr)\sin(2\pi k(2n-1)/p)\\ &=\frac{2p\sqrt{p}}{\pi^2}\sum_{n=1}^{\infty}\frac{(-1)^{n-1}}{(2n-1)^2}\biggl(\frac{2n-1}{p}\biggr)\\ &=\frac{2p\sqrt{p}}{\pi^2}L(2,\chi),\\ \end{align} where $\chi$ is the primitive Dirichlet character modulo $4p$ with order $2$.
Let $p=4r-1$, note that we may rewrite $S(p)$ as follows $$ S(p):=\frac{p}{2\pi}\sum_{k=0}^{p-1}\biggl(\frac{k}{p}\biggr)f(2\pi k/p)=\sum_{k=0}^{r-1} k\biggl(\frac{k}{p}\biggr)+\sum_{k=r}^{3r-1} (p/2-k)\biggl(\frac{k}{p}\biggr)+\sum_{k=3r}^{4r-2} (k-p)\biggl(\frac{k}{p}\biggr).\\ $$
Third step, we show that the sum in your question equals $S(p)$. \begin{align} \sum_{m=1}^{2r-1}\sum_{n=1}^{2r-1}\biggl(\frac{m+n-r}{p}\biggr)&=\sum_{k=2-r}^{r-1}\sum_{n=1}^{k+r-1}\biggl(\frac{k}{p}\biggr)+\sum_{k=r}^{3r-2}\sum_{n=k-r+1}^{2r-1}\biggl(\frac{k}{p}\biggr)\\ &=\sum_{k=2-r}^{r-1}(k+r-1)\biggl(\frac{k}{p}\biggr)+\sum_{k=r}^{3r-2}(3r-1-k)\biggl(\frac{k}{p}\biggr)\\ &=\sum_{k=0}^{r-1}(k+r-1)\biggl(\frac{k}{p}\biggr)+\sum_{k=r}^{3r-2}(3r-1-k)\biggl(\frac{k}{p}\biggr)+\sum_{k=2-r}^{-1}(k+r-1)\biggl(\frac{k}{p}\biggr)\\ &=\sum_{k=0}^{r-1}(k+r-1)\biggl(\frac{k}{p}\biggr)+\sum_{k=r}^{3r-2}(3r-1-k)\biggl(\frac{k}{p}\biggr)+\sum_{k=3r+1}^{4r-2}(k-3r)\biggl(\frac{k}{p}\biggr)\\ &=\sum_{k=0}^{r-1}(k+r-1)\biggl(\frac{k}{p}\biggr)+\sum_{k=r}^{3r-1}(3r-1-k)\biggl(\frac{k}{p}\biggr)+\sum_{k=3r}^{4r-2}(k-3r)\biggl(\frac{k}{p}\biggr)\\ &=S(p)+(r-1)\sum_{k=0}^{r-1}\biggl(\frac{k}{p}\biggr)+\frac{2r-1}{2}\sum_{k=r}^{3r-1}\biggl(\frac{k}{p}\biggr)+(r-1)\sum_{k=3r}^{4r-2}\biggl(\frac{k}{p}\biggr)\\ &=S(p)+(r-1)\sum_{k=0}^{p-1}\biggl(\frac{k}{p}\biggr)+\frac{1}{2}\sum_{k=r}^{3r-1}\biggl(\frac{k}{p}\biggr)\\ %&=S(p)+\frac{r-1}{2}\biggl(\sum_{k=0}^{p-1}\biggl(\frac{k}{p}\biggr)-\sum_{k=0}^{p-1}\biggl(\frac{-k}{p}\biggr)\biggr)+\frac{1}{4}\biggl(\sum_{k=r}^{3r-1}\biggl(\frac{k}{p}\biggr)-\sum_{k=r}^{3r-1}\biggl(\frac{-k}{p}\biggr)\biggr)\\ &=S(p).\\ \end{align}
For the last step, we use the estimations $$L(2,\chi)\leqslant\sum_{n=1}^\infty\frac{1}{(2n-1)^2}=(1-2^{-2})\zeta(2)=\frac{\pi^2}{8},$$ $$L(2,\chi)=\prod_{\text{prime}~p}(1-\chi(p)p^{-2})^{-1}\geqslant\prod_{\text{prime}~p>2}(1+p^{-2})^{-1}=(1+2^{-2})\frac{\zeta(4)}{\zeta(2)}=\frac{\pi^2}{12}.$$ Then we get $$\frac{p\sqrt{p}}{6}\leqslant S(p)\leqslant\frac{p\sqrt{p}}{4}.$$ That is exactly what you conjectured.
I believe there are more general and deep theories on values of Dirichlet L-functions which may provide much better solution. However, I have little knowledge of that, so I just mention this field as a possibility for you.