Let us consider two vectors $\mathbf{a}=(a_1,...,a_n)$ and $\mathbf{b}=(b_1,...,b_m)$ so that each quantity is an **integer** $a_i,b_j \in \mathbb{N}$. It represents for example supply and demand. Let $\mathbf{C}=(c_{ij})_{ij} \in \mathbb{R}_{+}^{n \times m}$ be a matrix representing the cost of transporting goods from $i$ to $j$.

My question is the following: can we guarantee that a solution of the Optimal transport problem between $\mathbf{a}$ and $\mathbf{b}$ is a matrix whose elements **are integer**?

More precisely I am interested in the solutions of the following Optimal Transport problem: \begin{equation} \underset{\boldsymbol\pi \in \Pi(\mathbf{a},\mathbf{b})}{\min} \sum_{i=1,j=1}^{n,m}c_{ij}\pi_{ij}=\underset{\pi \in \Pi(\mathbf{a},\mathbf{b})}{\min} \langle \mathbf{C},\boldsymbol\pi \rangle \end{equation}

where $\Pi(\mathbf{a},\mathbf{b})$ is the "transport" polytope defined by: \begin{equation} \Pi(\mathbf{a},\mathbf{b})=\{\boldsymbol\pi \in \mathbb{R}^{n \times m} \ | \pi_{ij}\geq 0; \ \forall j, \ \sum_{i=1}^{n} \pi_{ij}=b_j; \ \forall i, \ \sum_{j=1}^{m} \pi_{ij}=a_i\} \end{equation} In this way $\pi_{ij}$ is the amount of goods transferred from $i$ to $j$ so that the overall amount is preserved. (Remark: usually this problem is considered with histograms but here I am not normalizing the quantities.)

We can solve this problem using the simplex algorithm for example.

More formally my question is: can we guarantee that there exists an optimal transport plan $\boldsymbol\pi^{*}$ of the previous problem such that $\pi^{*}_{ij}$ is an integer?

- My hint:

The problem is a standard linear program (LP) as it can be written: \begin{equation} \underset{\begin{smallmatrix}\mathbf{p} \in \mathbb{R}_{+}^{nm} \\ \mathbf{A}\mathbf{p}=\begin{pmatrix} \mathbf{a} \\ \mathbf{b} \end{pmatrix} \end{smallmatrix}}{\min} \mathbf{c}^{\top} \mathbf{p} \end{equation}

by vectorizing all the quantities, *i.e.* $\mathbf{c}$ is $nm$ vectors equal to the stacked columns contains in $\mathbf{C}$ and same for $\mathbf{p}$ which corresponds to $\boldsymbol\pi$. In this form the constraint matrix writes:
\begin{equation}
\mathbf{A}=\begin{pmatrix}\mathbf{1}_n^{\top} \otimes I_m \\ I_n \otimes \mathbf{1}_m^{\top} \end{pmatrix} \in \mathbb{R}^{(n+m)\times nm}
\end{equation}

where $\otimes$ is the Kronecker product, $I_n$ the identity matrix of size $n$ and $\mathbf{1}_n$ is a vector of ones of size $n$. This problem is a relaxation of the Integer Linear Programming (ILP) problem where elements of $\mathbf{p}$ are constrained to be integers (which is want we want). Solving (ILP) is NP-hard in general and we can use the (LP) relaxation and round the entries but the solution may not be optimal, and even not feasible.

However as written in [1] (see also https://en.wikipedia.org/wiki/Integer_programming) elements of an optimal solution of the (LP) are integers when $\mathbf{A}$ is totally unimodular, that is every square submatrix of $\mathbf{A}$ has determinant $0,+1$ or $-1$. Is $\mathbf{A}$ totally unimodular in this case?

[1] Alexander Schrijver, Theory of Linear and Integer Programming, John Wiley & Sons, 1986