Given $M$ a real symmetric positive-definite matrix, I would like to characterise all matrices $A$ such that $A^\top M A=M$. Note that the question of finding $A$ solutions to $A^\top M A=M$ for all positive-definite matrices $M$ has already been answered here.
Necessarily $\det(A)^2=1$.
Reducing the problem to diagonal p.q. matrices $M$ is symmetric, so diagonalisable in an orthonormal basis: $M=Q^\top D Q$ with $D$ diagonal and $QQ^\top=I$. Then,
\begin{align} &A^\top Q^\top D Q A=Q^\top D Q \\ &(QA)^\top D(QA)=Q^\top D Q \\ &(QAQ^\top)^\top D (QAQ^\top)=D \end{align}
so I think the problem can be simplified in
$$\text{Find the matrices $B$ such that $B^\top D B=D$}$$ where $D=\operatorname{diag}(d_1,\dots,d_n)$, $d_i>0$.
and then recover $A$ from $A=Q^\top B Q$ ($Q$ and $D$ are known since $M$ is given). The $d_i>0$ stems for $M$ positive-definite.
Trying to solve this problem The $2^n$ matrices $B=\operatorname{diag}(\pm 1,\dots, \pm 1)$ are obvious solutions. What about the other ones? Writing $B^\top D B=D$ in components, the term $(i,j)$ is given by $$\sum_k B_{k,i}d_k B_{k,j}=d_i\delta_{ij}$$ In particular, $$\sum_k B_{k,i}^2d_k =d_i.$$ However, I don't see any clear conclusions using these equalities.
Question: How can to find/characterise the $B$ such that $B^\top D B=D$?
The matrix $D$ has a positive square root, name it $E$, which is a diagonal matrix with entries $\sqrt{d_i}>0$. Hence $E^{\top}=E$ and $(E^{-1})^{\top}=E^{-1}$.
From $B^{\top}DB=D$ we obtain $E^{-1}B^{\top}E=EB^{-1}E^{-1}$, hence $$ (EBE^{-1})^{\top}=(E^{-1})^{\top}B^{\top}E^{\top}=EB^{-1}E^{-1}=(EBE^{-1})^{-1}. $$ So $EBE^{-1}$ is an orthogonal (unitary) matrix.
Starting from any orthogonal (unitary) matrix $C$ we set $B:=E^{-1}CE$. Then $B^{\top}=E^{\top}C^{\top}(E^{-1})^{\top}=EC^{\top}E^{-1}$ and so $$ B^{\top}DB=B^{\top}EEB=EC^{\top}E^{-1}EEE^{-1}CE=EC^{\top}CE=E^2=D. $$