Create matrix containing values in [0,1] where sum of all diagonals and anti-diagonals is fixed

Question

The problem I am facing sounds at first glance pretty simple. However, as very often, it seems more complicated than I first assumed:

I want to calculate a matrix $P = (p_{j,k}) \in \mathbb{R}^{n \times n}$, $n\in\mathbb{N}$, satisfying the following constraints:

1. $p_{j,k} \in [0,1]$,
2. $\sum_{j,k} \ p_{j,k} = 1$,
3. the sums of all $2n-1$ diagonals are fixed by $b_{1-n},\ldots,b_{n-1}\in (0,1)$, e.g. $\sum_{j} p_{j,j} = b_0$ or $\sum_{j} p_{j,j+1} = b_{-1}$, and finally
4. the sums of all $2n-1$ antidiagonals are fixed by $a_1,\ldots,a_{2n-1}\in (0,1)$.

Obviously, the problem is easy to solve for $n=2$, since the "corner elements" of the matrix are directly given from 3. and 4. However, for $n>2$, it becomes more difficult.

A simple linear algebra approach for $n=3$ leads to the problem of "solving" $Mp=c$ with $$ M = \left( \begin{array}{ccccccccc} 1 & & & 0 & 0 & 0 & 0 & 0 & 0 \\ & 1 & & 1 & & & 0 & 0 & 0 \\ & & 1 & & 1 & & 1 & & \\ 0 & 0 & 0 & & & 1 & & 1 & \\ 0 & 0 & 0 & 0 & 0 & 0 & & & 1\\ 0 & 0 & 0 & 0 & 0 & 0 & 1 & &\\ 0 & 0 & 0 & 1 & & & & 1 & \\ 1 & & & & 1 & & & & 1\\ & 1 & & & & 1 & 0 & 0 & 0 \\ & & 1 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1& 1 & 1& 1 & 1 & 1 & 1& 1&1 \\ \end{array} \right), \quad p = \left( \begin{array}{c} p_{1,1} \\ p_{1,2} \\ p_{1,3} \\ p_{2,1} \\ p_{2,2} \\ p_{2,3} \\ p_{3,1} \\ p_{3,2} \\ p_{3,3} \end{array} \right), \quad c = \left( \begin{array}{c} a_1 \\ a_2 \\ a_3 \\ a_4 \\ a_5 \\ b_{-2} \\ b_{-1} \\ b_{0} \\ b_{1} \\ b_{2} \\ 1 \end{array} \right) $$ where $p \in \mathbb{R}^{n^2}$, $b \in \mathbb{R}^{4n-1}$ and $M$ is of size $(4n-1)\times n^2$ with entries in $\{0,1\}$. Note that $M$ is rank-deficient and the rows that specify the values for the "corner elements" are clearly visible. Dropping these rows leaves $$ \tilde{M} = \left( \begin{array}{ccccccccc} & 1 & & 1 & & & 0 & 0 & 0 \\ & & 1 & & 1 & & 1 & & \\ 0 & 0 & 0 & & & 1 & & 1 & \\ 0 & 0 & 0 & 1 & & & & 1 & \\ 1 & & & & 1 & & & & 1\\ & 1 & & & & 1 & 0 & 0 & 0 \\ 1& 1 & 1&1 &1 & 1 & 1 & 1 & 1 \\ \end{array} \right), \quad \tilde{p} = \left( \begin{array}{c} p_{1,2} \\ p_{2,1} \\ p_{2,2} \\ p_{2,3} \\ p_{3,2} \\ \end{array} \right), \quad \tilde{c} = \left( \begin{array}{c} a_2 \\ a_3 \\ a_4 \\ b_{-1} \\ b_{0} \\ b_{1} \\ 1-a_1-a_5-b_{-2}-b_2 \end{array} \right) . $$ The question is, if there exists a solution (for $n \in \mathbb{N}$) for this problem satisfying all properties and how to calculate it (least-squares, SVD,...)? Or is there an approach different from the one I chose which is suitable to compute a solution iteratively? I have the intuition, that due to the somewhat special constraints, there is a close connection to probability ("values are positive, sum equals one") or graph (adjacency matrices) theory.

Answer 1

The given constraints are a system of linear inequalities, so you can find a feasible solution (or prove that none exists) by feeding these constraints to a linear program (LP) solver with some arbitrarily chosen objective function.