Question
I would like to know the most general and optimal form that the $n\times n$ square matrices $A$ and $B$ must take s.t., the matrix $AB^T$ is symmetric. (e.g., $A\in GL(n)$, $B=\lambda I$ would work, but this is clearly not the most general).
We can assume that $A$ and/or $B$ are invertible if necessary.
I think that the optimal number of parameters that parameterize $A$ and $B$ would be $ p_{opt}:=2n^2 - n(n-1)/2 = 3n(n-1)/2 $, (i.e., $2n^2$ to parameterize an arbitrary $A$ and $B$, minus $n(n-1)/2$ for their product to be symmetric). So ideally, the solution would have this optimal number of free parameters.
.
.
.
Examples
Example 1a: Letting $A=C\Lambda$ and $B=C\Pi$, where $C$ is an arbitrary square matrix and $\Lambda$ and $\Pi$ are commuting matrices (e.g., diagonal matrices), then clearly $AB^T-BA^T =0 $. Say that $\Lambda$ and $\Pi$ are diagonal, then this parameterization has $n(n+2)<p_{opt}$ free parameters. This is fine, but probably not the most general parameterization.
Example 1b: Letting $\Lambda$ and $\Pi$ be symmetric commuting matrices means that they are simultaneously diagonalizable: $\Lambda = Q\, diag(a) \,Q^T$ and $\Pi = Q \,diag(b) \,Q^T$ for $Q\in O(n)$ orthogonal. I think that $Q$ can be parameterized by $n(n-1)/2$ parameters (as the Lie algebra of $O(n)$ are the skew-symmetric matrices, which are have $n(n-1)/2$ free parameters). Therefore this parameterization has $n^2 + 2n + n(n-1)/2 = 3n(n+1)/2>p_{opt}$. Meaning it is not optimal.
Example 2: Let $A = I + S_3S_2$ and $B = S_1 + S_3 + S_3S_2S_1$, where $S_i$ are arbitrary symmetric matrices. Then this satisfies $AB^T-BA^T=0$, but has $3n(n+1)/2>p_{opt}$. Again not optimal, but probably general. We could however assert that $S_i$ have zeros on the diagonal, making the number of free parameters $3n(n-1)/2=p_{opt}$, meaning this is the optimal number of parameters. But does it span the most general set of matrices that satisfy $AB^T=BA^T$? This is yet to be proven and I wouldn't know how to prove that. Maybe one could write this out each $S_i = \sum \alpha_{ij} E_j$, where $E_j$ are the matrix basis elements for symmetric matrices, then show independence of the $\alpha_{ij}$'s?
Remark: $p_{opt}=3n(n-1)/2$ is three times the number of free parameters of a skew-symmetric matrix, which suggests that possibly this could be parameterized by three skew symmetric matrices?