I am trying to find a nice and efficient way to approach the following problem: I need to solve certain type of equations (see below for an example) for a set of square matrices, which have at most one non-zero element in each row and the sum of the $j$-th row across all matrices is equal to a constant.
To give an example, consider a set $\mathcal{A}$ of four matrices labeled as $\mathcal{A}: \{A_1,A_2,A_3,A_4\}$ obeying the properties mentioned. I am looking to find a set of two matrices $\mathcal{X}:\{X_1,X_2\}$ satisfying $$ A_3 \otimes X_2 = A_2\otimes X_1 \,, $$ subject to the constraint that the $X_i$'s have at most one non-zero element in each line and that the sum $X_1+X_2$ is a row stochastic matrix. In fact, it is the last two constraints that I do not know how to enforce using Solve or Reduce.
In general, the set $\mathcal{A}$ will have $n^2$ matrices while the set $\mathcal{X}$ will have $n$, and the equations that need to be solved become more complicated (multiple terms on each side and multiple equations). The example mentioned is just one of the simplest cases and could be worked out by hand for low dimensional matrices but I'm looking for a way to write a code that can handle bigger sets (ideally analytically). Of course solutions need not exist or could be multiple, but I am just interested in finding one, if it exists.
Any ideas how to elegantly enforce such constraints? Can I tell Mathematica in advance that matrices $X_i$ can only have one non-zero element on each row?