# Linear matrix equation involving $\sum_i A_i X B_i$

I am dealing with a non-standard linear matrix equation: $$$$X \quad+\quad \sum_{i=1}^r \Big( (A_i X B_i) + (A_i X B_i)^T + (B_i X A_i) + (B_i X A_i)^T \Big) \quad=\quad C$$$$ The symbols in capital letters, i.e. $$X, A_i, B_i, C$$, denote $$n \times n$$ real matrices. $$X$$ is the unknown that we wish to solve for. Additional considerations:

1. The number of terms $$r$$ is large.
2. $$C$$ is symmetric
3. $$A_i$$'s are all rank one matrices (not sure if this information is helpful for expressing $$X$$ though).

Is there an elegant way to express the solution form $$X$$?

Each summand in the summation sign is symmetrised and $$C$$ is symmetric. Hence $$X$$ must be symmetric and the equation is equivalent to $$X+\sum_i\left(A_iXB_i+B_i^TXA_i^T+B_iXA_i+A_i^TXB_i^T\right)=C,$$ which can be rewritten as (see Wikipedia) $$\left[I+\sum_i\left(B_i^T\otimes A_i+A_i\otimes B_i^T+A_i^T\otimes B_i+B_i\otimes A_i^T\right)\right]\operatorname{vec}(X)=\operatorname{vec}(C).$$ Call the matrix inside the pair of square brackets $$M$$. The equation $$M\operatorname{vec}(X)=\operatorname{vec}(C)$$ is solvable if and only if $$MM^+\operatorname{vec}(C)=\operatorname{vec}(C)$$, where $$M^+$$ denotes the Moore-Penrose pseudoinverse of $$M$$. In case it is solvable, $$\operatorname{vec}(X)=M^+\operatorname{vec}(C)$$ is always a solution.