# Finding Hall blockers by matrix scaling

Koyo Hayashi111Department of Information and System Engineering, Chuo
University, Tokyo 112-8551, Japan<EMAIL_ADDRESS>Hiroshi
Hirai222Corresponding author, Department of Mathematical Informatics, Graduate
School of Information Science and Technology, The University of Tokyo, Tokyo,
113-8656, Japan<EMAIL_ADDRESS>
###### Abstract

For a given nonnegative matrix $A=(A_{ij})$, the matrix scaling problem asks
whether $A$ can be scaled to a doubly stochastic matrix $XAY$ for some
positive diagonal matrices $X,Y$. The Sinkhorn algorithm is a simple iterative
algorithm, which repeats row-normalization $A_{ij}\leftarrow
A_{ij}/\sum_{j}A_{ij}$ and column-normalization $A_{ij}\leftarrow
A_{ij}/\sum_{i}A_{ij}$ alternatively. By this algorithm, $A$ converges to a
doubly stochastic matrix in limit if and only if the bipartite graph
associated with $A$ has a perfect matching. This property can decide the
existence of a perfect matching in a given bipartite graph $G$, which is
identified with the $0,1$-matrix $A_{G}$. Linial, Samorodnitsky, and Wigderson
showed that a polynomial number of the Sinkhorn iterations for $A_{G}$ decides
whether $G$ has a perfect matching.

In this paper, we show an extension of this result: If $G$ has no perfect
matching, then a polynomial number of the Sinkhorn iterations identifies a
Hall blocker—a certificate of the nonexistence of a perfect matching. Our
analysis is based on an interpretation of the Sinkhorn algorithm as
alternating KL-divergence minimization (Csiszár and Tusnády 1984, Gietl and
Reffel 2013) and its limiting behavior for a nonscalable matrix (Aas 2014). We
also relate the Sinkhorn limit with parametric network flow, principal
partition of polymatroids, and the Dulmage-Mendelsohn decomposition of a
bipartite graph.

Keywords: Matrix scaling, perfect matching, Hall blockers, alternating
minimization, KL-divergence, polymatroids

##  1 Introduction

For a given nonnegative matrix $A=(A_{ij})$, the (doubly stochastic) matrix
scaling problem [22] asks whether $A$ can be scaled to a doubly stochastic
matrix $XAY$ for some positive diagonal matrices $X,Y$. The Sinkhorn algorithm
(also called RAS algorithm or IFP procedure) is a simple iterative algorithm,
which repeats row normalization $A_{ij}\leftarrow A_{ij}/\sum_{j}A_{ij}$ and
column normalization $A_{ij}\leftarrow A_{ij}/\sum_{i}A_{ij}$ alternatively.
See Idel’s survey [16] for rich literature of matrix scaling and the Sinkhorn
algorithm. By this algorithm, $A$ converges to a doubly stochastic matrix in
limit if and only if the bipartite graph $G$ associated with $A$ has a perfect
matching [23]. Here $G$ is defined to have row and column indices $A$ as
vertices and have edges $ij$ with $A_{ij}>0$. This convergence property can
decide whether a given bipartite graph $G$ has a perfect matching. Consider
0,1-matrix $A_{G}$ defined by $(A_{G})_{ij}:=1$ for $ij\in E$ and zero for
other indices. Linial, Samorodnitsky, and Wigderson [17] showed that a
polynomial number of iterations of the Sinkhorn algorithm applied to $A=A_{G}$
decides whether $G$ has a perfect matching. The same result is independently
obtained by Gurvits and Yianilos [15]. Specifically, within $O(n^{2}\log n)$
iterations, if $\|{\bf 1}-A{\bf 1}\|_{2}^{2}<1/n$ holds after column
normalization, then there is a perfect matching. Otherwise, there is no
perfect matching. Here ${\bf 1}$ denotes the all-one vector, and $n$ is the
number of vertices in one of the color classes of $G$. Although the time
complexity of this algorithm is inferior than the standard augmenting path
algorithm, it is interesting in its simplicity of the algorithm description as
well as its conceptual difference from the standard one. Also it is connected
to the recent development on operator scaling; see the end of introduction.

In this paper, we address an extension of this result. Hall’s marriage theorem
says that $G$ has a perfect matching if and only if there is no Hall blocker—a
subset $X$ of vertices in one color class having neighbors $\Gamma(X)$ with
$|X|>|\Gamma(X)|$. The above algorithm outputs neither a perfect matching nor
a Hall blocker. The main result of this paper is to show that a polynomial
number of the Sinkhorn iterations can identify a Hall blocker if $G$ has no
perfect matching. Particularly, it determines a subset $X$ with maximum
$|X|-|\Gamma(X)|$, and hence determines the maximum matching number
$n-|X|+|\Gamma(X)|$ by König-Egervary theorem.

Our analysis is based on an interpretation of the Sinkhorn algorithm as an
alternating minimization for convex optimization. There are such two
interpretations. The first one is alternating minimization for capacity
$\inf_{x_{1},x_{2},\ldots,x_{n}>0}\prod_{i}\sum_{j}A_{ij}x_{j}/(x_{1}x_{2}\ldots
x_{n})$. The second one is alternating minimization for the minimum KL-
divergence between the row-normalized and column-normalized spaces of
matrices. The analysis of the above algorithm [15, 17] is based on the first
interpretation. If $G$ has no perfect matching, then the convex optimization
formulation of the (logarithmic) capacity has unbounded optimal value
$-\infty$. So it seems difficult from this analysis to retrieve information of
Hall blockers and obtain an iteration bound.

Our analysis is based on the second interpretation. This alternating
minimization approach and its convergence property can be described by a
framework of Csiszár and Tusnády [5]. This fact is explicitly noticed by Gietl
and Reffel [13]. In the case where $A$ is not (approximately) scalable, the
Sinkhorn algorithm does not converge. However it is oscillating in limit
between two matrices attaining the minimum KL-divergence between the row-
normalized space and column-normalized space; see [13]. Then Aas [1] revealed
a block-diagonalized structure of the oscillating limit and gave a polynomial
time procedure to determine to this block structure.

It turns out that this block structure includes information of Hall blockers.
We provide a sophisticated analysis on the oscillating limit based on the idea
of principal partition of polymatroids [9, 24]. Particularly, we reveal that
the block structure of the limit matrix equals a refined version of the
Dulmage-Mendelsohn decomposition [7] and is obtained via parametric network
flow. Also we give an explicit formula of the limit of the row-marginal
$A^{*}{\bf 1}$ of the column-normalized limit $A^{*}$. The block-diagonal form
is obtained by arranging $A$ with respect to the decreasing order of
$A^{*}{\bf 1}$. We provide an estimate of convergence to $A^{*}{\bf 1}$. After
a polynomial number of the Sinkhorn iterations, the row-marginal $A{\bf 1}$
after column-normalization is sufficiently close to $A{\bf 1}^{*}$. By sorting
$A{\bf 1}$ we can identify a vertex subset $X$ with maximum $|X|-|\Gamma(X)|$,
which is a Hall blocker if $G$ has no perfect matching. The estimated number
$O(n^{7}\log n)$ of required iterations is not good, since our convergence
analysis seems rather loose. The convergence rate of the Sinkhorn algorithm
for the scalable case (convergence case) is well-studied; see [16, Section 7].
We could not find any existing work on the convergence analysis of the
Sinkhorn algorithm for the nonscalable case. We expect that further
convergence analysis will improve the iteration bound.

It is an interesting future research to extend the present results to operator
scaling [6, 11, 12, 14]— a quantum generalization of matrix scaling to totally
positive operators. The Sinkhorn algorithm is generalized to the operator
Sinkhorn algorithm (or Gurvitz algorithm) [14]. A recent important result [11]
is that this algorithm can decide the (doubly stochastic) scalability in
polynomial time, which is viewed as a generalization of the above-mentioned
perfect matching testing by matrix scaling. A Hall blocker, a certificate for
nonscalability in this setting, is a certain vector subspace (called an shrunk
subspace). It is a natural question whether the operator Sinkhorn algorithm
can find a Hall blocker for the nonscalable case. This may need to
characterize the limits of the operator Sinkhorn algorithm, which is raised as
an open problem in [12]. Although the operator Sinkhorn algorithm is
interpreted as alternating minimization for the capacity in this setting, a
divergence minimization interpretation is unknown; see [19] for a progress in
this direction. We hope that our results will be useful for tackling these
problems.

### Organization.

In Section 2, we summarize basic facts on matrix scaling. In Section 3, we
provide a polymatroid-based analysis on the limit of the Sinkhorn algorithm.
In Section 4, we present our algorithm and prove its validity.

### Notation.

For an positive integer $n$, let $[n]:=\\{1,2,\ldots,n\\}$. Let ${\mathbb{R}}$
denote the set of real numbers. For a vector $p\in{\mathbb{R}}^{n}$ and subset
$X\subseteq[n]$, let $p(X):=\sum_{i\in X}p_{i}$. For a subset $X\subseteq[n]$,
let ${\bf 1}_{X}\in{\mathbb{R}}^{n}$ denote the 0,1-vector defined by $({\bf
1}_{X})_{i}=1\Leftrightarrow i\in X$. The all one vector ${\bf 1}_{[n]}$ is
simply denoted by ${\bf 1}$.

For an $n\times m$ matrix $A=(A_{ij})\in{\mathbb{R}}^{n\times m}$,
$I\subseteq[n]$, and $J\subseteq[m]$, let $A[I,J]$ denote the submatrix of $A$
consisting of $A_{ij}$ for $i\in I,j\in J$. If $A$ is a block diagonal matrix
with diagonals $D_{1},D_{2},\ldots,D_{k}$, then $A$ is written as
$A=\bigoplus_{i\in[k]}D_{i}$.

For two nonnegative vectors $p,q\in{\mathbb{R}}^{n}$, let $D(p\|q)$ denote the
KL-divergence between $p$ and $q$:

$D(p\|q):=\sum_{i=1}^{n}p_{i}\log p_{i}/q_{i},$

where $0\log 0/x:=0$ for $x\geq 0$ and $x\log x/0:=\infty$ for $x>0$.

When $p([n])=q([n])$, the KL-divergence relates the $\ell_{1}$-distance via
Pinsker’s inequality:

$D(p\|q)\geq\frac{1}{2p([n])}\|p-q\|_{1}^{2}.$ (1.1)

##  2 Matrix scaling

In this section, we summarize basic results on matrix scaling particularly
from the KL-divergence interpretation. Let $A=(A_{ij})$ be an $n\times m$
nonnegative matrix. We assume that $A$ has neither zero rows nor no zero
columns. A scaling of $A$ is a matrix $\tilde{A}=(\tilde{A}_{ij})$ written as
$\tilde{A}_{ij}=A_{ij}x_{i}y_{j}$ for some positive vectors
$x\in{\mathbb{R}}^{n}$ and $y\in{\mathbb{R}}^{m}$.

Let $r$ and $c$ be $n$-dimensional and $m$-dimensional positive vectors,
respectively. An $(r,c)$-scaling of $A$ is a scaling $\tilde{A}$ of $A$ such
that $\tilde{A}{\bf 1}=r$ and $\tilde{A}^{\top}{\bf 1}=c$. The matrix scaling
problem asks an $(r,c)$-scaling of $A$. If an $(r,c)$-scaling of $A$ exists,
then $A$ is said to be $(r,c)$-scalable. If there is a scaling $\tilde{A}$ of
$A$ such that $(\tilde{A}{\bf 1},\tilde{A}^{\top}{\bf 1})$ is arbitrarily
close $(r,c)$ (under some norm of ${\mathbb{R}}^{n+m}$), then $A$ is
approximately $(r,c)$-scalable. Let $R:=r([n])$ and $C:=c([m])$. We do not
assume $R=C$, which is an obvious necessary condition for scalability. The
case of $(r,c)=({\bf 1},{\bf 1})$ corresponds to the doubly stochastic scaling
mentioned in introduction, where “nonscalable” is meant as “not approximately
$(r,c)$-scalable.”

We next explain combinatorial characterizations of scalability by Rothblum and
Schneider [21]. Let $G(A)=([n]\sqcup[m],E)$ denote the bipartite graph
representing the nonzero pattern of $A$, Specifically, vertices $i\in[n]$ and
$j\in[m]$ are joined by an edge $ij$ in $G(A)$ if and only if $A_{ij}>0$. The
notation $X\sqcup Y$ for a subset of vertices is meant as the disjoint union
of $X\subseteq[n]$ and $Y\subseteq[m]$. A subset $X\sqcup Y$ of vertices is
said to be stable if there is no edge $ij$ with $i\in X$ and $j\in Y$. Let
${\cal S}$ denote the set of all stable set.

######  Theorem 2.1 ([21]).

  * (1)

$A$ is approximately $(r,c)$-scalable if and only if $R=C$ and $r(X)+c(Y)\leq
C$ for $X\sqcup Y\in{\cal S}$.

  * (2)

$A$ is $(r,c)$-scalable if and only if $R=C$ and $r(X)+c(Y)\leq C$ for
$X\sqcup Y\in{\cal S}$, and the equality $r(X)+c(Y)=C$ for $X\sqcup Y\in{\cal
S}$ implies $A[[n]\setminus X,[m]\setminus Y]=O$.

  * (3)

An $(r,c)$-scaling of $A$ is unique if it exists.

The Sinkhorn algorithm is a simple method to obtain an approximate
$(r,c)$-scaling. The row-normalization $R(A)$ of $A$ is a scaling of $A$
defined by $R(A)_{ij}:=(r_{i}/p_{i})A_{ij}$ for $p_{i}:=(A{\bf 1})_{i}$.
Similarly, the column-normalization $C(A)$ of $A$ is a scaling of $A$ defined
by $C(A)_{ij}:=(c_{j}/q_{j})A_{ij}$ for $q_{j}:=(A^{\top}{\bf 1})_{j}$. The
Sinkhorn algorithm is to repeat row-normalization $A\leftarrow R(A)$ and
column-normalization $A\leftarrow C(A)$ alternatively.

We mainly concern the situation where $A$ is not approximately
$(r,c)$-scalable or the Sinkhorn algorithm is not convergent. So we define a
sequence $\\{(M_{k},N_{k})\\}_{k=0,1,2,\ldots}$ generated by the Sinkhorn
algorithm in the following way:

$\displaystyle N_{0}$ $\displaystyle:=$ $\displaystyle C(A),$ $\displaystyle
M_{k}$ $\displaystyle:=$ $\displaystyle R(N_{k}),$ $\displaystyle N_{k+1}$
$\displaystyle:=$ $\displaystyle C(M_{k}).$ (2.1)

Define the spaces ${\cal M},{\cal R},{\cal C}$ of $n\times m$ matrices by

$\displaystyle{\cal M}$ $\displaystyle:=$
$\displaystyle\\{M\in{\mathbb{R}}^{n\times m}\mid M_{ij}=0\
(i,j:A_{ij}=0),M_{ij}\geq 0\ (i,j:A_{ij}>0)\\},$ $\displaystyle{\cal R}$
$\displaystyle:=$ $\displaystyle\\{M\in{\cal M}\mid M{\bf 1}=r\\},$
$\displaystyle{\cal C}$ $\displaystyle:=$ $\displaystyle\\{M\in{\cal M}\mid
M^{\top}{\bf 1}=c\\}.$

Consider the convex optimization for the minimum KL-divergence between ${\cal
R}$ and ${\cal C}$:

${\rm inf}.\quad D(M\|N)=\sum_{i,j}M_{ij}\log M_{ij}/N_{ij}\quad{\rm
s.t.}\quad M\in{\cal R},N\in{\cal C}.$ (2.2)

Then the Sinkhorn algorithm is viewed as alternating minimization of this
problem:

######  Lemma 2.2 (see [13]).

The following hold:

  * (1)

$D(M\|C(M))=\inf_{N\in{\cal C}}D(M\|N)$ $(M\in{\cal M})$.

  * (2)

$D(R(N)\|N)=\inf_{M\in{\cal R}}D(M\|N)$ $(N\in{\cal M})$.

The convergence of the sequence $(M_{k},N_{k})$ can be proved by the following
properties.

######  Proposition 2.3 ([5]; see also [13]).

For $M\in{\cal R},N\in{\cal C}$, the following hold:

(3-point property)

$D(M_{k}\|N_{k})+D(M\|M_{k})=D(M\|N_{k})$.

(4-point property)

$D(M\|N_{k+1})\leq D(M\|M_{k})+D(M\|N)$.

(5-point property)

$D(M_{k}\|N_{k})+D(M\|N_{k+1})\leq D(M\|N)+D(M\|N_{k})$.

For completeness, the proofs of Lemma 2.2 and Proposition 2.3 are given in
Appendix.

Choose an arbitrary minimizer $(M^{*},N^{*})$ of (2.2). It actually exists,
since $D$ is a lower semi-continuous convex function on compact set ${\cal
R}\times{\cal C}$. From the 5-point property, we have

$D(M_{k}\|N_{k})-D(M^{*}\|N^{*})\leq
D(M^{*}\|N_{k})-D(M^{*}\|N_{k+1})\quad(k=0,1,2,\ldots).$

Adding them from $k=0$ to $\ell$, we obtain

$\sum_{k=1}^{\ell}D(M_{k}\|N_{k})-D(M^{*}\|N^{*})\leq
D(M^{*}\|N_{0})-D(M^{*}\|N_{\ell+1}).$

Then it holds $\ell(D(M_{\ell}\|N_{\ell})-D(M^{*}\|N^{*}))\leq$ LHS since

$D(M_{k}\|N_{k})\geq D(M_{k}\|N_{k+1})\geq D(M_{k+1}\|N_{k+1})\geq
D(M^{*}\|N^{*}).$

Hence we obtain:

######  Lemma 2.4.

For a positive integer $\ell>0$, it holds

$D(M_{\ell}\|N_{\ell})-D(M^{*}\|N^{*})\leq\frac{D(M^{*}\|N_{0})}{\ell}.$ (2.3)

Thus the sequence $D(M_{k}\|N_{k})$ converges to the infimum of (2.2).
Particularly, $\\{(M_{k},N_{k})\\}_{k}$ has an accumulating point
$(M^{*},N^{*})$ in compact set ${\cal R}\times{\cal C}$, which is also a
minimizer of (2.2). By adding 3-point property for $k+1$ and 4-point property
for $k$, we obtain $D(M^{*}\|M_{k+1})\leq D(M^{*}\|M_{k})$. Necessarily
$D(M^{*}\|M_{k})$ converges to zero. By Pinsker’s inequality (1.1), we have:

######  Theorem 2.5 ([5]; see [13]).

The sequence $(M_{k},N_{k})$ converges to a minimizer $(M^{*},N^{*})$ of
(2.2).

We call $(M^{*},N^{*})$ the Sinkhorn limit of $A$. If ${\cal R}\cap{\cal
C}\neq\emptyset$, then $D(M^{*}\|N^{*})=0$ and $M^{*}=N^{*}$; this is
precisely the case where $A$ is approximately $(r,c)$-scalable.

######  Theorem 2.6 ([3, 23]).

Suppose that $A$ is approximately $(r,c)$-scalable. $M_{k}$ and $N_{k}$
converge to an $(r,c)$-scaling $M^{*}=N^{*}$.

It is known [4] that the Sinkhorn limit $M^{*}=N^{*}$ in this case is the
unique minimizer of $D(M\|A)$ over $M\in{\cal R}\cap{\cal C}$.

Note that the approximate scalability is equivalent to ${\cal R}\cap{\cal
C}\neq\emptyset$. This is a network flow feasibility condition. Then Theorem
2.1 (1) can be deduced from max-flow min-cut theorem; see the beginning of the
next section.

######  Remark 2.7.

Alternating minimization of this kind is studied in the literature of first
order methods—a current active area of optimization; see [2, Chapter 14].
Since KL-divergence does not satisfy the $L$-smoothness assumption, we could
not find from the literature results giving better estimates.

##  3 Polymatroidal analysis on the Sinkhorn limit

In this section, we study the Sinkhorn limit $(M^{*},N^{*})$ in more detail,
particularly from polymatroidal viewpoint. Consider the set $P:=P(A,c)$ of all
vectors $p\in{\mathbb{R}}^{n}$ for which $A$ is approximately
$(p,c)$-scalable. Then it holds

$\displaystyle P$ $\displaystyle=$ $\displaystyle\\{p\in{\mathbb{R}}^{n}\mid
p(X)+c(Y)\leq C\ (X\sqcup Y\in{\cal S}),p([n])=C\\}$ (3.1) $\displaystyle=$
$\displaystyle\\{p\in{\mathbb{R}}^{n}\mid p=N{\bf 1},N\in{\cal C}\\}.$

The first equality follows from Theorem 2.1 (1). One can see the second
equality from a standard argument of network flow: Regard each edge $ij$ in
$G$ a directed edge from $i$ to $j$ having infinite capacity. Add source $s$
and sink $t$ together with directed edges $si$ $(i\in[n])$ and $jt$
$(j\in[m])$. The capacity of $si$ is defined as $p_{i}$ and the capacity of
$jt$ is defined as $c_{i}$. Let $\vec{G}(A,p,c)$ denote the resulting network.
From $N$ in (3.1). we obtain an $(s,t)$-flow in $\vec{G}(A,p,c)$ such that the
flow-values of edges $si,ij,jt$ are $p_{i},N_{ij},c_{j}$, respectively. In
this way, $N$ can be identified with a maximum flow of the flow-value
$C=p([n])$. Also an stable set $X\sqcup Y\in{\cal S}$ is identified with an
$(s,t)$-cut $\\{s\\}\cup X\cup([m]\setminus Y)$ with finite capacity
$p([n])-p(X)+C-c(Y)$. Then, by the max-flow min-cut theorem, the condition
$p(X)+c(Y)\leq C$ $(X\sqcup Y\in{\cal S})$, $p([n])=C$ is necessarily and
sufficient for $p$ to be represented as $p=N{\bf 1}$.

######  Lemma 3.1.

  * (1)

$D(R(N)\|N)=D(r\|N{\bf 1})$ for $N\in{\cal M}$.

  * (2)

$\inf_{M\in{\cal R},N\in{\cal C}}D(M\|N)=\inf_{p\in P}D(r\|p)$.

###### Proof.

(1). Substitute $R(N)_{ij}=(r_{i}/N{\bf 1})_{i}N_{ij}$ and
$r_{i}=\sum_{j}R(N)_{ij}$ for $D(R(N)\|N)=\sum_{i,j}R(N)_{ij}\log
R(N)_{ij}/N_{ij}$.

(2). We have $\leq$ from (1). Conversely, choose $(M^{*},N^{*})\in{\cal
R}\times{\cal C}$ attaining the infimum. By Lemma 2.2 (2) and the strict
convexity of $M\mapsto D(M\|N^{*})$, it must hold $M^{*}=R(N^{*})$. Then
$D(M^{*}\|N^{*})=D(R(N^{*})\|N^{*})=D(r\|p^{*})$, where $p^{*}=N^{*}{\bf 1}$.
∎

Therefore, the value $D(M^{*}\|N^{*})$ is also given by optimization problem:

${\rm inf}.\quad D(r\|p)\quad{\rm s.t.}\quad p\in P.$ (3.2)

In particular, by $D(M_{k}\|N_{k})=D(r\|N_{k}{\bf 1})$, the sequence
$N_{k}{\bf 1}$ converges to an optimal solution $p^{*}$ of the problem.

We are going to give an explicit formula of $p^{*}$. Let $\Gamma:2^{[n]}\to
2^{[m]}$ be defined by

$\Gamma(X):=\\{j\in[m]\mid\mbox{$j$ is adjacent to a node in $X$ in
$G(A)$}\\}.$ (3.3)

Since every stable set $X\sqcup Y$ is contained in a larger stable set
$X\sqcup([m]\setminus\Gamma(X))$. the inequality system for $P$ can be written
as

$p(X)\leq c(\Gamma(X))\ (X\subseteq[n]),\quad p([n])=c(\Gamma([n]))(=C).$
(3.4)

Observe that $X\mapsto c(\Gamma(X))$ is a monotone submodular function.
Therefore, $P$ is the base polytope of the polymatroid; see [8].

Consider a map ${\cal S}\to{\mathbb{R}}^{2}$ by

$X\sqcup Y\mapsto(r(X),c(Y))\quad(X\sqcup Y\in{\cal S}),$ (3.5)

and consider the convex hull, denoted by $\mathop{\rm Conv}{\cal S}$, of
points $(r(X),c(Y))$ over all stable sets $X\sqcup Y\in{\cal S}$. The convex
polygon $\mathop{\rm Conv}{\cal S}$ is contained by box $[R,0]\times[0,C]$ and
contains $(0,0)$, $(R,0)$, $(0,C)$ as extreme points. A stable set $X\sqcup Y$
is called extreme if $(r(X),c(Y))$ is a nonzero extreme point of $\mathop{\rm
Conv}{\cal S}$. See the left of Figure 1.

Figure 1: $\mathop{\rm Conv}{\cal S}$ in $(y,x)$-plane (left) and the
associated block-triangularization of $A$ (right)

######  Lemma 3.2.

For extreme stable sets $X\sqcup Y$, $X^{\prime}\sqcup Y^{\prime}$, if
$r(X)\geq r(X^{\prime})$ and $c(Y)\leq c(Y^{\prime})$, then $X\supseteq
X^{\prime}$ and $Y\subseteq Y^{\prime}$. In particular, the map (3.5) is
injective on the family of extreme stable sets.

###### Proof.

It suffices to consider the case where points $(r(X),c(Y))$ and
$(r(X^{\prime}),c(Y^{\prime}))$ are equal or adjacent extreme points. Observe
that $(X\cap X^{\prime})\sqcup(Y\cup Y^{\prime})$ and $(X\cup
X^{\prime})\sqcup(Y\cap Y^{\prime})$ are both stable, and that

$(r(X),c(Y))+(r(X^{\prime}),c(Y^{\prime}))=(r(X\cap X^{\prime}),c(Y\cup
Y^{\prime}))+(r(X\cup X^{\prime}),c(Y\cap Y^{\prime})).$

This means that the point of the RHS divided by $2$ is equal to the midpoint
of the edge between $(r(X),c(Y))$ and $(r(X^{\prime}),c(Y^{\prime}))$.
Necessarily $(r(X),c(Y))=(r(X\cup X^{\prime}),c(Y\cap Y^{\prime}))$ must hold
otherwise $(r(X\cap X^{\prime}),c(Y\cup Y^{\prime}))$ or $(r(X\cup
X^{\prime}),c(Y\cap Y^{\prime}))$ goes outside of $\mathop{\rm Conv}{\cal S}$.
Since $r$ and $c$ are positive vectors, it must hold $X=X\cup X^{\prime}$ and
$Y=Y\cap Y^{\prime}$, as required. ∎

Consider all extreme stable sets $X_{\kappa}\sqcup Y_{\kappa}$
$(\kappa=0,1,2,\ldots,\theta)$. Since $c$ is a positive vector, it necessarily
holds

$Y_{\kappa}=[m]\setminus\Gamma(X_{\kappa}).$ (3.6)

By Lemma 3.2, they can be ordered as

$\displaystyle[n]=X_{0}\supset X_{1}\supset\cdots\supset
X_{\theta}=\emptyset,$ $\displaystyle\emptyset=Y_{0}\subset
Y_{1}\subset\cdots\subset Y_{\theta}=[m].$

Let $(I_{1},I_{2},\ldots,I_{\theta})$ and $(J_{1},J_{2},\ldots,J_{\theta})$ be
partitions of $[n]$ and $[m]$, respectively, defined by

$I_{\kappa}:=X_{\kappa-1}/X_{\kappa},\quad
J_{\kappa}:=Y_{\kappa}/Y_{\kappa-1}\quad(\kappa=1,2,\ldots,\theta).$ (3.7)

Notice that $A[X,Y]=O$ if and only if $X\sqcup Y\in{\cal S}$. We can arrange
$A$ so that these submatrices $A[I_{\kappa},J_{\kappa}]$ are placed in
diagonal positions and their lower left blocks are all zero matrices. See the
right of Figure 1.

Let $r^{\kappa}$ denote the $|I_{\kappa}|$-dimensional vector obtained by
restricting $r$ to the indices in $I_{\kappa}$. Similarly, let $c^{\kappa}$
denote the restriction of $c$ to the index set $J_{\kappa}$. Let
$R_{\kappa}:=r(I_{\kappa})$ and $C_{\kappa}:=c(J_{\kappa})$ denote the total
sums of $r^{\kappa}$ and $c^{\kappa}$, respectively.

Define $p^{*}\in{\mathbb{R}}^{n}$ by $p^{*}_{i}:=(C_{\kappa}/R_{\kappa})r_{i}$
for $i\in I_{\kappa}$. By an appropriate ordering, $p^{*}$ is written as

$p^{*}=\left(\begin{array}[]{c}(C_{1}/R_{1})r^{1}\\\ (C_{2}/R_{2})r^{2}\\\
\vdots\\\ (C_{\theta}/R_{\theta})r^{\theta}\end{array}\right).$ (3.8)

######  Lemma 3.3.

  * (1)

Each submatrix $A[I_{\kappa},J_{\kappa}]$ is approximately
$(r^{\kappa}/R_{\kappa},c^{\kappa}/C_{\kappa})$-scalable.

  * (2)

It holds $R_{1}/C_{1}<R_{2}/C_{2}<\cdots<R_{\theta}/C_{\theta}$.

###### Proof.

(1). Suppose otherwise. Then there is a stable set $X\sqcup Y\subseteq
I_{\kappa}\sqcup J_{\kappa}$ such that
$(1/R_{\kappa})r^{\kappa}(X)+(1/C_{\kappa})c^{\kappa}(Y)>1$. Notice that the
equation of the line through extreme points
$(r(X_{\kappa-1}),c(Y_{\kappa-1}))$ and $(r(X_{\kappa}),c(Y_{\kappa}))$ in the
$(x,y)$-plane ${\mathbb{R}}^{2}$ is given by
$(1/R_{\kappa})(x-r(X_{\kappa}))+(1/C_{\kappa})(y-c(Y_{\kappa-1}))=1$. Then
$(X\cup X_{\kappa})\sqcup(Y\cup Y_{\kappa-1})$ is a stable set in the original
graph $G$ such that the point $(r(X\cup X_{\kappa}),c(Y\cup Y_{\kappa-1}))$ is
outside of $\mathop{\rm Conv}{\cal S}$. This is a contradiction.

(2). It is clear from the observation that $R_{\kappa}/C_{\kappa}$ is a slope
of the edge between extreme points $(r(X_{\kappa-1}),c(Y_{\kappa-1}))$ and
$(r(X_{\kappa}),c(Y_{\kappa}))$ in the convex polygon $\mathop{\rm Conv}{\cal
S}$. ∎

######  Proposition 3.4.

$p^{*}$ is the unique minimizer of (3.2).

###### Proof.

We first verify that $p^{*}$ is a point in $P$. By Lemma 3.3 (1) and (3.1),
for each $\kappa$, there is $|I_{\kappa}|\times|J_{\kappa}|$ matrix
$B_{\kappa}$ such that $(B_{\kappa})_{ij}=0$ for each $i,j$ with
$(A[I_{\kappa},J_{\kappa}])_{ij}=0$, $B_{\kappa}{\bf
1}=r^{\kappa}/R_{\kappa}$, and $B_{\kappa}^{\top}{\bf
1}=c^{\kappa}/C_{\kappa}$. Let $B=\bigoplus_{\kappa}C_{\kappa}B_{\kappa}$ be
the block diagonal matrix with diagonals $C_{\kappa}B_{\kappa}$. Then
$B\in{\cal C}$, and $p^{*}=B{\bf 1}$.

We next verify that $p^{*}$ satisfies the KKT-condition for (3.2). The
gradient $\nabla D(r\|p)=(\partial D(r\|p)/\partial p_{i})_{i}$ at $p^{*}$ is
given by

$\displaystyle\nabla D(r\|p^{*})=-\sum_{i}\frac{r_{i}}{p^{*}_{i}}{\bf
1}_{\\{i\\}}=-\sum_{\kappa=1}^{\theta}\frac{R_{\kappa}}{C_{\kappa}}({\bf
1}_{X_{\kappa-1}}-{\bf 1}_{X_{\kappa}})$
$\displaystyle=-\frac{R_{1}}{C_{1}}{\bf
1}-\sum_{\kappa=1}^{\theta-1}\left(\frac{R_{\kappa+1}}{C_{\kappa+1}}-\frac{R_{\kappa}}{C_{\kappa}}\right){\bf
1}_{X_{\kappa}}.$ (3.9)

Also $p^{*}$ satisfies the following inequalities from (3.4) in equality:

$({\bf
1}_{X_{\kappa}})^{\top}p^{*}=\sum_{\lambda=\kappa+1}^{\theta}p^{*}(I_{\lambda})=\sum_{\lambda=\kappa+1}^{\theta}(C_{\lambda}/R_{\lambda})R_{\lambda}=\sum_{\lambda=\kappa+1}^{\theta}c(J_{\lambda})=c(\Gamma(X_{\kappa})),$

where we use $[m]\setminus Y_{\kappa}=\Gamma(X_{\kappa})$ in the last
equality. Also $P$ has the equality constraint ${\bf 1}^{\top}p=C$. By Lemma
3.3 (2), it holds $R_{\kappa+1}/C_{\kappa+1}-R_{\kappa}/C_{\kappa}>0$. Thus,
(3.9) says the KKT-condition that $-\nabla D(r\|p^{*})$ is a nonnegative
combination of ${\bf 1}_{X_{\kappa}}$. By strict convexity of $-\log$, it is a
unique minimizer. ∎

Thus we have the following, which is implicit in [1].

######  Theorem 3.5.

$N^{*}{\bf 1}=p^{*}$.

In particular, it holds ${\bf 1}_{I_{\kappa}}^{\top}N^{*}{\bf
1}_{Y_{\kappa-1}}=p^{*}(I_{\kappa})=C_{\kappa}={\bf 1}_{[n]\setminus
X_{\kappa}}^{\top}N^{*}{\bf 1}_{J_{\kappa}}$. This means that all off-diagonal
blocks of $N^{*}$ are zero, i.e.,
$N^{*}=\bigoplus_{\kappa}N^{*}[I_{\kappa},J_{\kappa}]$, where the diagonal
block $N^{*}[I_{\kappa},J_{\kappa}]$ has marginal
$((C_{\kappa}/R_{\kappa})r^{\kappa},c^{\kappa})$. Accordingly,
$M^{*}=R(N^{*})$ is also a diagonal matrix
$M^{*}=\bigoplus_{\kappa}M^{*}[I_{\kappa},J_{\kappa}]$ such that
$M^{*}[I_{\kappa},J_{\kappa}]$ has marginal
$(r^{\kappa},(R_{\kappa}/C_{\kappa})c^{\kappa})$. A sharper version is given
in Theorem 3.8 below.

Next we give an estimate that relates $D(r\|p)-D(r\|p^{*})$, the distance
between $p*$ and $p=N{\bf 1}$, and the magnitude of the off-diagonal blocks of
$N$. Suppose that $p=N{\bf 1}$ for $N\in{\cal C}$. For
$1\leq\kappa<\lambda\leq\theta$, define $\Delta_{\kappa\lambda}$ by

$\Delta_{\kappa\lambda}:=\sum_{i\in I_{\kappa},j\in J_{\lambda}}N_{ij}.$
(3.10)

This is the total sum of entries in the off-diagonal block
$N[I_{\kappa},J_{\lambda}]$. Also define $\Delta$ and $\Delta_{\kappa}$ by

$\displaystyle\Delta$ $\displaystyle:=$
$\displaystyle\sum_{1\leq\kappa<\lambda<\theta}\Delta_{\kappa\lambda},$
$\displaystyle\Delta_{\kappa}$ $\displaystyle:=$
$\displaystyle\sum_{\gamma:\kappa<\gamma}\Delta_{\kappa\gamma}-\sum_{\gamma:\gamma<\kappa}\Delta_{\gamma\kappa}.$

Observe the relation:

$p(I_{\kappa})=C_{\kappa}+\Delta_{\kappa}=p^{*}(I_{\kappa})+\Delta_{\kappa}.$
(3.11)

######  Lemma 3.6.

$D(r\|p)-D(r\|p^{*})\geq\sum_{\kappa=1}^{\theta}\frac{R_{\kappa}}{2(C_{\kappa})^{2}}\left\|p^{*\kappa}-\frac{C_{\kappa}}{C_{\kappa}+\Delta_{\kappa}}p^{\kappa}\right\|_{1}^{2}+\sum_{1\leq\kappa<\lambda\leq\theta}\left(\frac{R_{\lambda}}{C_{\lambda}}-\frac{R_{\kappa}}{C_{\kappa}}\right)\Delta_{\kappa\lambda}.$

###### Proof.

Let $\tilde{p}^{\kappa}:=p^{\kappa}C_{\kappa}/(C_{\kappa}+\Delta_{\kappa})$.
The LHS is given by

$\displaystyle\sum_{i}r_{i}\log\frac{r_{i}}{p_{i}}-\sum_{i}r_{i}\log\frac{r_{i}}{p^{*}_{i}}=\sum_{i}r_{i}\log\frac{p^{*}_{i}}{p_{i}}=\sum_{\kappa=1}^{\theta}\frac{R_{\kappa}}{C_{\kappa}}\sum_{i\in
I_{\kappa}}p^{*}_{i}\log\frac{p^{*}_{i}}{p_{i}}$
$\displaystyle=\sum_{\kappa=1}^{\theta}\frac{R_{\kappa}}{C_{\kappa}}\left\\{D(p^{*\kappa}\|\tilde{p}^{\kappa})-\sum_{i\in
I_{\kappa}}p^{*}_{i}\log\left(1+\frac{\Delta_{\kappa}}{C_{\kappa}}\right)\right\\}$
$\displaystyle\geq\sum_{\kappa=1}^{\theta}\frac{R_{\kappa}}{C_{\kappa}}\left\\{\frac{\|p^{*\kappa}-\tilde{p}^{\kappa}\|_{1}^{2}}{2C_{\kappa}}-C^{\kappa}\frac{\Delta_{\kappa}}{C_{\kappa}}\right\\}$
$\displaystyle=\sum_{\kappa=1}^{\theta}\frac{R_{\kappa}}{2(C_{\kappa})^{2}}\|p^{*\kappa}-\tilde{p}^{\kappa}\|_{1}^{2}-\sum_{\kappa=1}^{\theta}\frac{R_{\kappa}}{C_{\kappa}}\Delta_{\kappa}$
$\displaystyle=\sum_{\kappa=1}^{\theta}\frac{R_{\kappa}}{2(C_{\kappa})^{2}}\|p^{*\kappa}-\tilde{p}^{\kappa}\|_{1}^{2}-\sum_{1\leq\kappa<\lambda\leq\theta}\left(\frac{R_{\kappa}}{C_{\kappa}}-\frac{R_{\lambda}}{C_{\lambda}}\right)\Delta_{\kappa\lambda},$

where the inequality follows from Pinsker’s inequality (1.1) and
$-\log(1+x)\geq-x$. ∎

Although this estimate is enough for our goal (Theorem 4.1), we continue to
study the Sinkhorn limit $(M^{*},N^{*})$ to give a final form. Extreme stable
sets are obtained via the following parametric optimization for parameter
$\lambda\in[0,1]$:

${\rm max}.\quad(1-\lambda)r(X)+\lambda c(Y)\quad{\rm s.t.}\quad X\sqcup
Y\in{\cal S}.$ (3.12)

This is a minimum $(s,t)$-cut problem on the network
$\vec{G}(A,(1-\lambda)r,\lambda c)$. Indeed, $X\sqcup Y$ corresponds to
$(s,t)$-cut $\\{s\\}\cup X\cup([m]\setminus Y)$ with capacity
$-(1-\lambda)r(X)-\lambda c(Y)+(1-\lambda)R+\lambda C$. Particularly, (3.12)
for all $\lambda$ can be efficiently solved by a parametric maximum flow
algorithm, e.g., [10].

An extreme stable set is precisely a unique maximizer of (3.12) for some
$\lambda\in[0,1]$. A parameter $\lambda$ is said to be critical if a minimizer
is not unique. The number of critical parameters, i.e., the number of slopes
in $\mathop{\rm Conv}{\cal S}$, equals $\theta$. Let
$\lambda_{1},\lambda_{2},\ldots\lambda_{\theta}$ be critical parameters
ordered as $0<\lambda_{1}<\lambda_{2}<\cdots<\lambda_{\theta}<1$. Let ${\cal
S}_{\kappa}$ denote the family of maximizers of (3.12) for
$\lambda=\lambda_{\kappa}$. Then it is not difficult to see

  * •

$X\sqcup Y,X^{\prime}\sqcup Y^{\prime}\in{\cal
S}_{\kappa}\Longrightarrow(X\cup X^{\prime})\sqcup(Y\cap Y^{\prime}),(X\cap
X^{\prime})\sqcup(Y\cup Y^{\prime})\in{\cal S}_{\kappa}$.

  * •

$X_{\kappa-1}\sqcup Y_{\kappa-1},X_{\kappa}\sqcup Y_{\kappa}\in{\cal
S}_{\kappa}$.

  * •

$X\sqcup Y\in{\cal S}_{\kappa}\Longrightarrow X_{\kappa-1}\supseteq X\supseteq
X_{\kappa},Y_{\kappa}\subseteq Y\subseteq Y_{\kappa-1}$.

Namely, ${\cal S}_{\kappa}$ admits a distributive lattice structure, where
$(X\sqcup Y)\wedge(X^{\prime}\sqcup Y^{\prime}):=(X\cup
X^{\prime})\sqcup(Y\cap Y^{\prime})$ and $(X\sqcup Y)\vee(X^{\prime}\sqcup
Y^{\prime}):=(X\cap X^{\prime})\sqcup(Y\cup Y^{\prime})$. This is essentially
a well-known fact that the family of minimum $(s,t)$-cuts forms a distributive
lattice. Choose a maximal set (chain) $\\{X_{\kappa,\alpha}\sqcup
Y_{\kappa,\alpha}\\}_{\alpha=0,1,2,\ldots,\alpha_{\kappa}}\subseteq{\cal
S}_{\kappa}$ with the property:

$\displaystyle X_{\kappa}=X_{\kappa,0}\supset X_{\kappa,1}\supset\cdots\supset
X_{\kappa,\alpha_{\kappa}}=X_{\kappa+1},$ $\displaystyle
Y_{\kappa}=Y_{\kappa,0}\subset Y_{\kappa,1}\subset\cdots\subset
Y_{\kappa,\alpha_{\kappa}}=Y_{\kappa+1}.$

This is also efficiently obtained from the residual network of
$\vec{G}(A,(1-\lambda)r,\lambda c)$ with respect to any maximum flow.
Accordingly, define $I_{\kappa,\alpha}$ and $J_{\kappa,\beta}$ by

$I_{\kappa,\alpha}:=X_{\kappa,\alpha-1}/X_{\kappa,\alpha},\quad
J_{\kappa,\alpha}:=Y_{\kappa,\alpha}/Y_{\kappa,\alpha-1}\quad(\alpha=1,2,\ldots,\alpha_{\kappa}).$
(3.13)

Now $\\{I_{\kappa,\alpha}\\}_{\kappa,\alpha}$
and$\\{J_{\kappa,\alpha}\\}_{\kappa,\alpha}$ are partitions of $[n]$ and $[m]$
refining $\\{I_{\kappa}\\}_{\kappa}$ and$\\{J_{\kappa}\\}_{\kappa}$,
respectively. Again, we arrange $A$ so that
$A[I_{\kappa,\alpha},J_{\kappa,\alpha}]$ are diagonal blocks and their lower
left blocks are zero. This is a finer block-triangularization of $A$.

The point $(r(X_{\kappa,\alpha}),c(Y_{\kappa,\alpha}))$ lies on the segment
between extreme points $(R_{\kappa},C_{\kappa})$ and
$(R_{\kappa+1},C_{\kappa+1})$ in $\mathop{\rm Conv}{\cal S}$. By the
maximality of the chain, there is no stable set $X\sqcup Y$ in
$I_{\kappa,\alpha}\sqcup J_{\kappa,\alpha}$ such that the stable set
$(X_{\kappa,\alpha+1}\cup X)\sqcup(Y_{\kappa,\alpha}\cup Y)$ is mapped to the
segment between $(r(X_{\kappa,\alpha}),c(Y_{\kappa,\alpha}))$ and
$(r(X_{\kappa,\alpha+1}),c(Y_{\kappa,\alpha+1}))$. This implies the condition
of Theorem 2.1 (2):

######  Lemma 3.7.

Each submatrix $A[I_{\kappa,\alpha},J_{\kappa,\alpha}]$ is
$(r^{\kappa,\alpha}/R_{\kappa},c^{\kappa,\alpha}/C_{\kappa})$-scalable.

Notice that $C_{\kappa}/R_{\kappa}=c(J_{\kappa,\alpha})/r(I_{\kappa,\alpha})$
for each $\alpha$. Again, all off-diagonal blocks in the limit
$M^{*}[I_{\kappa},J_{\kappa}]$ must be zero; see the argument after Theorem
3.5. Thus we have the following, which is implicit in [1].

######  Theorem 3.8.

The Sinkhorn limit $(M^{*},N^{*})$ is given by

$M^{*}=\bigoplus_{\kappa\in[\theta],\alpha\in[\alpha_{\alpha}]}R_{\kappa}B_{\kappa,\alpha},\quad
N^{*}=\bigoplus_{\kappa\in[\theta],\alpha\in[\alpha_{\alpha}]}C_{\kappa}B_{\kappa,\alpha},$
(3.14)

where $B_{\kappa,\alpha}$ denotes the (unique)
$(r^{\kappa,\alpha}/R_{\kappa},c^{\kappa,\alpha}/C_{\kappa})$-scaling of
$A[I_{\kappa,\alpha},J_{\kappa,\alpha}]$.

### Relation to principal partition and DM-decomposition.

The problem (3.12) is a parametric submodular optimization, since it is also
written as minimization of $X\mapsto tc(\Gamma(X))-(1-t)r(X)$. Such a
parametric problem has been studied in a general framework of the principal
partition of polymatroids [9, 24]. The above nested distributive lattice
structure and the associated partition are generalized to this setting, which
has a number of applications in combinatorial optimization; see the above
surveys.

In the case of $(r,c)=({\bf 1},{\bf 1})$, the block-triangularization of $A$
obtained from the chain $X_{\kappa,\alpha}\sqcup Y_{\kappa,\alpha}$ is a
refinement of Dulmage-Mendelsohn decomposition (DM-decomposition) [7]; see
[18, Section 4.3] and [20, Section 2.2.3]. The DM-decomposition uses a chain
of maximizers of (3.12) only for $\lambda=1/2$. This refined DM-decomposition
considering all parameters is due to N. Tomizawa (unpublished 1977); see [24].

##  4 Finding Hall blockers by Sinkhorn iteration

Let $G=(U\sqcup V,E)$ be a bipartite graph with color classes $U,V$ and edge
set $E$. We assume that there is no isolated node. Let $n_{1}:=|U|$,
$n_{2}:=|V|$, and $n:=n_{1}+n_{2}$. Define $n_{1}\times n_{2}$ matrix $A_{G}$
by $(A_{G})_{ij}=1$ if $ij\in E$ and zero otherwise, where rows and columns of
$A_{G}$ are indexed by $U$ and $V$, respectively. The following is well-known:

(Hall’s theorem)

$G$ has a perfect matching if and only if $|X|\leq\Gamma(X)$ for all
$X\subseteq[n]$.

(König-Egervary theorem)

The maximum cardinality of a matching in $G$ is equal to the minimum of
$|U|-|X|+\Gamma(X)$ over $X\subseteq U$.

Letting $(r,c):=({\bf 1},{\bf
1})\in{\mathbb{R}}^{n_{1}}\times{\mathbb{R}}^{n_{2}}$, consider a doubly
stochastic scaling for $A_{G}$, and consider the following algorithm for a
positive integer $\ell>0$:

Algorithm: Sinkhorn & Sorting

0:

$A\leftarrow C(A_{G})$.

1:

Repeat the row- and column-normalization $\ell$ times:

$A\leftarrow C(R(A))$.

2:

Let $p:=A{\bf 1}$, and sort $p$ as $p_{i_{1}}\leq p_{i_{2}}\leq\cdots\leq
p_{i_{n_{1}}}$.

3:

Output $U_{k}:=\\{i_{1},i_{2},\ldots,i_{k}\\}$ for $k=0,1,2,\ldots,n_{1}$,
where $U_{0}:=\emptyset$.

The goal of this section is to prove the following:

######  Theorem 4.1.

For $\ell\geq 64n^{7}\log n$, the output contains a subset $X\subseteq U$ with
maximum $|X|-\Gamma(X)$.

In particular, if $G$ has no perfect matching, then above $X$ is a Hall
blocker. Also the maximum matching number is obtained from $X$.

We start the proof. According to (2.1), define the sequence $(M_{k},N_{k})$
with the initial point $N_{0}:=C(A_{G})$. Consider the family
$\\{X_{\kappa}\sqcup Y_{\kappa}\\}_{\kappa=0,1,\ldots,\theta}$ of extreme
stable sets and the associated partitions $(I_{1},I_{2},\ldots,I_{\theta})$
and $(J_{1},J_{2},\ldots,J_{\theta})$ in (3.7). Then the limit
$p^{*}=N^{*}{\bf 1}$ is given by

$p^{*\kappa}=\frac{|J_{\kappa}|}{|I_{\kappa}|}{\bf
1}\quad(\kappa=1,2,\ldots,\theta).$

Observe from
$|J_{1}|/|I_{1}|>|J_{2}|/|I_{2}|>\cdots>|J_{\theta}|/|I_{\theta}|$ that
sorting $p^{*}$ recovers all $X_{\kappa}$.

In the step 2 of the algorithm, $A$ equals $N_{\ell}$. Let $p:=A{\bf 1}$.

######  Lemma 4.2.

If $\|p^{*}-p\|_{\infty}<1/(2n^{2})$, then the output contains $X_{\kappa}$
for all $\kappa=0,1,2,\ldots,\theta$.

###### Proof.

If $p_{i}>p_{j}$ holds for all $(i,j)\in I_{\kappa}\times I_{\lambda}$ with
$\kappa<\lambda$, then the output contains all $X_{\kappa}$. Since
$|J_{\kappa}|/|I_{\kappa}|-|J_{\lambda}|/|I_{\lambda}|\geq 1/n^{2}$ for
$\kappa<\lambda$, for $i\in I_{\kappa},j\in I_{\lambda}$, we have
$p_{i}-p_{j}>p_{i}^{*}-1/(2n^{*})-(p_{j}^{*}+1/(2n^{*}))>1/n^{2}-1/n^{2}=0$,
as required. ∎

######  Lemma 4.3.

$D(M^{*}\|N_{0})\leq n\log n$.

###### Proof.

Note that $(N_{0})_{ij}$ for $ij\in E$ equals $1/\deg(j)$, where $\deg(j)$ is
the degree of node $j\in V$. From $M^{*}_{ij}\leq 1$ and
$\sum_{i,j}M^{*}_{ij}=n_{1}$, we have $D(M^{*}\|N_{0})=\sum_{ij\in
E}M^{*}_{ij}\log M^{*}_{ij}+\sum_{ij\in E}M^{*}_{ij}\log\deg(j)\leq n\log n$.
∎

###### Proof of Theorem 4.1.

From $Y_{\kappa}=V\setminus\Gamma(X_{\kappa})$, any extreme stable set
$X_{\kappa}\sqcup Y_{\kappa}$ optimal to (3.12) with $\lambda=1/2$ gives a
maximizer $X_{\kappa}$ of $|X|-\Gamma(X)$. From Lemma 4.2, it suffices to show
that $\|p^{*}-p\|_{\infty}<1/(2n^{2})$ holds in step 2.

Let $\ell\geq Ln^{7}\log n$, where $L$ is a positive integer. By Lemmas 2.4,
3.6, and 4.3, we have

$\displaystyle\frac{1}{Ln^{6}}\geq\frac{D(M^{*}\|N_{0})}{\ell}\geq
D(R(A)\|A)-D(M^{*}\|N^{*})=D({\bf 1}\|p)-D({\bf 1}\|p^{*})$
$\displaystyle\geq\sum_{\kappa=1}^{\theta}\frac{|I_{\kappa}|}{2|J_{\kappa}|^{2}}\left\|p^{*\kappa}-\frac{|J_{\kappa}|}{|J_{\kappa}|+\Delta_{\kappa}}p^{\kappa}\right\|_{1}^{2}+\sum_{1\leq\kappa<\lambda\leq\theta}\left(\frac{|I_{\lambda}|}{|J_{\lambda}|}-\frac{|I_{\kappa}|}{|J_{\kappa}|}\right)\Delta_{\kappa\lambda},$
$\displaystyle\geq\sum_{\kappa=1}^{\theta}\frac{1}{2n^{2}(1+\Delta_{\kappa}/|J_{\kappa}|)^{2}}\left\|p^{*\kappa}-p^{\kappa}+{\bf
1}\Delta_{\kappa}/|I_{\kappa}|\right\|_{1}^{2}+\frac{1}{n^{2}}\Delta$

where we use $p^{*\kappa}={\bf 1}|J_{\kappa}|/|I_{\kappa}|$. Then
$\Delta<1/Ln^{4}$. Also we have

$\|p^{*\kappa}-p^{\kappa}+{\bf
1}\Delta_{\kappa}/|I_{\kappa}|\|^{2}_{1}<\frac{2n^{2}(1+\Delta_{\kappa}/|J_{\kappa}|)^{2}}{Ln^{6}}<\frac{8}{Ln^{4}},$

where we use
$|\Delta_{\kappa}|/|J_{\kappa}|\leq|\Delta_{\kappa}|\leq\Delta<1/Ln^{4}\leq
1$. From
$\|p^{*\kappa}-p^{\kappa}\|_{1}-|\Delta_{\kappa}|\leq\|p^{*\kappa}-p^{\kappa}+{\bf
1}\Delta_{\kappa}/|I_{\kappa}|\|_{1}$, we have

$\|p^{*\kappa}-p^{\kappa}\|_{1}<\frac{\sqrt{8}}{L^{1/2}n^{2}}+|\Delta_{\kappa}|\leq\frac{\sqrt{8}}{L^{1/2}n^{2}}+\frac{1}{Ln^{4}}<\frac{\sqrt{8}+1}{L^{1/2}n^{2}}<\frac{4}{L^{1/2}n^{2}}.$

For $L=64$, we have $|p^{*}_{i}-p_{i}|<1/2n^{2}$, as required. ∎

## Acknowledgments

We thank Satoru Fujishige for bibliographical information on principal
partition. The first author was supported by Grant-in-Aid for JSPS Research
Fellow, Grant No. JP19J22605, Japan. The second author was supported by JST
PRESTO Grant Number JPMJPR192A, Japan.

## References

  * [1] E. Aas, Limit points of the iterative scaling procedure. Annals of Operations Research 215 (2014), 15–23.
  * [2] A. Beck, First-Order Methods in Optimization. Society for Industrial and Applied Mathematics, 2017.
  * [3] L. Bregman, The relaxation method of finding the common point of convex sets and its application to the solution of problems in convex programming. USSR Computational Mathematics and Mathematical Physics 7 (1967) 200–217.
  * [4] I. Csiszár, I-divergence geometry of probability distributions and minimization problems. The Annals of Probability 3 (1975) 146–158.
  * [5] I. Csiszár and G. Tusnády, Information geometry and alternating minimization procedures. Statistics & Decisions, Supplement Issue No.1, (1984), 205–237.
  * [6] C. Franks, Operator scaling with specified marginals, preprint, (2018), arXiv:1801.01412, the conference version in STOC 2018.
  * [7] A. L. Dulmage and N. S. Mendelsohn, Coverings of bipartite graphs. Canadian Journal of Mathematics 10 (1958), 517–534.
  * [8] S. Fujishige, Submodular Functions and Optimization, 2nd Edition. Elsevier, Amsterdam, 2005.
  * [9] S. Fujishige, Theory of principal partitions revisited. In: W. J. Cook, L. Lovász, and J. Vygen, eds, Research Trends in Combinatorial Optimization, Springer, Berlin, 2009, pp.127–162.
  * [10] G. Gallo, M. D. Grigoriadis, and R. E. Tarjan, A fast parametric maximum flow algorithm and applications. SIAM Journal on Computing 18 (1989), 30–55.
  * [11] A. Garg, L. Gurvits, R. Oliveira, and A. Wigderson, Operator scaling: theory and applications. Foundations of Computational Mathematics (2019).
  * [12] A. Garg and R. Oliveira, Recent progress on scaling algorithms and applications preprint, (2018), arXiv:1808.09669.
  * [13] C. Gietl and F. P. Reffel, Accumulation points of the iterative proportional fitting procedure. Metrika 76 (2013), 783–798.
  * [14] L. Gurvits, Classical complexity and quantum entanglement. Journal of Computer and System Sciences 69 (2004), 448–484.
  * [15] L. Gurvits and P. N. Yianilos, The deflation–inflation method for certain semidefinite programming and maximum determinant completion problems. Technical report, NEC Research Institute, (1998).
  * [16] M. Idel, A review of matrix scaling and Sinkhorn’s normal form for matrices and positive maps, preprint, (2016), arXiv:1609.06349.
  * [17] N. Linial, A. Samorodnitsky, and A. Wigderson, A deterministic strongly polynomial algorithm for matrix scaling and approximate permanents. Combinatorica 20 (2000), 545–568.
  * [18] L. Lovász and M. Plummer, Matching Theory. North-Holland, Amsterdam, 1986.
  * [19] T. Matsuda and T. Soma, Information geometry of operator scaling, preprint, (2020), arXiv:2005.01453.
  * [20] K. Murota, Matrices and Matroids for Systems Analysis. Springer-Verlag, Berlin, 2000.
  * [21] U. G. Rothblum and H. Schneider, Scalings of matrices which have prespecified row sums and column sums via optimization. Linear Algebra and Its Applications 114/115 (1989), 737–764.
  * [22] R. Sinkhorn, A relationship between arbitrary positive matrices and doubly stochastic matrices. Annal of Mathematics Statistics 35 (1964), 876–879.
  * [23] R. Sinkhorn and P. Knopp, Concerning nonnegative matrices and doubly stochastic matrices. Pacific Journal of Mathematics 21 (1967), 343–348.
  * [24] N. Tomizawa and S. Fujishige, Historical survey of extension of the concept of principal partition and their unifying generalization ot hypermatroids. System Science Research Report No. 5, Department of Systems Science, Graduate School of Science and Engineering, Tokyo Institute of Technology, 1982.

##  Appendix A Appendix

### Proof of Lemma 2.2.

(1). Let $q_{j}:=\sum_{i}M_{ij}$. For any $N\in{\cal C}$, we have

$\displaystyle D(M\|N)=\sum_{i,j}M_{ij}\log
M_{ij}/N_{ij}\geq\sum_{i,j}M_{ij}\log\sum_{i}M_{ij}/\sum_{i}N_{ij}$
$\displaystyle=\sum_{i,j}M_{ij}\log q_{j}/c_{j}=\sum_{i,j}M_{ij}\log
M_{ij}/C(M)_{ij}=D(M\|C(M)),$

where we use the convexity of $D$ and $C(M)_{ij}=(c_{j}/q_{j})M_{ij}$.

(2). Let $p_{i}:=\sum_{j}N_{ij}$. For any $M\in{\cal R}$, we have

$\displaystyle D(M\|N)\geq\sum_{i,j}M_{ij}\log\sum_{j}M_{ij}/\sum_{j}N_{ij}$
$\displaystyle=\sum_{i}r_{i}\log r_{i}/p_{i}=\sum_{i,j}R(N)_{ij}\log
R(N)_{ij}/N_{ij}=D(R(N)\|N),$

where $R(N)_{ij}=(r_{i}/p_{i})N_{ij}$.

### Proof of Proposition 2.3.

The 5-point property is obtained by adding the 3-point and 4-point properties.
The 3-point property follows from

$\displaystyle D(M\|N_{k})=\sum_{i,j}M_{ij}\log\frac{M_{ij}}{(N_{k}{\bf
1})_{i}(M_{k})_{ij}/r_{i}}$ $\displaystyle=D(M\|M_{k})+D(r\|N_{k}{\bf
1})=D(M\|M_{k})+D(M_{k}\|N_{k}).$

The 4-point property follows from

$\displaystyle
D(M\|N_{k+1})-D(M\|M_{k})=\sum_{i,j}M_{ij}\log\frac{(M_{k})_{ij}}{(N_{k+1})_{ij}}=\sum_{j}\sum_{i}M_{ij}\log\frac{(M_{k})_{ij}}{C(M_{k})_{ij}}$
$\displaystyle=\sum_{j}(M^{\top}{\bf 1})_{j}\log\frac{(M_{k}^{\top}{\bf
1})_{j}}{(N^{\top}{\bf 1})_{j}}=D(M^{\top}{\bf 1}\|N^{\top}{\bf
1})-D(M^{\top}{\bf 1}\|M_{k}^{\top}{\bf 1})$ $\displaystyle\leq D(M^{\top}{\bf
1}\|N^{\top}{\bf 1})\leq D(M\|N),$

where we use the log-sum inequality.