I have two symmetric $3\times 3$ matrices $A, B$. I am interested in solving the system
$$ A X = X B $$
Is there a way this is usually done? The matrices are not necessarily non-singular.
I have two symmetric $3\times 3$ matrices $A, B$. I am interested in solving the system
$$ A X = X B $$
Is there a way this is usually done? The matrices are not necessarily non-singular.
You can use vectorization and Kronecker products.
Assume we have the equation $$\bf AXB=C$$
$$({\bf B}^T \otimes {\bf A})\text{vec}({\bf X}) = \text{vec}({\bf AXB}) = \text{vec}({\bf C})$$
Now how can we modify it to get what we want?
As stated in the comment from Christoph, here is a direct approach solving $AX=XB$ given both $A$ and $B$ are square matrices of order $k$.
Proposition. Let $A=(a_{ij})_{k\times k}$, $B=(b_{ij})_{k\times k}$ and the unknown $X=(x_{ij})_{k\times k}$. We also let $$Y=(x_{11},x_{21},\cdots,x_{k1},x_{12},x_{22},\cdots,x_{k2},\cdots,x_{1k},x_{2k},\cdots,x_{kk})^{T}$$ be the column vector with dimension $k^2$. Then, the solution of $AX=XB$ is given by the solution of the homogeneous equation $QY=0$, where $Q$ is the block matrix $$Q=\begin{bmatrix} Q_{11}&Q_{12}&\cdots&Q_{1k}\\ Q_{21}&Q_{22}&\cdots&Q_{2k}\\ \cdots&\cdots&\cdots&\cdots\\ Q_{k1}&Q_{k2}&\cdots&Q_{kk}\\ \end{bmatrix}$$ with $Q_{ii}=A-b_{ii}I$ ($i=1,2,\dots,k$ and $I$ is the identity matrix of order $k$) and $Q_{ij}=-b_{ji}I$ ($i=1,2,\dots,k$).
The proof follows by direct multiplication of matrices. We get the system of linear equations about every elements in $X$ and write the equation as the form $QY=0$. Since there always has solutions for the homogeneous equation $QY=0$, we get all the solutions for $AX=XB$.
Example. Solve $AX=XB$, where $A=\begin{bmatrix}-4&1\\-9&2\end{bmatrix}$ and $B=\begin{bmatrix}-2&1\\-1&0\end{bmatrix}$.
From the proposition, we have $$Q=\begin{bmatrix} -2&1&1&0\\ -9&4&0&1\\ -1&0&-4&1\\ 0&-1&-9&2\\ \end{bmatrix}\xrightarrow{\text{row operation}} \begin{bmatrix} 1&0&4&-1\\ 0&1&9&-2\\ 0&0&0&0\\ 0&0&0&0\\ \end{bmatrix}=P.$$ Since $PY=0$ has the same solution to $QY=0$, we get $$Y=(-4u+v,-9u+2v,u,v)^T,\qquad \forall u,v\in\Bbb C.$$ Therefore, $$X=\begin{bmatrix} -4u+v&u\\ -9u+2v&v\\ \end{bmatrix}= u\begin{bmatrix}-4&1\\-9&0\end{bmatrix}+v\begin{bmatrix}1&0\\2&1\end{bmatrix},\qquad \forall u,v\in\Bbb C.$$
From the example we can see the steps are
express $Q$ in terms of $A$ and $B$, apply row operations to obtain $P$,
solve $Y$ for $PY=0$,
write down the elements in $X$ (corresponding to the elements in $Y$).
This is a partial answer, but for a generalized setting where $A,B\in M_n(\Bbb{K})$ for some field $\Bbb{K}$. Suppose that $A$ and $B$ are diagonalizable over $\Bbb K$. Let $a_1,a_2,\ldots,a_n$ be a complete set of eigenvectors of $A$ with eigenvalues $\lambda_1,\lambda_2,\ldots,\lambda_n$, and $b_1,b_2,\ldots,b_n$ a complete set of eigenvectors of $B^t$ with eigenvealues $\mu_1,\mu_2,\ldots,\mu_n$. Then, the matrices $C_{i,j}=a_ib_j^t$ form a basis of $M_n(\Bbb{K})$.
Consider the linear map $T:M_n(\Bbb{K})\to M_n(\Bbb{K})$ given by $T(X)=AX-XB$. If $X=\sum_{i,j}x_{i,j}C_{i,j}$, then $$T(X)=\sum_{i,j}x_{i,j}(\lambda_i-\mu_j)C_{i,j}.$$ In particular, $T(X)=0$ if and only if $x_{i,j}=0$ whenever $\lambda_i\neq \mu_j$. Therefore, $$\ker T=\operatorname{span}\big\{C_{i,j}:\lambda_i=\mu_j\big\}.$$
In your case, the matrices $A$ and $B$ are real symmetric (you didn't specify the field, so I guess the field is the reals), and so they are diagonalizable. Therefore, you can use my approach.
However, if the field is $\Bbb{C}$, then a symmetric matrix may not be diagonalizable. I do not know the answer if $A$ or $B$ is not diagonalizable. But the wiki link given by Tianlalu may help.
After I made some calculations, I found this result. Assume that the characteristic polynomials of $A$ and $B$ split into linear factors over $\Bbb{K}$. Suppose that $A$ has $k$ Jordan blocks of sizes $r_1,r_2,\ldots,r_k$ with eigenvalues $\lambda_1,\lambda_2,\ldots,\lambda_k$, and suppose that $B$ has $l$ Jordan blocks of sizes $s_1,s_2,\ldots,s_l$ with eigenvalues $\mu_1,\mu_2,\ldots,\mu_l$. That is, we have the following proposition.
Proposition: Let $T:M_n(\Bbb{K})\to M_n(\Bbb{K})$ be defined by $T(X)=AX-XB$, where $A$ and $B$ are fixed elements of $M_n(\Bbb{K})$. Then, the linear map $T$ is diagonalizable over $\Bbb K$ if and only if both $A$ and $B$ are diagonalizable over $\Bbb{K}$.
Write $(a_i^1,a_i^2,\ldots,a_i^{r_i})$ for a generalized eigenvector sequence of $A$ in the $i$th block, that is,
$$Aa_i^p=\lambda_ia_i^p+a_i^{p+1}$$
for $p=1,2,\ldots,r_i$, with $a_i^{r_i+1}=0$. Similarly, write $(b_j^1,b_j^2,\ldots,b_j^{s_j})$ for a generalized eigenvector sequence of $B$ in the $j$th block, that is,
$$Bb_j^q=\mu_jb_j^q+b_j^{q+1}$$
for $q=1,2,\ldots,s_j$, with $b_j^{s_j+1}=0$. Then, for fixed $i=1,2,\ldots,k$ and $j=1,2,\ldots,l$, the span of $C_{i,j}^{p,q}=a_i^p(b_j^q)^t$ with $p$ and $q$ ranging from $\{1,2,\ldots,r_i\}$ and $\{1,2,\ldots,s_j\}$ is a direct sum of generalized eigenspaces of $T$ with eigenvalue $\lambda_i-\mu_j$.
This is because
$$\big(T-(\lambda_i-\mu_j)\big)^hC_{i,j}^{p,q}=\sum_{r=0}^h(-1)^r\binom{h}{r}C_{i,j}^{p+h-r,q+r},$$
where $C_{i,j}^{p,q}=0$ if $p>r_i$ or $q>s_j$. Therefore, if $h>r_i-p+s_j-q$, then $$\big(T-(\lambda_i-\mu_j)\big)^hC_{i,j}^{p,q}=0.$$ In characteristic $0$, there are $m_{i,j}=\min\{r_i,s_j\}$ generalized eigenspaces and the dimensions of the generalized eigenspaces are
$$r_i+s_j-1,r_i+s_j-3,r_i+s_j-5,\ldots,r_i+s_j+1-2m_{i,j}.$$
This shows that $$\operatorname{span}\{C^{r_i,s_j}_{i,j}:\lambda_i=\mu_j\big\}\subseteq\ker T\subseteq \operatorname{span}\{C_{i,j}^{p,q}:\lambda_i=\mu_j\big\}.$$ For $X\in M_n(\Bbb{K})$, we write $X=\sum_{i,j,p,q}x_{i,j}^{p,q}C_{i,j}$. If $X\in \ker T$, then $x_{i,j}^{p,q}=0$ whenever $\lambda_i\neq \mu_j$, and so $$T(X)=\sum_{i,j,p,q}x_{i,j}^{p,q}(C_{i,j}^{p+1,q}-C_{i,j}^{p,q+1})=\sum_{i,j}^{p,q}(x_{i,j}^{p-1,q}-x_{i,j}^{p,q-1})C_{i,j}^{p,q},$$ with $x_{i,j}^{0,q}=x_{i,j}^{p,0}=0$. This is as far as I can go, but it tells me that, unless both $A$ and $B$ are diagonalizable, we have $$\operatorname{span}\{C^{r_i,s_j}_{i,j}:\lambda_i=\mu_j\big\}\subsetneq \ker T\subsetneq \operatorname{span}\{C_{i,j}^{p,q}:\lambda_i=\mu_j\big\}.$$