Regular subsets of $\text{PSL}(2, q)$

This is not a complete solution of the problem, but only an outline of a possible strategy, along with some numerical evidence: Let $S$ be a subset of $G=\operatorname{PSL}(2,q)$. For $g\in G$ let $P_g$ be the permutation matrix of $g$, and let $J$ be the all-$1$-matrix of the same size. Then regularity of $S$ is equivalent to $\sum_{g\in S}P_g=J$. Thus the existence of a regular subset of $G$ is equivalent to the existence of a solution of \begin{equation} \sum_{g\in G}x_gP_g=J \end{equation} for $x_g\in\{0,1\}$.

If there is a regular subset, then we may assume that it contains $1$, so the remaining elements are fixed-point-free. Let $G^\star$ be the set of fixed-point-free elements from $G$ together with $1$. Then the existence of a regular subset is equivalent to the $0$-$1$-solvability of \begin{equation} \sum_{g\in G^\star}x_gP_g=J. \end{equation} Now it seems to be the case that this system of linear equations in the $x_g$ isn't even solvable over $\mathbb F_2$. I have verified that for all $q\le 89$.

How could one prove that in general? Let $V$ be the space of $(q+1)\times(q+1)$ matrices over $\mathbb F_2$. Then $V\times V\to\mathbb F_2$, $(U,V)\mapsto\operatorname{trace}(UV)$ is a nondegenerate bilinear form on $V$.

From that we get the following: If the system is not solvable, then there is a matrix $M$ which is orthogonal to the space of matrices generated by the left hand side, but is not orthogonal to $J$. The converse holds of course too.

Thus we are left to find $M$ such that $\operatorname{trace}(P_gM)=0$ for all $g\in G^\star$, and $\operatorname{trace}(JM)=1$.

One may even make some additional assumptions on $M$. For $h\in G$ set $M^h=P_h^{-1}MP_h$. If $H$ is a subgroup of $G$ of odd order, and $M$ fulfills the assumptions from above, then $\sum_{h\in H}M^h$ fulfills this assumption too.

Thus in the search for $M$ we may restrict to those matrices for which $M^h=M$ for all $h\in H$. If $q$ is a power of the prime $p$, and $H$ is a $p$-Sylow subgroup of $G$, then it seems to be that there are exactly $8$ choices for such an $M$, no matter what $q$ is.

A working family of matrices $M$ seems to be as follows: Use $\mathbb F_q\cup\{\infty\}$ to number the rows and columns of $M$. Then $M[i,j]=1$ if and only if either $j=\infty$, $i\in\mathbb F_q$, or $i,j\in\mathbb F_q$ are distinct and $i-j$ is a square.

To prove that this $M$ works reduces to a concrete question about $\mathbb F_q$ which shouldn't be hard to prove.

I would be curious to learn the background or the motivation of this question.

Remark: It would be interesting to see other methods for proving non-existence of regular sets of permutations. Some years ago I worked with Theo Grundhöfer and Gabor Nagy on such questions. In most cases, the non-solvability of the first displayed equation modulo a prime was the most successful technique. In those cases, one could even choose $M$ as a rank $1$ matrix, which allowed for a direct and easy to use argument (we called it contradicting subsets lemma, see Lemma 2 in https://arxiv.org/abs/1005.1598). The special features of your question seem to be: (a) One really has to work with the second displayed equation, the first one does not give a contradiction, and (b) For $q\ge9$ it seems that there is no matrix $M$ as above (with no assumption on $H$-symmetry) of rank $1$, so the contradicting subsets lemma does not work here.

This is an attempt to understand the very nice argument by Peter Mueller (completed by OP, Sean Eberhard, in the comments). I am not completely happy with it, and would like to encourage others to look for more conceptual explanations.

Denote $\Omega:=P^1\mathbb{F}_q=\mathbb{F}_q\cup \infty$, so that $|\Omega|=q+1$ and the group $G=PSL(2,q)$ acts on $\Omega$ by projective transformations. Assume that $S\subset G=PSL(2,q)$, $|S|=q+1$ is chosen so that $$\sum_{s\in S} \mathbb{1}_{s(i)=j}=1,\quad\forall i,j\in \Omega.$$ Then the same holds for the set $g_0S$ for any $g_0\in G$, thus we may suppose that $id\in S$, therefore other elements of $S$ do not have fixed points. Denote by $G^\star$ the set of fixed-point-free elements from together with 1, we have $S\subset G^\star$. Consider the function $M(i,j)$ on $\Omega\times \Omega$ defined as $$ M(i,j)=\begin{cases}1,\, \text{if}\,\, i=\infty,\, j\ne \infty\\ 1,\, \text{if}\,\, i,j\in \mathbb{F}_q,\chi(i-j)=1\\ 0,\, \text{otherwise}.\end{cases} $$ Here $\chi$ is a quadratic character of $\mathbb{F}_q$ (Legendre symbol if $q$ is prime). We get $$ \sum_{i,j\in \Omega,s\in S} M(i,j)\mathbb{1}_{s(i)=j}=\sum_{i,j\in \Omega} M(i,j)=q+q(q-1)/2 $$ is odd. Thus to get a contradiction it suffices to prove that $$ \sum_{i,j\in \Omega} M(i,j)\mathbb{1}_{s(i)=j} $$ is even for any fixed element $s\in G^\star$. For $s=id$ all summands are just zeroes. If $s$ does not have fixed point, than $$ \sum_{i,j\in \Omega} M(i,j)\mathbb{1}_{s(i)=j}=1+\sum_{i\in \mathbb{F}_q} M(i,s(i))=1+\left|i\in \mathbb{F}_q:\chi(s(i)-i)=1\right|. $$ Note that for fixed $\alpha\in \mathbb{F}_q$, the equation $$s(x)-x=\alpha\quad \quad (1)$$ is quadratic with respect to $x$, thus it has even number of roots unless it has a double root. We prove that there exists unique quadratic residue $\alpha$ for which (1) has a double root, the result would follow.

It may be proved by straightforward calculations, but let me give an argument if not conceptual, but at least almost calculations-free.

First of all, we use the following characterization of $PSL(2,q)$ (probably well known, but I did not see it before.)

Lemma. A projective transformation $s(x)=\frac{ax+b}{cx+d}, c\ne 0$, belongs to $PSL(2,q)$ if and only if the equation $s'(x)=1$ has two roots in $\mathbb{F}_q$.

Proof. Note that both properties are preserved when we replace $s(x)$ to $s(x)+C$ and $s(x+C)$ for constant $C\in \mathbb{F}_q$. So we may suppose that $a=0$ (subtract constant $a/c$ from $s(x)$) and also that $d=0$ (replace $s(x)$ by $s(x-d/c)$). So $s(x)=M/x$ for some constant $M$, and both conditions say that $-M$ is a non-zero square.

Next observation concerns the discriminants of the quadratic equation and critical values of corresponding functions. For the equation $f(x)=-x^2+2ax-b=0$, we all know that its discriminant $a^2-b$ equals $f(x_0)$, where $x_0=a$ is the critical value. But we deal with a slightly different equation $h(x):=s(x)-x=0$. After shifting the variable we get a function of the form $h_0(x)=-A/x+B-x$. Its discriminant is $D=B^2-4A$, and the critical values are $x_1=\sqrt{A},x_2=-\sqrt{A}$, so in our situation $A$ is a square by Lemma. Also we get $$D=h_0(x_1)h_0(x_2).$$

Since $h(x)=s(x)-x$ does not have roots, $D$ is not a square. Thus exactly one of two values of $h$ at critical points is a square. This value is the unique appropriate $\alpha$, as desired.