Is there any available implementation of a fast method to solve a Sylvester equation:
$$AX + XB = C,$$
where, $A \in \mathbb{R}^{$n \times n$}$, $B \in \mathbb{R}^{$p \times p$}$ and $C \in \mathbb{R}^{$n \times p$}$.
In my case, $B$ is a symmetric and large matrix with $p >> n$.
I tried the solve_sylvester function from python (scipy.linalg), and the two functions sylvester and lyap from Matlab as well, but all of them are slow.
Thanks,