Differentiable extensions with rounding guarantees for combinatorial optimization over permutations
Authors: Robert R Nerem
Zhishang Luo
Akbar Rafiey
Yusu Wang
Venue: NeurIPS 2025 poster ()
Keywords: Combinatorial optimization
Permutations
Birkhoff-von-Neumann decomposition
matching
quadratic assignment
Primary area: optimization
TL;DR: We introduce the Birkhoff Extension, a polynomial-time, almost-everywhere differentiable relaxation of permutation objectives to doubly stochastic matrices, enabling efficient combinatorial optimization.

Abstract: Continuously extending combinatorial optimization objectives is a powerful technique commonly applied to the optimization of set functions. However, few such methods exist for extending functions on permutations, despite the fact that many combinatorial optimization problems, such as the quadratic assignment problem (QAP) and the traveling salesperson problem (TSP), are inherently optimization over permutations. We present Birkhoff Extension (BE), an almost-everywhere-differentiable continuous polytime-computable extension of any real-valued function on permutations to doubly stochastic matrices. Key to this construction is our introduction of a continuous variant of the well-known Birkhoff decomposition. 
Our extension has several nice properties making it appealing for optimization problems. 
First, BE provides a rounding guarantee, namely any solution to the extension can be efficiently rounded to a permutation without increasing the function value. Furthermore, an approximate solution in the relaxed case will give rise to an approximate solution in the space of permutations. 
Second, using BE, any real-valued optimization objective on permutations can be extended to an almost-everywhere-differentiable objective function over the space of doubly stochastic matrices. 
This makes our BE amenable to not only gradient-descent based optimization, but also unsupervised neural combinatorial optimization where training often requires a differentiable loss. 
Third, based on the above properties, we present a simple optimization procedure which can be readily combined with existing optimization approaches to offer local improvements (i.e., the quality of the final solution is no worse than the initial solution). 
Finally, we also adapt our extension to optimization problems over a class of trees, such as Steiner tree and optimization-based hierarchical clustering.
We present experimental results to verify our theoretical results on several combinatorial optimization problems related to permutations.

Differentiable Extensions with Rounding Guarantees for Combinatorial Optimization over Permutations
Robert R. Nerem UCSD FrneremQucsd.ced edu
Zhishang Luo UCSD zluo@uesd.cm edu
Akbar Rafiey NYU ar95300nyu.edu
Yusu Wang UCSD yusuwanguucdd.edu
Abstract
Continuously extending combinatorial optimization objectives a powerful tech- nique commonly applied to the optimization of set functions. However, few such methods exist for extending functions on permutations, despite the fact that many combinatorial optimization problems, such the quadratic assignment problem (QAP) and the traveling salesperson problem (TSP), are inherently optimization over permutations. We present Birkhoff Extension (BE), almost-everywhere- differentiable continuous polytine-computable extension of any real-valued func- tion permutations doubly stochastic matrices. Key this construction is introduction of a continuous variant the well-known Birkhoff decomposition. Our extension has several nice properties making it appealing for optimization prob- lems. First, BE provides a rounding guarantee, namely any solution to the extension can be efficiently rounded Furthermore, an approximate solution in the relaxed case will give rise to an ap- proximate solution the space of permutations. Second, using BE, any real-vall real-valued optimization objective on permutations can be extended to an almost-everywheet differentiable objective function over the space of doubly stochastic matrices. This makes our BE amenable to not only gradient-descent based optimization, but also unsupervised neural combinatorial optimization where training often requires a differentiable loss. Third, based on the above properties, we present simple optimization procedure which can be readily combined with existing optimization approaches offer local improvements (i.e., the quality of the final solution is no worse than the initial solution). Finally, we also adapt our extension to optimiza- tion problems hierarchical clustering. We present experimental results to verify our theoretical results several combinatorial optimization problems related to permutations.
permutation without increasing the function value.
class of trees, such Steiner tree and optimization-ben based
Introduction
Continuously extending combinatorial objectives is common technique in combinatorial optimiza- tion, e.g., relaxation through linear programming, which can offer efficient optimization algorithms. Continuous extensions are particularly useful if they are differentiable, making them amenable to gradient- optimization methods. However, it often non-trivial to develop continuous exten- sions with theoretical guarantees that relate the optimization of the extension to the optimization of the combinatorial objective. In this paper, we consider combinatorial optimization problems where the goal is to minimize real-valued functions permutations. We aim develop extensions for functions on permutations with theoretical guarantees that allow for gradient- based optimization of these combinatorial objectives.
Optimization of functions on permutations is a setting that encompasses many combinatorial opti- mization problems with important applications. One of the most famous permutation optimization
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

problems is the traveling salesperson problem (TSP), where the aim is to find the order in which set of cities should be visited in tour minimize the length of the tour. TSP is an example of vertex ordering problem, class that contains many permutation optimization problems, in which the goal is to find an order of vertices in a graph that minimizes some objective. Examples of such problems are: directed feedback are set (DFASP), graph cutwidth, and minimum linear arrangement (MLA). Another essential permutation optimization problem is the quadratic assignment problem (QAP), in which a bijection between facilities and locations is sought that minimizes a quadratic objective function (such bijections can be identified with permutations).
1.1 Our work
In this paper develop techniques that solve combinatorial optimization problems over permutations. Formally, let Pn. denote the set of n x n permutation matrices (w.loo.g. permutations are viewed as matrices). For a given real-value objective function f Pn + R2 the goal is find
min S(P) S.t. P
Searching over the discrete space of permutation matrices is computationally prohibitive, It therefore appealing to turn our focus to continuous space, such the convex hull of permutation matrices D,,. This space, known as the Birkhoff polytope, is the set of matrices that are doubly stochastic, i.e., have non-negative entries and rows and columns summing one. The seminal result of Birkhoff [8] states that any matrix in the Birkhoff polytope can be decomposed into a convex combination of permutation matrices. This provides a probabilistic interpretation: any matrix A in the Birkhoff polytope arises as distribution (A) over permutation matrices A = (Rp=((A) P. We can define our Birkhoff extension as F(A)= = IE [A)=(())((()) and solve
min F(A) S.t. AED,
The challenge is that a Birkhoff decomposition (the distribution v(A) is not unique and the resulting extension may not be continuous and differentiable. Indeed, previous techniques for performing Birkhoff decomposition 8 [9] are not continuous.
differentiable Birkhoff extension. Our first theoretical contribution a Birkhoff decomposition
that is continuous and almost everywhere (a.c.) differentiable (Thm. 2.4). which gives rise to an extension that is also continuous and a.c. differentiable (Property D Continuity achieved by utilizing an arbitrary but fixed total order of permutations and decomposing according to this order. However, computing the extension with respect an arbitrary ordering of permutations is intractable. To solve this, we introduce score matrix S, order permutations by their inner product with S, and show that the resulting continuous and differentiable Birkhoff decomposition / extension can be efficiently computed (Property 2)
Optimization-enablinn properties. Another appealing strength of Birkhoff extension is that minima of the extension F directly correspond to minima of the function f over permutations (Property Bh. Moreover, Birkhoff extension admits a scheme (Property (4) for rounding a doubly stochastic matrix A to a permutation P' that is guaranteed not to degrade the quality of the solution, i.c., f(P) < (A)) These properties ensure that optimizing (or approximating) the extension yields optimal (or approximate) solutions to f.
Different choices for S yield different continuous Birkhoff extensions, which is valuable flexibil- ity. Interestingly, we can choose any permutation P (e.g, an approximate solution produced by a comparatively fast algorithm) to produce score matrix (Property 4 4 which can then be combined with optimization and rounding to produce solutions that are at least as good as P. That is, given any existing solution to a combinatorial optimization problem, we can use it as the score matrix for Birkhoff extension to further improve it, thereby yielding a local improvement procedure.
Optimization algorithm. Given the a.e. differentiability of Birkhoff extension, gradient. However, gradient descent cannot be directly applied to optimization over the Birkhoff polytope, since after each step the resulting matrix may not be doubly stochastic. We propose Frank- Wolfe algorithm that overcomes this issue and preserves double stochasticity throughout optimization.
compute its
Birkhoff extension not necessarily convex and, thus, gradient- based optimization could converge to local minima. We alleviate this issue by changing the score matrix whenever the optimization

---

converges to local minimum. Changing the score matrix changes the extension being optimized, potentially changing the extension one where the current iterate is not at local minimum, allowing further optimization, We further show that for specific changes to the score matrix (based on Property we are guaranteed the quality of the rounded solution does not decrease.
Application to neural combinatorial optimization. Birkhoff' extension can also be used for unsu- pervised neural combinatorial optimization where training often requires that the objective function is differentiable. Neural approaches are promising new paradigm in combinatorial optimization as, unlike traditional techniques, they inherently leverage the distribution of problem instances being solved [7]. However, supervised neural combinatorial optimization is often prohibitively expensive as it requires computing exact solutions to create labels for training. Unsupervised learning, such as the set extension proposed in [29 (28). circumvents this issue by removing the need for labels. In
App. B similar to [29], we propose an unsupervised neural approach based on our Birkhoff extension.
The properties of Birkhoff extensions offer advantages for unsupervised learning in that rounding
guarantees ensure the minima sought in training correspond well with the combinatorial objective.
Beyond permutation functions. It is compelling to consider when similar techniques can be applied other combinatorial functions. In App.[C] we present analysis for applying Birkhoff extension to functions on rooted binary trees over a fixed set of leaves. Optimization of these tree functions arises in many combinatorial optimization problems such as Steiner tree problems and hierarchical clustering [17]
Experiments. In Section we perform experiments quadratic assignment (QAP), traveling salesperson (TSP), and directed feedback arc set (DFASP) problems, showing Birkhoff extension is an effective approach for optimizing permutation functions. Our method is especially effective for the quadratic assignment problem (QAP), where it outperforms a range of Gurobi implementations and common heuristics.
To summarize, our main contributions are:
A novel Birkhoff decomposition, which has continuous, a.e. differentiable, and efficiently com- putable coefficients.
A continuous a.e. differentiable extension of permutation functions to real-valued functions on the Birkhoff polytope, with properties, such as rounding guarantees, desirable for optimization.
A pooreticalll-justffiee permutation-functionn optimization procedure that combines this extension with gradient- optimization.
Combination of Birkhoff extension with an NN to yield an unsupervised neural optimizer (App. B
1.2 Related work
Extensions and optimization. The optimization literature often focuses on building extensions with desirable optimization properties, particularly convexity and concavity [15
38
48
A classical
approach extending a discrete set function f {0, 13" IR is by computing the convex closure, which is the point-wise supremum over linear functions that lower bound f 122 48 detailed comparison with this approach and ours is made in App. H prominent example of successful application of these methods to combinatorial optimization is submodular functions. The convex closure for a submodular function identical the Lovász extension [34), also known as the Choquet integral in decision theory 13. which leads to polynomial-tine algorithms for submodular minimization 255). series of works []] 53 introduced and studied multilinear extension of submodular functions which results in approximation algorithms for certain constrained submodular maximization problems. See [5] for further details on extensions submodular functions. Convex extensions have been applied to a broader class of set functions beyond submodular ones [21]. as
well as to combinatorial penaltics with structured sparsity 200 40).
Extensions and neural combinatorial optimization. Unsupervised learning for combinatorial
optimization problems has recently attracted great attention 28 46 (50. Many frameworks based on RL B1 B3
18, 6 [5 or supervised learning 52, 24 do not hold any special requirements on
the formulation of combinatorial problems. However, these approaches often suffer from dependence on labeled data or unstable training, respectively. In contrast, unsupervised learning for combinatorial optimization problems, where continuous relaxations of discrete objectives are utilized, is superior in its faster training, good generalization, and strong capability of dealing with large-scale problems.
3

---

The general idea is to use, loss function, function on a continuous domain that extends the discrete function. Notable examples of these types of work are 28 54 relaxation of discrete functions are used for the loss.
where a probabílistic
Neural set function extension. Karalias al. 29 propose several continuous and differentiable extensions for training neural networks to optimize set functions, such as for the Max Clique and Max Independent Set problems. These extensions for the efficient rounding schemes. However, many combinatorial optimization problems, such as the optimization of permutation functions, have no natural formulation as set function optimization. We use this convex decomposition framework as inspiration for our extension of permutation functions.
formed by convex combinations, which allows
2 Birkhoff Extension
In this section we introduce our continuous and a.e. differentiable Birkhoff decomposition. We then use this decomposition to construct extension of permutation functions to the Birkhoff polytope. Finally, we show this extension has many advantageous properties. Proofs of the claims made in this section are given in App.
2.1 Preliminaries
doubly stochastic n matrix is one with non non-negative entries where each row and column sums to 1. A permutation matrix is special doubly stochastic matrix with binary entries and single 1 in every row and column, The class of n doubly stochastic matrices is convex polytope known as the Birkhoff polytope Dn.- The Birkhoff polytope lies in an (n 112-ddimeniooss affine subspace of n² -i-dimensional Euclidean space defined by 2n independent linear constraints specifying that the row and column sums all equal 1. Let Pn, denote the set of
permutation matrices.
Theorem 2.1 (Birkhoff decomposition 8 Any doubly stochastic matrix A E Dn, can be decom- posed as A
Σ k=1 agPh where M < n2
1, as 0, Σ k Og and Ps E P..
To construct this decomposition we view as the biadjacency matrix of a bipartite graph with vertices [n] L [n] and edges (i,j) of weight A(i,j). In this graph consider the following set of permutations, ie. perfect matching in this graph, that do not have edges of weight zero
Definition 2.2. permutation matrix E Pa is a perfect matching of non-negative matrix A iff' P(i,j) = implies A(i,j) > 0. We denote the space of permutations that are perfect matchings of A as P(A).
The standard algorithm 8 for constructing such decomposition is given in Alg. [ This is an iterative algorithm that starts from Bo = A and, at each iteration k. takes the matrix Be resulting from the previous step, which is proportional to doubly stochastic matrix, and finds permutation P that a perfect matching of Bg. The existence such a matching is a consequence of B, being proportional to doubly stochastic matrix and Hall's marriage theorem. Furthermore. this matching P can be computed using a standard bipartite matching algorithm in Onn²) time. If P is a perfect matching of B₁ and is the value of the smallest entry (i, in B& such that P(i.j)
1, then B++1 = B& - a is a matrix proportional to a doubly stochastic matrix with one less non-zero entry than in Bg. Note that since P is matching of BR, have a > 0. This process is repeated until the resultant matrix is the zero matrix.
2.2 A continuous and a.e. differentiable Birkhoff decomposition
We extend function Pa R on permutations function F D,, + E Birkhoff polytope via the Birkhoff decomposition: F(A) = E. (Pk). However, Birkhoff decomposition non-unique; there may be many different ways to represent doubly stochastic matrix as a convex combination of permutations. This non-un -onuuniuuessss is evident at each step of the decomposition the multiple choices of which permutation matrix P E P(Bx) to subtract. We particular decomposition so that the coefficients ag and the extension F, are continuous functions of the matrix A being decomposed (with the standard La-induced topology on D.).
describe how to fix
The key insight for this construction is to fix an arbitrary total order over permutation matrices and, at each step in the decomposition, always pick the valid permutation that comes first in the order.

---

Algorithm 1 Classical Birkhoff decomposition [8]
Algorithm 2 Continuous Birkhoff decomposition
Require: AED Dn
Require: A E Dn. identifying score matrix S
Ensure: ((ax, ]=))==] s.t. A Σ 0 Pk. Ensure: P kr k= s.t. A= = Ex=1 axPi Ex az = 1, and 0% > O.
Σx a2 = 1, and 0%k > O.
+ Bg < A
h 1, By + +
while B3 $ 0 do
while B& # do
Ps f P PP(Bk)
Ps t argmaxpeP(B.) (P,S)
0₁ f min.y B (i. Pa(i.)) = =
a +- miney M BBx (i,j) Pa(i,.) = 1}
B&+1 < Bb Px
B&+1 < BB or Pk
k +++
++
end while
end while
M +- k return {(aR, P=)
M + K
return ((ax,,P)]
Previous decomposition algorithms fail to achieve continuity because small changes to the matrix A being decomposed could change which permutation is subtracted at each step, which alter the trajectory of the decomposition. By fixing the order in which permutations decomposition, we circumvent this issue. Below we introduce a continuous Birkhoff decomposition scheme, and prove its correctness (i.e., validity and continuity) in Thm. Th4.2.4
subtracted in the
Definition 2.3 (Continuous Birkhoff Decomposition), Given enumeration PPyee Pn (i.e, fix a total order of all permutations), and given "D, the continuous Birkhoff decomposition of A induced by {Pr) is (ag, Pryz= where the coefficients are defined recurrently from E = 1 to n! in order by
8-1
Oe min if
A(I.))
Σ
Om Pmo (1,5) Pe((.))=11
m=1
Theorem 2.4. Given enumeration (Pt) and given A E Dno the coefficients of the continuous Birkhoff decomposition (o Pr)t of A are (i) Lipschitz continuous functions from DD, to R. (ii) all non-negative and sum 1, and (iii) yield a valid decomposition of via A >n! aePt. Furthermore, (iv) there are at most n²² n + coefficients being non-zero.
The a.e. differentiability of the continuous Birkhoff decomposition follows from Lipschitz continuity and an application of Rademacher's s theorem 43
Theorem 2.5. The coefficients {av !u! of the continuous Birkhoff decomposition are almost every- where differentiable functions from D, to R.
Now that we have constructed continuous Birkhoff decomposition, the next questions are how to represent that total order of all permutation matrices, and how to compute this decomposition efficiently. Indeed, for an arbitrary ordering of permutations, efficient computation is not feasible as it requires referencing the order of rz! elements. We overcome this challenge by instead focusing on orderings of permutations that arise from an inner product.
Definition 2.6 Score-Induced Birkhoff Decompositions). Given an n matrix S. the score of permutation P is (S, P) = E.j P(i.j) We call S score matrix, and say that S is identifying if it assigns a unique score to every permutation, thereby inducing a total order on Pn.
Furthermore, given identifying score matrix S, the Birkhoff decomposition as specified in Def. 2.3 with respect to an ordering of permutations by their score {Pe) is called an S-induced Birkhoff decomposition.
A simple example of identifying score matrix is given by S(i.)) = (i+))) The score-matrix induced total order is particularly effective because for any A E Dn, the permutation P E P(A) that comes first in this order can be found efficiently by solving a maximum weight matching problem. Consequently, decompositions with respect to this order can be constructed efficiently:
Theorem 2.7. Given an identifving score matrix S, the S-induced Birkhoff decomposition can be computed in O Onr) time by Alg. 2
5

---

In practice, we may wish to use score matrices other than S(i.j) 2(++n)) The following theorem shows that random assignment of S is sufficient for be identifying.
Claim 2.8. If the entries S are independent absolutely continuous random variables S(i.)) e R then S is identifying almost surely.
2.3 Properties of Birkhoff extension
We present several properties that make our score-induced Birkhoff extension a desirable candidate for optimization. Proofs of these properties are given App. A
Definition 2.9. Given A E Dn and ordering of permutations (P& )ol let (ag, Ps) be the non-zero Birkhoff coefficients defined in Def. 2.3 For any f Pa - IR, the Birkhoff extension f is the function F: Dn R. where
F(A)= Σ aPf(Ps)
(4)
We say F is score induced or S-induced if the ordering of permutations is induced by S. We sometimes emphasize the dependence on by using F' to denote the S-induced Birkhoff extension.
Almost everywhere differentiability and continuity of F are essential for gradient-based optimization. These properties follow from continuity and a.e. differentiability of the coefficients {ar {a Further- more, we show that when computing the gradient of F one only needs to consider the non-zero terms the Birkhoff decomposition.
Property Birkhoff extensions are Lipschitz continuous and almost everywhere differentiable. Fur- thermore if L+ = {f E [n!]] az > 0}, then V (AA) =Ete,, (V ( almost everywhere.
Computing a Birkhoff extension reduces to computing the corresponding Birkhoff decomposition, thus, Alg. 2 gives efficient computation of score-induced Birkhoff extensions (see Thm. 2.7
Property 2. Score-induced Birkhoff extensions F can be computed in O(n²) time.
One concern when optimizing a continuous extension to a combinatorial function f is that the optimization reaches some minimum in the extended space that does not correspond with minima of the combinatorial function (which is our true goal). Property| Seneertyybeloo below shows that (global) minima of the extension F (i.e., over D,)) are related to those of over the permutations.
Property 3. Let F be a Birkhoff' extension of / Pn + R. Then (1) minpeP, S(P) = MiD F(A) and (2) argmin F(A) C Cany(argminpep._ ((P)..
The Birkhoff decomposition leads to a simple rounding strategy:
Definition 2.10. Given matrix Dn and score matrix that induces a Birkhoff decomposition (ow, Pr we define
rounds(A) 4)= prgminn.
(5)
Note that rounds(A) can be computed Onn²) time by computing Birkhoff decomposition of A. The rounding scheme is lossless in that can only improve solution quality; see Property [ below. Consequently, optimizing, f reduces to optimizing the Birkhoff extension F. as any minimum of F can be used to derive minimum of f. Furthermore, approximations to minAeD, (As(A) can be rounded to approximations for minpep. f(P).
Property 4. Let F's be a score-induced Birkhoff extension of f: Pn - R. For any E D,, then
1. (ro)nd < Fs(A). Furthermore, if A is C-approximation for minAeP, Fs(A), then round (A) a C-approximation for minpep, f(P).
2. If P* E Pn, with max;y) PP (i.j) - S(i.5) < 2n' then froonds $ f(P".
Another useful quality of this rounding scheme, Property 2 is that if sufficiently close some permutation P. then rounding always yields a solution at least as good as P. This holds independent of the matrix that being rounded and gives useful flexibility to Birkhoff extension optimization. In particular if the score S is close to an approximate solution to the combinatorial optimization problem, then rounding always produces solution at least as good as the approximation. For example, S P apppoxx Q where P Papprox is solution produced by a fast approximation algorithm and Q is

---

Algorithm 3 Static score Frank-Wolfe over Dn Algorithm 4 Dynamic score Frank- over
Require: P + R, random A E Dn, score D,
Require: P + R. random A E Dn. score S for 1..-
for 1 -= 11... 1....TT do
Pe 4 argmax PEP. (VPs(A)), P)
P: f- argmax (VPs(Ai);
A,+1 + (1 A. A, + Pe
A++1 f (1 - ) XPP
p 4 p U Birkhoff
end for return rounds(Ar)
if Muudate__coon
then Q ~ [0. _
P Q arrmin (PP)
S S
end if end for return rounds(Ar)
random noise matrix with entries in [0, 1//]/. (Adding Q ensures that S is almost surely identifying by Claim 2.8) Then, gradient-based optimization of the Birkhoff extension associated with S finds solutions that are guaranteed to be no worse than the approximation Paperox. This means that we can use Birkhoff extension as local improvement strategy to potentially improve any given solution. It important that P Papprox used as the score matrix, not as an initialization for the optimized matrix A E D,, In fact, initializing A at P Papprox yields weaker guarantees as optimization may produce solution of worse quality this case.
2.4 Optimization procedure with dynamic score
Consider the optimization problem min fP) where the goal is to minimize a function f : Pa - R on permutations. A natural relaxation for this optimization problem is to optimize the Birkhoff extension of f over the set of doubly stochastic matrices; namely, the constrained optimization problem of the form min AgD, F(A). Here we propose iterative first-order optimization algorithm that is concerned with optimizing an objective function over the Birkhoff polytope Dn,-
Frank-Wolfe One difficulty in gradient-based approaches to optimize constrained optimization problems is the risk of stepping outside of the feasible region. We address this difficulty by adapting the famous Frank-Wolfe approach (23). In contrast to projected gradient descent approaches, the idea to follow a direction of descent that is best aligned with the negative of the gradient for which we can also easily ensure feasibility. This is done via optimizing the negative of the gradient over the extreme vertices P.. C DD, and then taking the obtained permutation as an alternative direction of descent. The overall process is outlined in Alg. Moreover, Frank-Wolft reduces to iterated linear minimizations for which highly efficient combinatorial algorithms exist, whereas projection- baasd methods require solving quadratic programs over the Birkhoff polytope. In particular, argmin (VFs(A)) P) can be computed by finding a maximum weight matching in the bipartite graph G that has vertices [n] L [n] and an edge from vertex to vertex j of weight FFs (i. the weight of matching P in G is (V I's
(A,) (A),, P). Similar applications of Frank- Wolfe to continuous optimization over the
Birkhoff polytope have been discussed by Tewari al. 49 and Jaggi 27). however these works do not consider applications to combinatorial optimization, as we do.
Dynamic score. One issue with gradient-based optimization of Birkhoff extension is that the extension is not necessarily convex, and thus, an optimization algorithm may converge to local minima. Even if the function is convex, analyzing the convergence of Frank-Wolfs type algorithms for non-smooth functions is highly nontrivial [3] 42 4 and, in general, they may not converge at all, see [39] for a counterexample. To alleviate these issues, we propose an optimization scheme in which the score matrix S is frequently changed to help escape local minima. The key to the efficacy of this approach that the score matrix can be changed without decreasing the quality of the rounded solution. Property 2 gives conditions for this change to be made without harming the optimization. In particular, this occurs under the condition that the new score matrix S' sufficiently close to permutation P satisfying f(P) < (round (A)). We further theoretically validate this approach by showing at any D, there is score matrix such that is not a local minimum of Fs (App. G
Theorem 2.11 (Escaping Local Minima). Let f Pa. + R. any function on permutations, there exists a score matrix S such that A is not a local minimum of Fs.
7

---

As finding such a score matrix is computationally hard, we use the following procedure in prac- tice. Let Birkhoff, s (A) be the permutations with non-zero coefficients the S-induced Birkhoff decomposition of A and for each iterate A, let p : s We update the score to 5" = prgmin f(P) + 2n Q where Q is a matrix with uniform random entries in [0. 1]. By Property 4 this update satisfies f(rounds: (At)) S (rounds (At).. This procedure is outlined in Alg. where update_scone is flag determined externaly indicating if convergence has occurred and the score should be updated.
Adaptations to constraints and trees. We provide further extensions based on our Birkhoff decomposition. First, show how Birkhoff extension can be adapted class of rooted binary trees over fixed leaf set (App. C). Next, we show that Birkhoff extension can easily handle broad class of simple constraint (App. Ip. In particular, we consider constraints that can be incorporated into the corresponding matching problem.
Unsupervised learning. We apply Birkhoff extension train neural combinatorial optimization solvers (App. B). As Birkhoff extension is a.e. differentiable, can be used as a loss function in unsupervised learning. Using this loss, neural networks can be trained to output solutions to combinatorial problems. Our proof-of-concept experiments show that training a neural network using Birkhoff extension yields model that predicts good approximate solutions. In particular, the output of this NN can be used to initialize our Birkhoff extension optimization algorithm to provide substantial speedups. Furthermore, this generalizes to larger problem instances than seen in training.
3 Experiments
We carry out experiments on three different combinatorial optimization problems: the quadratic assignment problem (QAP), the (Euclidean) traveling salesman problem (TSP), and the directed feedback are set problem (DFASP). For QAP, we seek a bipartite matching that minimizes a quadratic loss. For TSP, the goal is to find the shortest possible tour visiting each city once and returning to the starting point. In DFASP, we aim to find vertex ordering that minimizes the number of edges directed against the order. Detailed definitions and linear integer programming formulations for each problem are given in App.D) D For each problem instance, we apply a variant of Alg. to optimize and compare its performances with basclines. In particular, we update the score matrix every 10 epochs (i.e.. update_ score = True in Alg. i if and only if E {10, 20,... .., T)). To overcome the ((n)) time complexity, we truncate the Birkhoff decomposition only the first k = 5 terms. Automatic differentiation is used to compute gradients. Additional experiments and ablations are given in App.
E
Summary of results. We first summarize our results and then provide specific discussion (full details in App. E}. Our optimization algorithm performs best at QAP, where it outperforms two different Gurobi implementations and two popular heuristic approaches (which are the default scipy heuristics). We suspect Birkhoff extension performs particularly well here since permutation matrices naturally model matchings in that the i matrix entry indicates a match between facility i and location j. Other permutation problems are less amenable to this representation; e.g. in TSP the i matrix entry indicates that city i occurs jth in the tour. Indeed, Gurobi is highly effective for TSP and outperforms our Birkhoff extension-based approach both in terms of accuracy and efficiency. Nevertheless, we show that Birkhoff extension can still be used as a local improvement to further improve solutions provided by the minimum spanning tree (MST) approximation algorithm or by quadratic programming. Gurobi does not perform well for DFASP. In particular, our approach provides superior solutions for large instances for DFASP.
QAP. We evaluate on the QAPLIB library [0], which contains N = 136 problem instances with sizes between 12 and 256. Some instances have confirmed optimal objective values, while others have the best-known objective values reported. For Gurobi, we use the Settimann-Breccee lincarization 30 for the Mixed-Integer Linear Programming (MILP), and we set same maximum time budget for Gurobi and our method as 2n seconds, and we report the actual run- auntimee on average for all methods. For the same time budget, the canonical ILP formulation [56] fails on 70% of instances. We also test the same benchmark with two other heuristics: Fast Approximate QAP [52) and 2-opt algorithm [16]. Our method outperforms all benchmarks (see results in Table IT)
TSP. We generate instances by uniformly sampling n vertices v E [0, 1J2. For each size n E {20, 30, 40, 50, 100), we generate N 50 instances. Comparisons are made with an MST-based
8

---

Method
QAP: Assignment Cost
Average Gap from the best known
Run-Time (sees)
Gurobi, BeeeekkknBrreckx
7.67%
90.3
Fast Approximate QAP 52
15.41 %
0.01
2-opt 16
10.38%
22.82
Random S Init. Alg.
6,30 %
23.01
Table 1: Performances of methods for QAP in terms of assignment costs. As QAPLIB contains the best known solutions for each problem instance, we report the average gap to these objective values.
approximation algorithm [14]. a quadratic programming (QP) relaxation, and Gurobi (which is optimal). For Alg. we three different score matrix S initializations: uniformly random in [0, 10. MST approximation, and the QP relaxation solution. The learning rate is rj = 0.01. The maximum number of optimization steps is T 10000, and early stopping triggers after 2000 steps without improvement. All optimization processes converge (see optimization curves for different heuristics in App. E}. Results for TSP are presented in Table Although Gurobi is highly effective and returns the optimal solution, our approach can still provide local improvements over QP and MST (see the last
in Table by using those solutions the input score matrix for Alg. 4
DFASP. We generate instances using a directed Erdös-Rénye model with p E {0.1, 0.5, 0.9). Again, for each problem size E {20, 50, 100}. we generate N = 50 instances. Note, Gurobi is not able to find an optimal solution within a reasonable time for DFASP. Hence, to compare the quality of Gurobi vs. our approach, we limit the runtime of both algorithms to Tú minutes for equal comparison. Consequently, the number of steps T varies for different problem instances. (Additional experiments (App. E show these optimization algorithms offer minimal improvement after the 19 minute time limit.) The learning rate is 17 = 0.005. Results, given in Table 3 show our method is competitive with Gurobi and shows superior performance for larger instances.
TSP: Tour Length +
Method
11
20
30
40
50
100
Gurobi
3.889
4.531
5.190
5.707
7.748
MST
4.746
5.784
6.835
7.410
10.288
QP
4.553
5,666
6.818
7.782
11,896
MCTS W. Rand. Simulation
4.558 (0.5s)
5.563 (1.7s)
6.706 (4.7s)
7.468 (9.8s)
11.781 (113.37s)
Random S Init. Alg. 4 MST S Init. Alg. (Improv.)
4.518
5.681
7.139
8,433
15.615
4.345 (8.33 %)
5.290 (8.53 %)
6.328 (7.42
6.892 (6.99
9.836 (4.60 %)
QP Init. Alg. IImprvv)) 4.405 (3.25 %)
5.633 (0.58 %) 6.7 (1.60 %)
7.670 (1.45 %)
11.782 (0.96 %)
Table 2: Performances of methods for TSP in terms of tour length. Best results are marked as bold; second bests are blue. If score is initialized to an approximate solution, percentages indicate proportional improvement Alg. 4 oveer initialization. See appendix (E optimization
DFASP: Feedback-Are Set Size +
Method
p
0.1
0.5
0.9
72
20
50
100
250
20
50
100
250
20
50
100
250
Gurobi
3.8
44.8
269.1
875.6
67.7
475.4
2229.9
6831.7
156.7
1028.0
43322.7
13652.4
Random S Init. Alg.
4.7
53.0
289.0
991.3
65.0
496.5
2123.4
6459.8
157.5
1039.0
4253.7
12198.7
Table 3: Performances of DFASP algorithms in terms of feedback-are set cardinality.

---

4 Concluding remarks
We present Birkhoff extension, a continuous a.c. differentiable extension of permutation functions to doubly stochastic matrices, which has rounding guarantees. Combining this extension with a gradient- based optimization algorithm, we develop an iterative optimization framework for permutation functions. We present experiments to validate approach for combinatorial optimization problems.
We propose neural optimizer based on our Birkhoff extension, however, this direction requires further exploration and far more extensive experiments. We leave this as a future direction to investigate. We also note that currently we use a simple strategy to update the score matrix. It would be interesting to explore more effective update strategies. Finally, computing Birkhoff extension is expensive (O(n²) time complexity). although in practice, the number of permutations is usually far fewer than n² It will be interesting to investigate how to improve the time complexity, or how to obtain an updated the Birkhoff decomposition efficiently as the input matrix changes, given that one needs to compute this decomposition many times within our optimization framework. One potential for speedup is to develop matching algorithms that are more amenable to GPU computation. We initiate discussion by introducing GPU friendly matching algorithm in App. L
Acknowledgments
This work is partially supported by NSF (National Science Foundation) by grant CCF-2112665.
References
[1] Ryan P. Adams and Richard S. Zemel. Ranking via Sinkhorn propagation. ArXiv, abs/ 106.1925. 2011.
[2] Saeed Amizadeh, Sergiy Matusevych, and Markus Weimer. Learning to solve cireuit-SAT: An unsupervised differentiable approach. In International Conference on Learning Representations, 2019.
[3] Kamiar Asgari and Michael J. Neely. Projection-free non-smooth convex programming. arXiv preprint arXiv:2208.05512 05127, 2022.
[4] Kamiar Asgari and Michael J. Neely. Nonsmooth projection-free optimization with functional constraints. Computational Optimization and Applications, 89:927-975, 2024.
[5] Francis Bach. Submodular functions: from discrete to continuous domains. Mathematical Programming, 175:419-459, 2019.
[6] Irwan Bello, Hieu Pham, Quoe V Le, Mohammad Norouzi, and Samy Bengio. Neural combi- natorial optimization with reinforcement learning. In International Conference on Learning Representations, 2017.
[7] Yoshua Bengio, Andrea Lodi, and Antoine Prouvost. Machine learning for combinatorial optimization: methodological tour d'horizon. In European Journal of Operational Research, 2018.
[8] Garrett Birkhoff. Three observations on linear algebra. Univ. Nac. Tucumánn Revista A., 5:147-151, 1946.
[9] Fanchen Bu, Hyeonsoo Jo, Soo Yong Lee, Sungsoo Ahn, and Kijung Shin. Tackling prevalent conditions in unsupervised combinatorial optimization: Cardinality, minimum, covering, and more. In Ruslan Salakhutdinov, Zico Kolter, Katherine Heller, Adrian Weller, Nuria Oliver, Jonathan Scarlett, and Felix Berkenkamp, editors, Proceedings of the 41st International Confer- ence on Machine Learning, volume 235 of Proceedings of Machine Learning Research, pages 4629-4729, PMLR, Jul 2024.
[10] Rainer E. Burkard, Stefan E. Karisch, and Franz Rendl. Qaplib - quadratic assignment problem library. Journal of Global Optimization, 10(4):391-403, 1997.
10

---

[11] Gruia Calinescu, Chandra Chekuri, Martin Pal, and Jan Vondrák. Maximizing monotone submodular function subject matroid constraint. SIAM Journal on Computing, 40(6):1740- 1766, 2011.
[12] Gengjie Chen and Evangeline F. Y. Young. Dim sum: Light clock tree by small diameter sum. In Proceedings of the IEEE/ACM Design, Automation & Test in Europe Conference & Exhibition (DATE), pages 1 174-179, Florence, Italy, March 2019.
[13] Gustave Choquet. Theory of capacities. In Annales de l'institut Fourier, volume pages 131-295, 1954.
[14] Nicos Christofides. Worst-case analysis of new heuristic for the travelling salesman problem. Technical Report Report 388, Carnegie-Mellon University, Management Sciences Research Group, 1976.
[15] Yves Crama. Concave extensions for nonlinear 0-1 maximization problems. Matherattcall Programming, 61:53-60, 1993.
[16] G. A. Croes. A method for solving 812, 1958.
problems. Operations Research, 6(6):791-
[17] Sanjoy Dasgupta. cost function for similarity- hierarchical clustering. Proceedings of the forty-eighth annual ACM symposium on Theory of Computing, 2015.
[18] Arthur Delarue, Ross Anderson, and Christian Tjandraatmadja. Reinforcement learning with combinatorial actions: An application vehicle routing. Advances in Neural Information Processing Systems, 33:609-6200 -200, 2020.
[19] Fanny Dufossé and Bora Uçar. Notes Birkhoff-von Neumann decomposition of doubly stochastic matrices. Linear Algebra and its Applications, 497:108- -15, 2016.
[20] Marwa Halabi, Francis Bach, and Volkan Cevher. Combinatorial penalties: Which structures preserved by convex relaxations? In International Conference Artificial Intelligence and Statistics, pages 1551-1560. PMLR, 2018.
[21] Marwa El Halabi and Stefanic Jegelka. Optimal approximation for unconstrained non- submodular minimization. In International Conference on Machine Learning, pages 3961-3972. PMLR, 2020.
[22] James E Falk and Karla R Hoffman. A successive underestimation method for concave minimization problems. Mathematics of Operations Research, 1(3):251-259, 1976.
[23] Marguerite Frank and Philip Wolfe. An algorithm for quadratic programming. Naval research logistics quarterly, 3(1 (11-2):55-110, 1956.
[24] Maxime Gasse, Didier Chételat, Nicola Ferroni, Laurent Charlin, and Andrea Lodi. Exact combinatorial optimization with graph convolutional neural networks. Advances in Neural Information Processing Systems, 32, 2019.
[25] Martin Grötschel, László Lovász, and Alexander Schrijver. The ellipsoid method and its consequences in combinatorial optimization. Combinatorica, 1:169-197, 1981.
[26] Gurobi Optimization, LLC. Traveling salesman problem (TSP) example in python. https:// www. www.biurrb...om
httpr/aaaiiionuurreeeeeeeeeeae www.rreeteeeeapleeee/eee//e httpaapppessttppp.....lmll html 2024. Accessed:
2024-07-05.
[27] Martin Jaggi. Revisiting frank-wolfe: Projection-free sparse convex optimization. In Sanjoy Dasgupta and David McAllester, editors, Proceedings of the 30th International Conference on Machine Learning, volume 28 of Proceedings of Machine Learning Research, pages 427-435. PMLR, 2013.
[28] Nikolaos Karalias and Andreas Loukas. Erdos goes neural: an unsupervised learning framework for combinatorial optimization on graphs. Advances Neural Information Processing Systems, 1166699-667,,22 2020.
11

---

[29] Nikolaos Karalias, Joshua Robinson, Andreas Loukas, and Stefanie Jegelka. Neural set func- tion extensions: Learning with discrete functions in high dimensions. Advances in Neural Information Processing Systems, 153555388-155,, 2022.
[30] L.. Kaufman and F. Broeckx. An algorithm for the quadratic assignment problem using benders' decomposition. European Journal of Operational Research, 2(3):207-211, May 1978.
[31] Elias Khalil, Hanjun Dai, Yuyu Zhang, Bistra Dilkina, and Le Song. Learning combinatorial optimization algorithms over graphs. Advances in Neural Information Processing Systems, 30,
2017.
[32] Elias Khalil, Pierre Le Bodic, Le Song, George Nemhauser, and Bistra Dilkina. Learning to branch in mixed integer programming. In Proceedings of the AAAI Conference on Artificial Intelligence, volume 30, 2016.
[33] Wouter Kool, Herke Van Hoof, and Max Welling. Attention, learn to solve routing problems! arXiv preprint 10X3v:88030084/ 2018.
[34] László Lovász. Submodular functions and convexity. Mathematical Programming The State of the Art: Bonn 1982, pages 235-257, 1983.
[35] Haggai Maron, Heli Ben-Hamu, Nadav Shamir, and Yaron Lipman. Invariant and equivariant graph networks. In International Conference on Learning Representations, 2019.
[36] C. E. Miller, A. W. Tucker, and R. A. Zemlin. Integer programming formulation of traveling salesman problems. Journal of the ACM, 7(4):326-329, Oct 1960.
[37] Marcin Mucha and Piotr Sankowski. Maximum matchings via gaussian elimination. In 45th Annual IEEE Symposium on Foundations of Computer Science, pages 248-255. IEEE, 2004.
[38] Kazuo Murota. Discrete convex analysis. Mathematical Programming, 83:313-371, 1998.
[39] Y. Nesterov. Complexity bounds for primal-dual methods minimizing the model of objective function. Mathematical Programming, Series A, 171:311-330, 2018.
[40] Guillaume Obozinski and Francis Bach. Convex relaxation for combinatorial penalties, arXiv preprint arXiv: 1240, 2012.
[41] Michael O Rabin and Vijay V Vazirani. Maximum matchings in general graphs through randomization. Journal of algorithms, 10(4):557-567, 1989.
[42] Sathya N. Ravi, Maxwell D. Collins, and Vikas Singh. A deterministic nonsmooth frank wolfe algorithm with coreset guarantees. INFORMS Journal on Optimization, 1(2): 2019.
[43] Walter Rudin. Principles of Mathematical Analysis. 1964.
[44] Piotr Sankowski. Maximum weight bipartite matching in matrix multiplication time. Theoretical Computer Science, 410(44)44480 4488, 2009.
[45] Lex Schrijver. A Course in Combinatorial Optimization. 2003.
[46] Martin JA Schuetz, J Kyle Brubaker, and Helmut G Katzgraber. Combinatorial optimization with physics-inspired graph neural networks. Nature Machine Intelligence, 4(4):367-377, 2022.
[47] Hadar Serviansky, Nimrod Segol. Jonathan Shlomi, Kyle Cranmer, Eilam Gross, Haggai Maron, and Yaron Lipman. Set2graph: learning graphs from sets. Advances in Neural Information Processing Systems, 2020.
[48] Mohit Tawarmalani and Nikolaos V Sahinidis. Convex extensions and envelopes of lower semi-continuous functions. Mathematical Programming, 93(2):247-263, 2002.
[49] Ambuj Tewari, Pradeep K. Ravikumar, and Inderjit S. Dhillon. Greedy algorithms for struc- turally constrained high-dimensional problems. In Advances Neural Information Processing Systems, volume 24, pages 882-890, 2011.
12

---

[50] Jan Toenshoff, Martin Ritzert, Hinrikus Wolf, and Martin Grohe. Graph neural networks for maximum constraint satisfaction. Frontiers in Artificial Intelligncee 3:580607, 2021.
[51] Oriol Vinyals, Meire Fortunato, and Navdeep Jaitly. Pointer networks. Advances in Neural Information Processing Systems, 28, 2015.
[52] Joshua T. Vogelstein, James M. Conroy, Vince Lyzinski, Louis J. Podrazik, Steven G. Kratzer, Eric T. Harley, Carey E. Fishkind, R. Jacob Vogelstein, Cencheng Shen, Donniell E. Fishkind, Youngser Park, Michael S. Vogelstein, and Carey Priebe. Fast approximate quadratie program- ming for graph matching. PLoS ONE, 10(4)=01221022, 2015.
[53] Jan Vondrák. Optimal approximation for the submodular welfare problem in the value oracle model. In Proceedings of the fortieth annual ACM symposium on Theory of computing, pages 67-74, 2008
[54] Haoyu Peter Wang, Nan Wu, Hang Yang, Cong Hao, and Pan Li. Unsupervised learning for combinatorial optimization with principled objective relaxation, Advances in Neural Information Processing Systems, 35:31444 113314443114555 2022.
[55] Runzhong Wang. Zhigang Hua, Gan Liu, Jiayi Zhang, Junchi Yan, Feng Qi, Shuang Yang, Jun Zhou, and Xiaokang Yang. A bi-level framework for learning to solve combinatorial optimization on graphs. Advances in Neural Information Processing Systems, 34:21453-21466. 2021.
[56] Huizhen Zhang, Cesar Beltrán-Roye and Liang Ma. Solving the quadratic assignment problem by means general purpose mixed integer linear programming solvers. Technical Report 2622, Optimization Online, 2010. Technical report, April 19, 2010.
Table of Contents
Contents
Introduction
Our work
1.2
Related work
2
Birkhoff Extension
2.1
Preliminaries
2.2
continuous and a.e. differentiable Birkhoff decomposition
2.3
Properties of Birkhoff' extension
2.4
Optimization procedure with dynamic score
3
Experiments
Concluding remarks
10
Deferred Proofs
15
Unsupervised Neural optimizer
18
Extensions for Optimization over Trees
22
13

---

D Problem Details
22
D.I Traveling Salesperson Problem
23
D.2 DEASP
23
D.3 QAP
23
Experiments Details
23
Optimization Curves
23
E.2 Timing
24
Ablation Study
25
Effect of Truncation Depth k
25
F2
Effect of Score-Update Frequency
25
................................
F3
Initialization Strategy
................................
26
Optimization Method Comparison
26
Scalability on Large-Scale QAP Instances
26
Escaping Local Minima via Score Matrix Modification
26
H Comparison with Convex Closure
28
H.T Definition and Properties of Convex Closure
28
H.2 Computational Challenges and Advantages Our Approach
................................
Generalization to Constrained Permutations)
29
I.I Formal Setup
29
1.2 Key Properties
30
1.3 Compatible Constraint Types
30
Truncated Birkhoff Extension
31
Global Lipschitz bound for Birkhoff extension
32
GPU friendly maximum perfect matching
33
14

---

A Deferred Proofs
proof of Thm. 2.4 We show continuity (i) by induction. Note that Eq. 3 makes reference to the order { which we now induct on. For the base case, we have that æn is a Lipschitz continuous function of A as the sum in Eq. 3 disappears. Now, assume am for m < E is a Lipschitz continuous function of A. Then, by Eq. ar is Lipschitz continuous function of A, as min is Lipschitz continuous function.
Now for each l E [n|]], B = Em Om Pm- Clearly, have B B 0r Pe and ag min., { [ (i. (c( = 1}. Furthermore let Z' denote the set of indices (i,5) of non-zero entries in B"; that is, Z" = [(i.j) (i. 0}. First, note that since A and every P... are doubly stochastic, Be must be proportional to some doubly stochastic matrix for any l. Indeed, for any & the matrix Be has rows and columns that sum - Σm= amm- Next, we use induction to show that (a) Be has non-negative entries, (b) ag O or ag > 0, and (c) if oz > O, then 7 is strict subset of i... Z' C 2--1 Note that (c) implies that the number of non-zero entries in Bz is strictly smaller than that in B² whenever ag O. We then use these properties to prove statements (ii) - (iv) in the theorem.
First, it is easy to that properties (a) (c) hold for the base case 8 = 1. Now consider { > 1, and assume they hold for any r < l. Since all entries in BC-1 A -t=1 Om i are non-negative. namely Bi, A(i,j) Eme Om P.. (i, 2 O for any i,j E [m]. we have oz 2 0. proving property (b).
If ag = 0, then properties (a) and (c) trivially hold so assume ay > 0. Let , be indices that give rise to ay.. That is, 0g 1 (i", (i")))) and Pe("",j*) = 1. By definition of ag, for any other i, such that (r
1. we have that or S (i,)) Since P. is binary matrix, this means that for all i.j E [m]. (ePe S (i.(i..). Hence all entries in B = B -aPt are non-negative, proving property (a).
Furthermore, note that by construction, B((r, BC- - ar Pi(i", = 0. while B (i*) = a > 0. As each entry in Be is less than the corresponding entry in Bc-1, but are still non-negative, we can conclude Z' zt-l, In particular, (i", (²' E z but not in Z This proves property (c).
Hence by induction, properties (a), (b) and (c) hold for all l E [n!]. We now ready to prove statements (ii)- - (iv).
To prove (iii) we aim to show that D = O. We proceed by contradiction so suppose Dz! * 0. Then Bd is proportional to a doubly stochastic matrix and must have at least one matching, say Pe., by Hall's Marriage theorem. However, as all entries in B are non-decreasing as l increases, the set of non-zero entries in B -1 is super set of those in Ba! and, thus, super set of those in P:-. Therefore we must have Org- > O. Furthermore, let (i* ;j*)) be the pair of indices giving rise to ai-,, i.e., ag-- = B -f- (i) and Pa(i",j*) 1. Then Be (""j")= = Be" -1--1pp "i"","") - ax-
(i*))" = O. As all entries can only decrease, Bi' .j") O. However, since Pe- (i" = 1, this means that Pt- cannot be a matching for Bnd which is contradiction. Hence our assumption that Bdl is non-zero cannot be true, and we must have that Bdd O, meaning A Σn! atPt. This proves statement (iii).
To prove (ii), note that by (b), all of 0. Furthermore, as we mentioned earlier, for any E the matrix Be has rows and columns summing to E Om Since Bdd is the zero-matrix, follows that 1 Emi 4m = 0. proving (ii).
Finally, we bound the number of non-zero coefficients used in the decomposition and prove statement (iv). Consider an e" for which there are n² - n + non-zero coefficients O'm with m A 8* If no such t" exists we are done. If such ex exists, then B = A -A D 0'm Pm is either the zero matrix or is proportional to some doubly stochastic matrix. If B is the zero matrix, then are done. So assume BE proportional to some doubly stochastic matrix. Note that any doubly stochastic matrix necessarily has at least n non-zero entries (as each row-sum needs to be 1). Hence Be has at least non-zero entries. However, as the set of indices for non-zero entries zmh strictly decreases each time am > 0 (property (c) proved above), and there are n² n + 1 such am with < l", this means that
15

---

the number of non-zero entries in the original matrix A is at least n² - -n+1+n= - which is not possible (as there are only n² entries in n x n matrix), Hence B must be a zero matrix, and there cannot be more than n² statement (iv).
non-zero coefficients {arlé=1- This completes the proof of
proof of Thm. 2.3
D, lies in (n - 1) (n-1) dimensional affine subspace of RT defined by all n x n matrices with rows and columns summing to 1. We can parameterize this subspace with the linear map
Φ : R
R
(6)
a1.1
a1n-1
j=
a1
01,1
:
4
0n-1,1
0n-1n-1
0,1
1n-1,n-1
-ーー 04,1
1- -E (1 (1= and)
Note that X = (D,, is the subset of R(n-1)x(x-1) containing matrices with rows and column sums in [0. This is a closed set whose boundary has measure
since it has dimension less than (n - 1) x (n - 1). The Lipschitz continuity of Φ and ae gives that ag o Φ is Lipschitz continuous. Applying Rademacher's theorem [43 to the interior int(X) of yields that ay o is an a.e. differentiable functiion from ix) + R. Since X \ int(X) ) has measure zero, we have that ag & is an a.e. differentiable function on X. Let N C X be the of points for which ag o not differentiable. Then <(N) is the set of points such that ae is not differentiable. Furthermore, (NN) has Lebesgue measure zero in D,, since N has measure zero and ø is surjective linear map to D,.
Consider two distinct permutations P, E P Let I and I' contain the pairs (i.j) of indices such that P(i.j) 1 and pi( respectively. The score of these permutations (P. = E(.eer S(i,j) and (P', S) = E((,)r S(i,j), are equal if and only if
proof of Thm. 2.8
Σ S(i
Σ S(i,)).
(8)
(AAE
The left-hand side and the right-hand side of this equation independently distributed absolutely continuous random variables, so they are equal with probability zero. Since any pair of permutations have different scores almost surely, all permutations have different scores almost surely by the union bound.
proof of Property | Lipschitz continuity and a.c. diffcrentiability of F' follow from the Lipschitz continuity and a.e. differentiability of { (theorems 2.4 and 2.5
Now we show that the gradient V F is equivalent to
VAF(A)= E( k=1
We start with a dimensionality argument showing there are never more than n² + 2 terms in a Birkhoff decomposition. First we show the set MM of permutations with positive coefficients is linearly independent.
Suppose some P3 (PI...., PM3 can be written as linear combination Pt = EkE(M](O) C&PR for E R. We first show all with & < must be zero. The proof is by induction. We have e = 0 as there must be some i.j E [m] such that P (i,
but Ps (i, = for > Now 0. Then es 0 as there is a a E [n]
consider some A < l and suppose but for > k, ji)
Ck 0. We can then conclude e _... Gz-
such that P)(i, Finally, since there is i [r]] such that Pt(i.)) 1 but (i 0 for k > l we can conclude
16

---

Pt = EMEMMMe C&P = Ek=4 c&P cannot hold. We have, thus, shown that the set of permutations with positive coefficients is linearly independent.
Since this set of permutations is linearly independent, the maximum number of permutations with positive coefficient is one more than the dimension of D, which is n² - 2n + 2. Suppose that has full Birkhoff decomposition, that there nº 2n positive terms in its Birkhoff decomposition. Then, by continuity of {az Jo there is an open ball in Dn containing A such that all decompositions of points the ball have the same positive coefficients. Thus, the gradient the zero coefficients at A is zero.
It now remains to show that a.c. A D, has full Birkhoff decomposition. For each subset E C P, with [F = n² 2n + 1, note that the convex hull Conv(E) (E) of F has dimension n² 2n and, therefore, has measure zero in the space Dn, which has dimension n² 2n 1. Consider the union of such convex hulls
Conv
(10)
ECP. 2n4+
The space also has measure zero. Suppose A does not have full Birkhoff decomposition. Then E as can represented the convex combination of at most n² 2n permutation matrices. We can then conclude that a.e. A E D, has full Birkhoff decomposition.
proof of Property 3 (1.) Note that for E Pn we have f(P) F(P) since all Birkhoff co- efficients of P have only one term. Thus, minpeP, f(P) 2 minAED, F(A). Now suppose minpep,, f(P) minAED, F(A) and let A E aeg, F(A). Note F(A) is a convex combination F(A) EE (P)(PP) and E 0g = 1. We can then conclude that since minpep, f(P) EM Ek=1 (Ps) there must be some Pa such that f(Ps) S minpep, f(P),
a contradiction.
(2.) Suppose A minimizes ((A) over DR. Then F(A) minpep, S(P), which occurs only if for each Ps the convex combination F(A)= = Ei= axf(Pk) we have f(P)) = minper, f(P).
Since A EN k== Ork: Pk and Pk E pep, the claim argmin AED. F(A) S
Conv PeP, (()) holds.
proof of Property 4 (1) frounds(A)). F(A) then xf F(A) so for each Pk in the decomposition F(A) = Σx=1 f (Pk) we have f(Px) > F(A). However, since D ag 1 this implies Σ az (Ppk F(A), contradiction.
(2) Let the decomposition of F be F(A) = E (PP).. Recall by Thm. 2.7 that Alg. pro- duces this decomposition. Through Alg. 2
M
we have P1 argmaxrep, (P,S). Next we show argmaxpep, (P,S) p* Since each entry of S is within 1/2n of p", have
(S,
(11)
1
n1
(12)
2
Also, any o P* must differ from Pa by at least one entry so (P', Pd) S (n 1) and the inner product with S is
(S, < (n 1)
(13)
2n
(14)
where the second term accounts for the entry-wise differences between and P*. We have then shown P; argmaXpep, (P,S) = P* which implies (rounds I f(P)).
D
17

---

B. proof of Thm. 2.7 which returns only the non-zero Birkhoff coefficients. This algorithm is the same as Alg. except at each step, with B being the matrix to be decomposed at this step, we subtract off the permutation of maximum score that a matching of B, as opposed to an arbitrary matching of Recall that B is proportional to a doubly stochastic matrix, so it either has matching or is the zero matrix (in which case the algorithm terminates). To compute argmax (BP(B) (P. S) we first construct bipartite graph G that has an edge from vertex to vertex j with weight S j) if and only if ((i.j) > 0. It is easy to that (i) matchings of G correspond to the matchings of the scaled doubly stochastic matrix B: and (ii) for any matching P in this graph, its weight is exactly (P.S) which is the score of permutation P. Hence can compute argmax (gmaapeP(B)) (BP( (P, simply by computing the maximum-wett weight matching of this bipartite graph G. This computation takes (nn)) time using the Hungarian algorithm, and since there are at most (1²)) matchings to compute, the total time complexity is ( ).
Our algorithm for constructing the continuous Birkhoff decomposition is Alg.
We show the correctness of Alg.
2 Recall our algorithm returns a collection of permutations
PI..... PM- First, let Ω = [=] denote the total ordering of all permutations induced by the score matrix S. Now let A = Ei=1 âePr denote the Birkhoff decomposition of A w.r.t. the total order Ω as defined in Def. 2.3 Note the slight change in notation so that å and Pa represent the permutations and coefficients in Def.
2.33.2.3amm and or and Ps represent the permutations and coefficients returned by Alg. [ Let 1 < i2 ... < iR denote the set of indices whose corresponding coefficients Oi₃ are positive. That is, ignoring all zero coefficients in the decomposition, we have R Ex=1 die P,, where each 0sg > O. Our goal is for each i E [M] to show (cond-A): that 0g = ai and Px = P.z. We do so via induction on the index E [M].
A
We begin by showing the property that Pe is a matching of B6-1 = A Eee Fm Pm if and only if
0e > 0. Call this property (+). This property holds since ãg > 0 if and only if every element in the minimum defining āi in Def. 2.3 is non-zero, which means for each i E [n]] with (i.) = 1 we have ( > O, i.e., Pe is a matching of BC-1
The first permutation P1 returned by Alg. 2 is argmaxpeP(A)) (P,S) which is the matching of A with largest score. That is, Pi is the matching of A with smallest index in the total order f2. Furthermore, by (*). the permutation P., must be a matching of 1 A since ai, > O. For j < i we have ā, O, and again by (+), each P3 not matching of B A. We have shown P.. is the matching of A with smallest index in Z and, therefore, (cond-A) holds for the base case k = 1.
i2 Now assume (cond-A) A) holds for all m < k; aim to show that it holds for k. Let
... < i be the indices for the previous k- - non-zero coefficients. In the kth iteration of Alg. 2
B = A - Σx=1 O. Pm Bix 1, and Ps argmaXpep(n)) (P. S), which is the first permutation matrix in the total order 2 that is matching for B. Let j be the index of P& in the total order Ω that
is, Pk- = P,. Notice that j & {i....., } as this would contradict P%: being a matching of B. We claim that j > ik-1. If not, then by (+), āj > 0, a contradiction to our inductive hypothesis that i1y.... ..-.....--111 are the first k indices whose coefficients are non-zero in the Birkhoff decomposition of A. Hence, P3 is the first permutation in the list P. 1h.-tt,..... Pu that is matching of B.
On the other hand, since in is the first index in i 11 - - ..... nd such that ai, > 0, by (x), the index i is the first in this list such that the corresponding permutation is a matching of 136- We can then conclude is = and P., Px. Furthermore, Ois or since B B Bi- so the definition of og in Alg. og = ming { (i,, j) (i.)..) = 1} is equivalent to the definition of ai, in Eq. 3 This finishes the proof of the inductive step. Combining the base case with the inductive step, we have that (cond-A) holds for all & E [M]. hence the set {(a,, Pb) returned by Alg. 2 exactly corresponds to those terms in the Birkhoff decomposition (as computed by Def. 2.3 with non-zero coefficients.
B Unsupervised Neural optimizer
Since Birkhoff extensions are a.e. differentiable, they can be useful for training neural networks for
unsupervised neural combinatorial optimization. particular, similar to [29] we can train a neural network Ne with parameters 0 that maps an instance I of a problem to a doubly stochastic matrix A,,
18

---

which we aim to train for the optimization of the combinatorial objective f. See the illustration in Figure
For example, for TSP in the Fuclidean space Rd the instance is a set of cities, represented by a vector X, IRnd representing n points (21,....; ; } in Rd, The output of No is n doubly stochastic matrix A = (X(r,), and the neural network is trained in an unsupervised manner to optimize F(Ar). Note that once trained, when a new instance I' is given with input r,, we can simply return rounds ( (Xr))) as the TSP tour. Essentially, No can be viewed as an neural optimizer for the given optimization problem over the extended space Dn,- Once solution in D,, identified, it can be rounded to permutation without lowering the quality of the solution (Property (
Having a differentiable Birkhoff extension allows us to train such neural network model in an unsupervised manner. In particular, first, suppose we have a score matrix 5 - this score matrix can be chosen simply as random stochastic matrix; it can also be a canonical choice depending on the input problem instance. For example, in the case of TSP, we can choose 5 to be a perturbation of the permutation derived from the MST. With this choice of S, let A (Nrr be the output of the neural network. We have that (computed by Alg. 2
M
Fs(A) = Σ (A))((()))))) 1=1
(15)
Here, note that both and Ps depend on A, and A itself depends on the parameters @ of the neural network Ng. We simply minimize Fs (A) w.r.t. the parameters via backpropagation. More precisely,
computing @FelA) boils down to computing day Bog 0A for each positive Birkhoff coefficient a'k-
06
a0 = aλ 30
The above description is for training No only for a single instance. Usually one wishes to train No over family of instances I, that once trained, it can be used to produce solutions to new instances. In particular, during training, the loss is Erez 1
[ΣT (Nx(( Once trained, given a new instance I, we can simply compute A (Xr) and return the permutation rounds(AA) as a candidate solution. In practice we found that for a test instance it makes sense to optimize No for a few more iterations testing (as a fine-tuning) to further improve the quality A₁.
We also note that for the case where we have a score matrix that depends on the problem instance (e.g, using MST to induce score matrix for the TSP problem), it is beneficial to also take this score matrix Sy as input to the neural network Ng to better inform the output matrix A, = No St). This input is optional.
Finally, in the simplest form, No: Rnd + Rax= maps an n-vector X with each entry from Rd to a n X n matrix A. In particular, again using TSP as an example, here we represent a set of points as a vector X Rnd which assumes an ordering of these points. The output matrix assumes the same ordering of input points. In other words, this map No needs to be permutation equivariant, namely, if we permute the order of input points in X, then the output should permute in the same way. Mathematically, this means that No satisfies over n elements. Such a permutation equivariant neural network can be implemented using models such as the set2graph neural network of (47) and the equivariant-graph network of (35]
= PN₂ (NN())) PT for any permutation matrix P
We train a neural network No for solving TSP using a loss comprised of the Birkhoff extension F and penalty term X Σ. 4E,
+ Σ, С C, A.g -1)² that ensures double stochasticity.
where À is the weight. We use the Adam optimizer [ to train the model. We generate a training dataset with N 6000 instances and with mixture of instance sizes, n = 20, 30, and 40. The input to No. for each instance, is the vector X2 and the MST-derived score matrix S.. We trained N₂ for T= 100 epochs, and selected learning rate of ' 0.001 using a hyperparameter search of the set (0.01, 0.05, 0.001, 0.0005). >. At testing we optimize the trained neural network No for a few more iterations a fine-tuning. We compare the performance of this fine-tuned model (labeled MST 5 Init. NN w. Alg. 2 with the corresponding untrained model (labeled MST S Init. Alg. 4) in Table 4. We show that most of the cases, the trained neural network model can achieve similar solution quality to the untrained model with much less runtime. This result holds even for problem instances that are larger than the instances seen in training
19

---

Time
< Is
967s
1318s
24.56s
5.42s
n = 100
Obj-
10.288
13.560
9.836
10.239
10.220
Table 4: Performances trained neural network pure optimization both with Birkhoff extension on TSP.
Time
671s
640s
9.72s
3.02s
11 == 50
<
7.410
8.433
6.892
7.333
7.170
Obj
TSP: Tour Length J
Time
Is
612s
571s
8.29s
2.59s
11 == 40
Obj.
6.835
7.139
6.328
6.529
6.219
Time
<s
259%
214s
5.30s
2.15s
n= 30
Obj.
5.784
3.681
5.290
5.512
5.206
Time
1s
192s
156s
2.71s
1.75s
n 20
Obj.
4.746
4.518
4.345
4.541
4.180
MST S Init. NN W. Alg. 2 (Fine-tumed)
Random Init. Alg. MST 5 Init. Alg. 4 (Improv.)
MST S Init. Alg. (100 steps)
Method
MST
20

---

Penally = 1) + Σ Σ 1/- 119)
Σ Σ
Au
X
Matrix A fapproximately doubly stochastic)
Neural network N₁
Lose: (s(A) + penalty
First Few Birkhoff Extension Terms
Fs((A) Σ ((AA(P(()))
Score Matrix 5
Gradient
Forward Flow
Figure The pipeline of training neural network No for a single instance. For a given problem instance J, we have its representation X, and a score matrix S₁ as input to the neural network Birkhoff
No. The output of the neural network a doubly stochastic matrix A extensions are used to compute the loss (A) E (PP) and we minimize via backpropagation. In the figure above, for example. M and rounding produces the permutation P2 = rounds, (Ar), highlighted in red.
21

---

C Extensions for Optimization over Trees
In this section consider how a parallel framework can be applied to optimization over rooted binary trees with n labeled leaves. Problems such as the Steiner minimum tree problem and optimization- based hierarchical clustering reduce optimization over this space. For instance, given a fixed rooted binary tree with leaves the set of terminals, the optimal Steiner tree with this topology can be computed efficiently. Therefore, computing the Steiner minimum tree reduces optimizing the topology.
We begin by representing the space of rooted binary trees matrices. Definition C.1. Let Tn, be the space of rooted binary trees with leaf set [ni], which we represent by a directed graph with edges directed away from the root. Let W2n-2 Dan-2 2 be the space of (2n - 2) x (2n- 2) doubly stochastic matrices W that satisfy Wii, = 0 if either
1. i > and i S j. or
2. i S and i + (n j.
Let B2n-2 C W2n-2 2 be matrices in Wzn that have binary entries. We relate these spaces by the map
B2n-2 + Tn
(16)
B B+
(17)
where T the tree with vertices [2n- - 11, leaves [n], and for n the free T has an edge (1,5) iff B(i,j) = 1 or B(i - (n 1),5) =
<is S 2n and 2n
Lemma C.2. T well-defined and surjective.
Proof. First, we claim that for any B E B2n-2 the image 7 rB) ) is indeed tree in Tn Note that since the columns of B sum to I, each vertex of TB)) other than the root, which is vertex 2n 1, has in-degree 1. Similarly, since the row sums of B are each of the internal vertices (non-leaf' vertices), which are the vertices (n +..... 2n - 1}. has out-degree 2. (For each internal vertex, there two rows in B that give its children.) Furthermore, r(B) has edges (i,j) where i S j as the only entries in B that can give rise to such edges are zero by (1.) and (2.) Def. C.I We can then conclude that 7 r(B) ) is directed acyclic graph. Furthermore, rBB) has cycles (regardless of direction) since for such a cycle to not be a directed cycle must contain a vertex with in-degree greater than 1.
For surjectivity, consider an arbitrary E Tn- Let Φ {n + 1,..., 2n 13 (n 1,..., be an enumeration of the internal vertices of T' that respects topological order, i... Φ(i) > o(j) implies j is not a descendent and (i +
For each internal vertex i in T with children j. let B(i,j) = 1, = 1. Since respects the topological order, we guarantee there are no entries B 1) > with either < and i > n 1 or + (n 1) < j and i S n 1, so (1.) and (2.) of Def. C.I are satisfied and B B2n-2. This B satisfies 7BB) T. thus, is surjective.
Matrices the space W2n-2 can be decomposed using Birkhoff decomposition as W2n-2 C D2n-2- Additionally, if W E W2, then its Birkhoff decomposition only contains permutations in B C P2n-2. By this fact, we are then free to apply Birkhoff extension to extend any function on Ban-2 to function F: W R. even if is not defined for P B2n-1-
Additionally, we can extend functions on trees Tn to functions on W2n 2; the procedure as follows. First f is composed with 7 to yield a function for B2n-2 + R. Then Birkhoff extension is used to extend this function to F W - RR. If the extension W E W2n-2 the Birkhoff extension rounding scheme can be used to find B E B2n-2 such that f(r(B) < F(W). ). Here, T TB)) is a tree T satisfying f(T) < F(W), , so this procedure can be used to optimize f.
optimized to find some solution
D Problem Details
In this section we give more detailed problem definitions and show the integer linear programs we employ for optimization using Gurobi.
22

---

D.1 Traveling Salesperson Problem
Definition D.I. Given a set of cities {1, (1,2,...,, and a distance di,j between each puir of cities and j. the traveling salesperson problem (TSP) is to find the permutation [mx] -> [n]] that minimizes
(k), (mod n))
(18)
TSP can be formulated as an integer linear program using subtour elimination constraints such as in the Miller-Tucker- Beeker- formulation [36]. For optimization with Gurobi we use the formulation given in 26
Our experiments use the integer quadratic program minpep, Ei.1 " ("ODP(+)+PP"0DDPP
and its relaxation minpeD, ΣE= pT (i) + + Pd (P) for TSP. Here P(i) denotes
the i-th column the square matrix P.
D.2 DFASP
be directed graph where V the set of vertices and E. the set of directed edges. A feedback are set in G is subset of edges F C E. such that the subgraph G' = F) is acyclic. The Directed Feedback Are Set Problem (DFASP) to find the minimum feedback are set, i.e., the smallest subset of edges F' whose removal makes the graph G acyclic.
Definition D.2. Let G =
(V, =
Alternatively, DFASP can formulated as a vertex ordering (i.e. permutation) problem, where the goal is find an enumeration of vertices Jn that minimizes the cardinality of the set of backward edges. Here we give an LP formulation. Let Zig be a binary variable that equals I the edge (i,j) E
backward edge, and otherwise. Let yi be an integer variable for each vertex i V representing the position of vertex in topological ordering. DFASP can be formulated using the following integer linear program of [? J.
minimize
Σ Tij
(K))E
subject to
Tij
{0,1}
v(i.)) E E,
Vi
ViEV, V, v(i.))
+ 1 S yj + | Taj
D.3 QAP
Definition D.3. Given two non-negative n. X n matrices D and L find a permutation x [n] + [n] that minimizes
Σ
((.))
(19)
We the Nuufman-Broecke integer linear program formulation 30 for the Gurobi implementation.
E Experiments Details
In this section we provide additional experimental results for the TSP, DFASP and QAP.
E.1 Optimization Curves
Here, show the optimization curves for each experiment.
Below are plots of the averaged optimization curves for TSP problem instances of different scales.
The r--axi gives the number of steps / and the y-axis gives solution tour length, averaged over the N 50 instances.
23

---

n-20, N-50, TSP
D-30, N-50, TSP
n-40, N-50. TSP
Fandom 5 IIE. Alg.
Random Ale.
Bandom is Ini. Alg. Alg.
HET Int, 4
32
MIT Aug
31
GP Init Alg. 4 MS1
CP S lnit Alp. MB1T
MST
20-
Gurebi
Tour Loogehh
Gurcbi
Tour
0.115
3,002
5.15
2000
4000
5030
5000
10000
2000
4000
5000
8000
10000
2000
6000
5000
80000
10000
Stopa
Stops
Steps
n-
5-50,
n=100, No:0. TSP
22.5
Hardkam Alg.
MST Alg. 4
20.0
Alg.
oP Int. Alc. 4
17.5
75
HST
Ou obi
12.5
10.0
7:11
7.248
5.
2000
4200
60000
0020
10200
2000
4000
6000
5000
10000
Stegs
Stops
Below are plots of the optimization curves for DFASP. The a-axis gives the number of steps I and the y-axis gives the average cardinality of the feedback are set returned by the optimization algorithm. The runtime limit for both methods is a minutes (c.g. 2 minutes for n = 20, 5 minutes for n = 50 and 10 minutes for
TO
100). Results are averaged across N = 50 instances.
p=0.1. 5- Directed Feedbark Are Set Problem
n= p=0.5, N--00 Directed Feedbackk Arc Set Problem
n=20, p=0.5. N=500 instances. Directed Feedbrack Arc Set Problem
of boduaee Set
Rarcon Aly
Karkens kil. 4..4
170
Kaxion Irl, Aly
t0
100
161
100
Stepe
Steps
Stepe
ni p0001. NO instances,
n-:50, pii... N--50 Instances.
n=:5. p=0.5. N-:50 Instances.
Directed Feedback Arc Set Problem
Directed Feedbackk Arc Set Problem
Directed Feedbrack Arc Set Preblem
100
Karsorn Aly
600
Korkens kil. 4
Aly
Candrality of fwwtterrr.
1090
200
560
1000
10/0
1090
500
1000
(E0
200
350
Stepe
Steps
Steps
1=000 pm001. N-50 instances.
n=000, pmo0.s. N--50 instances.
n=100, p=O.S. N-55 instances.
Directed Feedback Arc Set Problem
Directed feedbback Arc Set Problem
Directed FeedBrack Arc Set Preblem
475
Kencom Aly
Nondoms Iill 4
4150
4425
24000
4400
400-
275-
2200
2250
105.
42000
375.
4250
420
1003
403
Sipps
Stps
Steps
E.2 Timing
We are using CentOS Linux 7 system on a Intel Xeon CPU E5-2650L processor. Gurobi is run on a single thread per process. We report average runtime for each method for different scales of TSP and QAP problems, see Table and Table ?? below. DFASP not included as we fix runtime for this problem
24

---

TSP
Method
n 20
30
40
n 50
100
Gurobi
< Is
< Is
< 1s
< 1s
< 1s
MST
< 1s
< 1s
< 1s
< 1s
< 1s
QP
36s
50s
69s
95s
240s
Random S Init. Alg.
192s
239s
612s
671s
967s
MST S Init. Alg. 4
156s
2114s
571s
640s
1318s
QP S Init. Alg.
135s
249s
532s
615s
953s
Table 5: Average runtime of TSP algorithms in seconds.
Ablation Study
This section provides additional ablation studies and extended experiments to support the main results. We examine the effects of truncation depth k, score-update frequency m, initialization strategy. optimization method, and instance size. Each subsection specifies the experimental setup used to isolate the respective factor.
F1
Effect of Truncation Depth
To analyze the influence of the truncation depth k in the Birkhoff decomposition, conduct controlled experiments on a subset 14 QAP instances selected from QAPLIB, with problem sizes ranging from 12 to
35. All experiments use identical hyperparameters: learning rate = 0.001, total runtime t = 2n seconds, and optimization via the Frank Wolfe (FW) algorithm with random initialization. We vary & {3, 5. 10, 15, 20) and report the average gap from the optimal objective. Figure 2 below shows that performance improves as l: increases up to k=10, then the performance decreases. This indicates that moderate truncation depth provides good balance between efficiency and solution quality.
Gap from Optimal by Truncation Size
Gap from Optimal by Update Frequency
4.8
Gap from optimal (%)
Gap from optimal (%)
4.6
5.50
4.4
5.00
4.2
4.50
2.5
5.0
7.5
10.0
12.5
15.0
17.5
20.0
10
20
30
40
50
Truncation size (k)
Update frequency (every m epochs)
Figure 2: Gap vs. truncation size
Figure 3: Effect of score-update frequency
F.2 Effect of Score-Update Frequency
We study the impact of the score-update frequency m. which controls how often the scoring matrix is recomputed during optimization. The same 14 QAPLIB instances and identical optimization settings are used in the truncation experiments, fixing k=5. We test m (1,5, 10, 20, 50), where m=1 means the score matrix is refreshed at every iteration. Figure show that more frequent updates consistently yield smaller optimality gaps, demonstrating that dynamic score re-weighting effectively helps the optimizer escape poor local minima.
25

---

F33 Initialization Strategy
We compare random initialization and barycenter initialization understand their effects on opti- mization stability. Both configurations are evaluated on the same 14 QAPLIB instances with identical hyperparameters (k=5, 1700.001, runtime t=2n seconds). In the barycenter case, optimization starts from the uniform doubly stochastic matrix, while the random case samples each entry uniformly from [0, 1] before projecting onto the Birkhoff polytope via Sinkhorn iterations. Both use the same Frank Wolfe optimizer. The comparison shows that random initialization achieves slightly smaller average gap to optimal 14% vs. 4.42%), see results below Table suggesting that stochasticity in initialization provides beneficial diversity in the optimization trajectory.
Initialization Method
Avg. Gap from Optimal (%) 1
Random Initialization
4.14
BaryCenter Initialization
4.42
Table 6: Comparison of initialization strategies on 14 QAPLIB instances (n = 12-35). Both methods
use identical hyperparameters (k 5, r = 0.001, runtime 2n s). Random initialization achieves slightly lower average gap from optimal than barycenter initialization.
F,4 Optimization Method Comparison
We compare two continuous optimization schemes over the Birkhoff polytope: (1) gradient descent with Sinkhorn projection and (2) the Frank- Wolfe (FW) algorithm proposed in this work. Both methods operate on the same 14 QAPLIB instances, with fixed k=5, learning rate 1z=00001, and runtime t=2n seconds. For the gradient- based variant, we apply 10 Sinkhorn iterations after each gradient step to enforce double stochasticity, FW achieves a lower average optimality gap (4%1 than Gradient+Sinkhom (4.55%), see results in Table 7 below.
Optimization Method
Avg. Gap from Optimal (
Gradient + Sinkhom Projection
4.55
Frank-Wolfe (ours)
4.14
Table 7: Comparison between Frank-Woll Wolfe and Gradient + Sinkhorn optimization over the Birkhoff polytope 14 QAPLIB instances (sizes n= Frank-Wolfe achieves lower average optimality gap. All experiments use identical hyperparameters (k=5, 1000001, runtime 2n s).
F.5 Scalability Large-Scale QAP Instances
To assess scalability, we perform large-scale experiments on 20 synthetic random QAP instances of size 1=1000. Because QAPLIB includes no instances larger than n=256, we generate our dataset as follows: we uniformly sample 1000 facility points and 1000 location points in [0, 1 The distance matrix D is constructed from pairwise Euclidean distances among locations, and the flow matrix F is defined as the inverse of pairwise distances among facilities. We evaluate BE using both Frank-Wolfe and Gradient-Sinkhom optimizers under random and FAQ initialization. Each configu- ration is run with runtime limits of 2000 and 4000 seconds, corresponding to approximately 10, 10,000 and 20, 20,000 optimization epochs, respectively. We also include comparisons with Fast Approximate QAP (FAQ) and 2-opt baselines, see the results below in Table
G Escaping Local Minima via Score Matrix Modification
In this section, we provide theoretical justification for Algorithm 4's strategy of modifying score matrices escape local minima. We first establish key lemma about the linearity of our Birkhoff extension along certain directions, then show that when the current solution is local (but not global) minimum, there always exists a score matrix that allows escape from this local minimum.
Lemma G.1 (Linearity Along Highest-Score Direction). Let f : Pr. - R be any function on permutations, and let Fs be the Birkhoff extension of f induced by score matrix S. Let P be permutation matrix that has the highest score under S. Then for any A E D, and B E [0, 1]:
Fs((1 B)A + BP) (1 (1--3)E5 (BFFFss()) + 3F()
26

---

Method
Run-time (secs)
of epochs
Avg. Objective :
Gap from best %
Random assignment
138,005.401
77.37%
2-opt (Restart=3)
2,000
119,216,266
53. 1.3%
2-opt (Restart=5)
2,000
109,965,065
41.28%
EAQ
345
77.894,862
0.10%
Ours (BE Grud+Sink, Random Init.)
2,000
10.133
100,894.88 812
29.64%
(BE Grad+Sink, Random Init.)
4,000
20,064
90,244,871
15.97%
Ours (BE w/ FW, Random Init.)
2,000
10,024
102,393,212
31.57%
Ours (BE w/ FW, Random Init.)
4,000
20,091
89. 15,
14.40%
Ours (BE w/ FW, FAQ Init.)
2,000
10,102
77,819,625
0.00%
Table 8: Comparison of QAP solution methods on large-scale instances.
Proof. Let A =K== agPk be the Birkhoff decomposition of A under score matrix S. Consider (1- 3) -PP
We consider two cases:
Case 1: P + Pi (i.e., P is not the first permutation in the decomposition of A), then the Birkhoff
decomposition of (1 -3) + 3P is simply 3P (1 3) E agP This follows because:
appears first with coefficient since it has highest score
The remaining terms proportional (1 (--)) 33 times the original decomposition
Therefore in this case:
Fs((1 3P) (F(P) (1-3))ssA
Case 2: If = Pi (i.e., P is the first permutation in the decomposition of A), then for any entry
(i.j) where P(i,j) 1:
M
((1 A (() = (a +
Σ 2=
a Px(l.))) +
M
(1 - 3) 3 + (1 3)
Σ on (i
2=
By definition of the Birkhoff decomposition, there exists at least one entry (i", where P(i".j")
1 such that Σk=2 (i" j²t) 0 (this is what determines on in the original decomposition). For this entry:
((1 - 3)4 (T) ("))")))"+)))"" B)an +
This means that in the Birkhoff decomposition of (1 -3)A + BP:
P must appear first (since it has highest score) with coefficient exactly (1 1B))xx + + (since this the minimum value at any position where P has 1)
After subtracting ((1 - 3)a + P, the remaining matrix is proportional (1 - 1(9) times the matrix obtained after removing P from A's decomposition
27

---

Therefore:
M
()(1--)AA-00 (1- (1- 3)
Σ
an
((P)
k=2
((1 -(1- (L=)))a)((P))++3())) (1- 3)
a (PP)
(1 ()(nJ(()) + Σ o
k=2
(1 (A)s(A)+BB()))
In both cases we obtain Fs((1 3)A + 3P) (1 (L))))
Proof of Thm 2.77 We choose S as S PM where P is some minimizer of f. This ensures that P has strictly higher score than any other permutation under S. For any permutation matrix o P the score of P under S is n (since each permutation has exactly 1 s)) while the score of P under equal the number of positions where P and p agree, which is strictly less than n. Therefore, in the Birkhoff decomposition induced by S, P appears as the first term.
Recall that vector is called descent direction point if moving infinitesimally in direction decreases the function value. This can be verified by checking if the directional derivative is negative. For our extension F, the directional derivative in direction p is:
(A + c(P'
A)) Fs(A)
(V (VFs(A), P A)
lim
c00
Fs((1 e) + cp (A)
lim
c-
(1-c)Fs(A) + (((P)) Ps(A)
lim
(by Lem. G.1
c-+0
e(f(P-) - -Fs(A))
lim
c+0
f(P") (A()
The third equality holds by Lem. G.1 since construction of S ensures P* appears first in any Birkhoff decomposition.
If f(P*) - (A) = 0 then A must be global minimum of Fs. If f(P)) - (A( < 0 then p A descent direction for Fs at either case, A is not local minimum of Fs.
This theorem provides theoretical justification for Algorithm 4's strategy of modifying score matrices when stuck at local minima. shows that as long as better solutions exist (i.e., permutations with lower objective value), there always exists matrix that enables escape from the current local minimum.
H Comparison with Convex Closure
In this section, compare Birkhoff extension with the convex closure, well-known theoretical concept in convex analysis. While the convex closure provides certain theoretical guarantees, our Birkhoff extension offers practical computability, making it more suitable for real-world optimization problems.
H.1 Definition and Properties Convex Closure
The convex closure of function ʃ P,, + R is defined as the greatest convex function that is pointwise less than or equal to f on the domain of permutation matrices. For our purposes, since we extending from P., to D,,, the convex closure can be expressed as:
28

---

Definition H.1 (Convex Closure). Given a function f : Pn + R, its convex closure Conv (f) Dr - defined as:
k
(o))(f)(X)= = inf Σ f(P( X Σ
a.P., 2 0, a =1
1=1
i=1
where the infunum is taken over all possible Birkhoff decompasitions of X.
The convex closure is convex by definition, agrees with f all permutation matrices, and for minimization problems, minimizing Conv(ʃ) over Dn gives the same optimal value as minimizing f over Pn-
H.2
Computational Challenges and Advantages of Our Approach
Despite its theoretical appeal, computing the convex closure for a general function f and doubly stochastic matrix X is computationally intractable. Determining the convex closure requires solving combinatorial optimization problem to find the optimal set permutation matrices and their weights, which is generally NP-hard for large
In contrast, our score-based Birkhoff extension can be computed in O(n)) time for any doubly stochastic matrix. It is continuous and almost everywhere differentiable, enabling the use of gradient- based optimization methods. While not necessarily convex, our extension preserves the property that minimizing over Dn and rounding gives the same value minimizing directly over Pr (Property 3). The flexibility provided by different score matrices allows for exploration of the solution space and escape from local minima.
The computational efficiency of our Birkhoff extension makes it applicable to practical problems where the convex closure would be intractable to compute. For optimization problems over permuta- tions, we can use gradient-bas based methods efficiently, integrate with machine learning approaches like neural networks, and handle larger problem instances than would be possible with exact convex closure computations.
In summary, while the convex closure provides a theoretically optimal convex relaxation, our Birkhoff extension offers practical, computable alternative that preserves many desirable properties while enabling efficient optimization.
I Generalization to Constrained Permutations
The framework presented in this paper be naturally extended consider only subset of permutations that satisfy certain constraints. This generalization is particularly useful in applications where not all permutations work with any subset permutations for which we can efficiently solve the corresponding matching problem under the given constraints.
valid or desirable solutions. The key insight is that framework
I.1 Formal Setup
Let C. S. Pn a subset of permutation matrices that satisfy some constraints. We define the constrained Birkhoff polytope (C) as the convex hull of C,,:
(C, A AERaX =A=
Σ apP, Σ ap==11 app.Dap=
1,ap 2
(20)
PEC,
For this framework to be practical, we require two key properties of the constraint set C,,
Efficient Matching: We must able to efficiently solve the constrained matching problem argmaxree. (P,S) for any score matrix S. This same matching problem is used in two contexts our framework: (1) finding the highest scoring permutation in the decomposition algorithm, and (2) finding the permutation with greatest inner product with the gradient in the Frank- Wolfe optimization step. The matching problem must be modified respect the given constraints while still finding the permutation that maximizes the inner product.
29

---

2. Efficient Initialization: We must be able to efficiently generate points within the con- strained polytope (C).). This is crucial because optimization methods like Frank Wolffe require feasible starting points. For many constraint types, finding such initial points is straightforward.
The continuous decomposition framework can be adapted to work with D, (C.) as follows:
Algorithm 5 Constrained Continuous Birkhoff decomposition
Require: A (C), identifying score matrix $ k= M.M
Ensure: P k= s.t. A = Ek=n PkPk,
Σx ag 1, 0% 0, and Pk Cn.
+ B <-
while B + do
Ps f argmax PEC, SP,S)
0kk f minej(B(i,)) Pa(i, = 1}
B+B- B-a.Pa +++
end while
M t- k return {lak, P1)
For optimization over the constrained polytope, we use the Frank- Wolfe algorithm with constrained matching:
Algorithm 6 Constrained Frank-Wolfe over (C)
Require: Initial point Ag E (C), score matrix S, step sizes f% (7)71 Ensure: Final doubly stochastic matrix Ar E (C) An <- initializeeD,.(())
for = to
Pt f argmax (P, (Fs(A-)))
A; f (1 -
end for return Ar
1.2 Key Properties
The constrained version maintains two key properties of the original framework:
1. Continuity: The decomposition remains continuous long as the constraint set C, is fixed and the matching problem under constraints can be solved efficiently. This follows from the same argument in the unconstrained case, where we fix an ordering over the constrained permutations. The continuity of the decomposition implies that the extension F(A) E. ou f(Ps) of any function / C, ->
also continuous.
2.
Rounding: The constrained version inherits all the rounding properties proved for the origi- nal Birkhoff extension, as the the extension still a convex combination permutations.
I.3 Compatible Constraint Types
The framework can handle any constraint set C, for which the matching problem argmax pec., (P, S) can solved efficiently. A key class of constraints that satisfy this requirement are those whose linear programming relaxation has integral solutions. This means that when we relax the permutation constraints to allow for doubly stochastic matrices, the optimal solution of the constrained matching problem is still a permutation matrix. This property ensures that the matching problem can be solved efficiently using standard linear programming algorithms, and the solution will automatically satisfy the permutation constraints.
30

---

Most notably this includes LPs with unimodular constraint matrix and an integer right-hand sides 45 The integrality of the LP relaxation powerful property that guarantees the existence of efficient algorithms for solving the matching problem, making these constraints particularly well-suited for our framework.
A simple but important example is the case of "don't It match" constraints, where certain pairs of elements are forbidden from being matched. This can be expressed by setting the corresponding entries in the score matrix S to -0,, effectively making these matches impossible. The LP relaxation of this constrained matching problem has integral solutions, as the optimal solution will never select forbidden match. This basic case demonstrates how the framework can handle even simple constraints while maintaining its key properties.
Truncated Birkhoff Extension
In practice, computing the full S-induced Birkhoff decomposition can be expensive. We introduce a truncated and normalized variant that uses only the first K non-zero terms of the decomposition while preserving the conexx-combination interpretation.
Definition J.1 (Truncated, normalized Birkhoff extension). Let S be an identifying score matrix and let (or (A), Ps (A)) be the non-zero terms in the S-induced Birkhoff decomposition of A E D,,. ordered as in Alg. 2 For any K E {1,...,M, M). define the truncated normalization factor
K
(x( = Σ (A),
(21)
R=1
and the normalized coefficients
(A)
0g
(K) (A) (A)
(A)
if S K,
(22)
ifk> K.
The truncated Birkhoff extension of f:f:fo. Pn +
R at level K is
K
F5 (K)
(A)
Σ
a k=1
(K)
(A)
(23)
which satisfies E 0t
(A) (A) by construction.
Theorem J.2 (Escaping Local Minima for Truncated Extension). Let f Pn - R and let (() be the truncated, normalized Birkhoff extension at level K > There exists score matrix S such that any A E D, is not a local minimum of (K)
Proof. The proof parallels that of Thm. 2.1T Choose an identifying score matrix S and choose it so that some P & arg minp f(P) is the unigue top-scoring permutation under S. Let the S-ordered decomposition A be A = Σi21 a,P, Define the path A(B) = (1 3)A + 3P* .
[0, 1].
Let J be the set of the top-K permutations (by S) B after inserting pa if needed: if pa & (PI,..., Px) then J
PK- otherwise J = -P Px.. Write
S
Σ ap, PEJ
Ts Σ Capf(P S(P),
F^(A) T3/S,,
PEJ
Since S is identifying, the S-order is strict; hence for all sufficiently small 3 O. the top-K K set along A(3) remains and all coefficients terms in J are scaled by (1 3) while p gains +B. Therefore
(
F (K) (A(3)) (1 (1-B)S)+3
Differentiating 1 0+ gives
d
F (K) (A(3))
-(Pt)- (s(A)
< O.
d3
S
3=
31

---

since f(P*) minp f(P) S fA(A). If the inequality is strict, P - A is a descent direction and A is not local minimum. If equality holds, then fA( = f(P), so every E already attains the global value and F
(A) = f(P"), ie., A already achieves the global minimum of FS F(K), This proves the claim for K 2
(K) F(K)
Remark J.3 (Truncated preservation of properties). All properties established for acore-induced Birkhoff extensions continue to hold for the truncated, normalized extension F(k);
Efficient computation (Prop. 2) Running Alg. for K iterations yields the truncated decomposition and F{K) in O(K r²) time (via K maximum- weight matchings), which is S 2O(re) in the full case.
Minima correspondence (Prop. 3 For any permutation P, FS (P) = f(P). For any A,
(x(K)
(A) is a convex combination of (F(Pr(A)))&Ck, so F5 (A) 2 minper. S(P). Thus F
mInAED. (A) minpep, (P), and any minimizer belongs Conv(argminper. (PP)
by the same argument as Prop. 3
Rounding and approximation (Prop. 4. Using the K permutations returned by the truncated decomposition, define round's K) (A) argmin SK
(PP() Since F (A) is an average over
have fround (A)) (A)) F (K) (A). The C-approximation transfer follows identically, For the score-update guarantee (Prop. 22. if is sufficiently close to P then P is the top-scoring permutation under S' and thus appears the first term, in particular among the first K terms for any K 2 1; the same argument goes through.
(A))
K. (S.,,
Dynamic score and escaping local minima. Thm. J.2 provides the truncated analogue of Thm. 2.17 the proof follows the same linearityallog- direction combined with normalization.
K
Global Lipschitz bound for Birkhoff extension
We now provide an explicit Lipschitz bound for the (non- (non-truncated) score induued Birkhoff extension that is independent of the choice of score matrix.
Theorem K.1 (Score-matrix-indeeeeeee
Lipschitz bound). Let f Pn + IR and define the oscillation و = maxp ((P) 100--mmnn minp ((P). For any identifying score matrix and any A, A' E D,,
[Fs(A) - Fs(AA) < + A A'll
(24)
In particular, the Lipschitz constant in 24 is independent of S.
Proof. By Prop. 1 Fig is a.e. differentiable on Dn with
VAFs(A) Σ (VAae(AA) [(Pt), L4 = (e: > 0}.
EEL+
From Def. 2.3 (A) is the minimum over the entries on the support of P: after subtracting earlier terms. Hence, at points of differentiability, V V is supported on a single entry (the active minimum), and different have disjoint supports along the decomposition path. Moreover, since Ee (e(A) for all A, we have E, (A))(A) = 0 entrywise. Indeed, the map A + E ae(A) function 1, so directional derivative any direction H zero: 0= O= E Aar((A, H). As this holds for every H, it follows that the matrix Ee V AaxeAA) must be identically zero. Using this, for any scalar e,
A) is the constant
V AsFs(A) = Σ (FPP) -).
FET.
Choosing C maxe ((P)+mmnn f(P) each entry of V (AF((A) receives contribution from at most one E (disjoint supports), it follows that
centers the values that (f( < A for all É. Because
V Fs())l < A/-
Along the segment A, = (1 4)4 +) + t.A', the fundamental theorem of calculus gives
Fs(A') - -s(A) = [ (V)Fs(A), A A) dt,
32

---

where (.-) denotes the Frobenius inner product. Then (X, S Oo |||| yields
for matrices X, Y
S ! ||V dt || - A|l - _ || - A|i,,
where we have applied the uniform bound on Fsllo a.e. along the segment.
L GPU friendly maximum perfect matching
We present a suboptimal randomized algorithm for the Maximum Weight Perfect Matching problem in bipartite graphs, with a runtime of ((n²). While this problem can solved deterministically in (n²) time, and more efficiently in O(n") time using randomized algorithms with low failure probability [37 (44, our goal different: we aim design an algorithm that is simple to implement GPU. Our algorithm runs in Onn²) time and returns perfect matching with weight
and efficient
at least half the optimal.
Let G = VU, V, E') be a bipartite graph, where [0] = V n. U (nn V = {(I..... ) Let A(G),,, be random number from {1, .... R} if (u,, u.) E E, otherwise A(G),.) O. Here R is chosen to be large enough e.g.,
n (11) and the matrix G. The matrix A(G') has nice properties some of which we use here to design our algorithm. First the rank of A(G) at most equal to the size of the maximum matching and the equality holds with probability at least - n/ Rt. More importantly, if G has perfect matching then with high probability det(A(G)) # 0, and so A(G) is invertible. Rabin and Vazirani
known as random adjacency matrix of
4] showed that
Theorem L.1. With high probability, A(G) # 0 if and only if graph G - {u;, 2j} has perfect matching.
In particular for edge (u,, uy). A(G), * O if and only if e is allowed, i.c. is contained in a perfect matching. We use this result and the techniques developed by [37] to design our greedy algorithm for finding maximum weight perfect matching. We note that in the worst case algorithm returns a perfect matching with weight at least 1/2 the optimal one.
Algorithm 7 2-approximation maximum weight perfect matching
Require: Bipartite graph G
B = AGG)
M 0
for c - n. do
Find maximum weight edge (ue, u,) that is allowed edge in G M i.e. row r, not yet eliminated, and such that B, f 0 and A(G),.r of O. eliminate the r-th row and the e-th column of B
Add (ue, Ur) to M
end for return M
Theorem L.2. Algorithm 7 runs in time O(n²) and with high probability returns perfect matching M such that the weight of M is at least 1/2 times the weight of the optimal matching.
Proof. The proof that Algorithm 7 runs in time (nn3) and with high probability returns a perfect matching follows from [37]. The approximation guarantee follows from our greedy selection among the feasible edges at each iteration.
NeurIPS Paper Checklist
Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
33

---

Justification: We claim to provide an extension (Sec. 2.3 optimization algorithm (Sec. 2.4 and experiments (Sec. all of which are clearly given the paper.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or
NA answer to this question will not perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results
be expected to generalize other settings.
It is fine to include aspirational goals as motivation long it is clear that these goals are not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification:
The limitations are discussed in the concluding section of the paper (Sec.
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create a separate "Limitations" section their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on a few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect on the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution is low images are taken in low lighting. Or speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness.
While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, worse outcome might be that reviewers discover limitations that aren't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an impor- tant role in developing norms that preserve the integrity of the community, Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory assumptions and prool's
Question: For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [Yes]
Justification: Assumptions are stated in theorem statements and all proofs are provided in the appendix.
Guidelines:
The answer NA means that the paper does not include theoretical results. the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
34

---

All assumptions should be clearly stated or referenced in the statement of any theorems.
The proof' can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide a short proof sketch provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proof' provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental result reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Answer: [Yes]
Justification: We have made code and data public (or reference the data souree where
appropriate). It is available at
Optimization
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution be necessary to either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how replicate the results, access a hosted model (e.g., in the case of a large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
specific model and empirical evaluation, may
While NeurIPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily a new algorithm, the paper should make it clear how reproduce that algorithm. (b) If the contribution is primarily new model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution is a new model (c.g., a large language model), then there should either be way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome describe the particular way they provide for reproducibility. In the case of closed-souree models, it may be that access to the model is limited in some way (c.g., to registered users), but it should be possible for other researchers to have some path to reproducing verifying the results.
5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, as described in supplemental material?
Answer: [Yes]
Justification: All data and code are public. The link to the repository is given in the previous answer.
35

---

Guidelines:
reproducible, they The answer NA means that paper does not include experiments requiring code.
Please see the NeurIPS code and data submission guidelines (https://nips.cr cc/
public/guido Soatles//CooSSuuumsiioooeee
for more details.
While we encourage the release of code and data, we understand that this might not be possible, so "No" is an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurlPS code and data submission guidelines (https for more details.
7 //nips. cc/ httee/ubliccaee gwwcc//gi/eeeeeeeeeeee wwwooeeueisssoooiiooooooooo
The authors should provide instructions on data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. If only subset of experiments should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information possible in supplemental material (appended the paper) is recommended, but including URLs to data and code is permitted.
6. Experimental setting/details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: This information is reported in sections 3 and E In particular, the choice of matrix is often complex and variable and we provide full details for this.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core the paper to a level of detail that is necessary appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, as supplemental material.
7. Experiment statistical significance
Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Question:
Answer: [No]
Justification: We do not report error bars statistical significance due to computational and time constraints; as is common in combinatorial optimization, we report mean results over multiple problem instances.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by error bars, confi- dence intervals, statistical significance tests, least for the experiments that support the main claims of the paper.
The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call a library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
should be clear whether the error bar is the standard deviation the standard error of the mean.
36

---

It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of enrors is not verified.
For asymmetric distributions, the authors should be careful figures symmetric error bars that would yield results that are out of range (e.g. negative error rates).
show in tables or
error bars are reported in tables or plots, The authors should explain the text how they were calculated and reference the corresponding figures tables in the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on com- puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [Yes]
Justification: These details are given in Sec.
Guidelines:
The answer NA means that the paper does include experiments.
The paper should indicate the type of compute workers CPU or GPU, internal cluster, cloud provider, including relevant memory and storage.
The paper should provide the amount compute required for each of individual experimental runs as well estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make into the paper).
9. Code ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurlPS Code of Ethics
Answer: [Yes]
Justification: Given the focus on foundations, we do not have human subjects or sensitive data and we do not perceive risks of harm misuse.
Guidelines:
The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., there a special consid- eration due to laws or regulations in their jurisdiction).
10. Broader impacts
Question: Does the paper discuss both potential positive societal impacts and negative socictal impacts of the work performed?
Answer: [NA]
Justification: The societal impact of improved combinatorial optimization or neural combi- natorial optimization techniques are minimal.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
If the authors answer NA or No. they should explain why their work has no socictal impact why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
37

---

The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology.
If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition attacks, mechanisms for monitoring misuse, mechanisms to monitor how system learns from feedback over time, improving the efficiency and accessibility ML).
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: The paper does not contain data or models that have a high risk for misuse.
Guidelines:
The answer NA means that the paper poses no such risks.
Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Question: Are the creators or original owners of assets (e.g.. code, data, models), used in the paper. properly credited and are the license and terms of use explicitly mentioned and properly respected?
Answer: [Yes]
Justification: We provide citations to all data and code we use. All material is licensed for such use.
Guidelines:
The answer NA means that the paper does not use existing assets.
The authors should cite the original paper that produced the code package dataset.
The authors should state which version of the asset is used and, if possible, include URL.
The name of the license (c.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided.
If assets are released, the license, copyright information, and terms of use in the Datasets.papersstt paperswi thcode. cww/ datasets package should be provided. For popular datasets,
has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided.
38

---

If this information is not available online, the authors are encouraged to reach out to the asset's creators.
13. New assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification: There are no new assets released in this paper.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/mmoel as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent obtained from people whose asset is used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsoureing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well details about compensation (if any)?
Answer: [NA]
Justification: There were no human subjects in this work.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included the main paper.
According the NeurIPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional review board (IRB) approvals or equivalent for research with human subjects
Question:
Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?
Answer: [NA]
Justification:
No human subjects or study participants were used in this research.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Depending the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and expect authors to adhere to the NeurlPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such the institution conducting the review.
16. Declaration of LLM usage
39

---

Question: Does the paper describe the usage of LLMs if it is important, original, or non-standard component of the core methods this research? Note that if the LLM used only for writing, editing, or formatting purposes and does not impact the core methodology, scientific rigorousness, or originality of the research, declaration is not required.
Answer: [NA]
Justification: LLMs were not used in the core methods of this research.
Guidelines:
The answer NA means that the core method development in this research does not involve LLMs any important, original, non-standard components.
Please refer to our LLM policy (https://neure) http://pp...pc/ooeeeeeee//////// for what should or should not be described.
40


Reviews:
[{"review_id": "583hua3LwP", "reviewer": "Reviewer_6Krw", "summary": "This paper presents a continuous and almost everywhere differentiable extension from functions over permutations to functions over doubly stochastic matrices. The starting point is the Birkhoff decomposition, which allows us to decompose any doubly stochastic matrix into a convex combination of a relatively small number of permutation matrices. If the doubly stochastic matrix $A$ is decomposed as $\\sum_k \\alpha_k P_k$, where the $P_k$ are permutation matrices, we might extend a function $f$ of permutation matrices to doubly stochastic matrices by defining $F(A) = \\sum_k \\alpha_k f(P_k)$. This would indeed be an extension, since when $A$ is a permutation matrix the only valid way to express it as a convex combination of permutation matrix is to put all the weight on $A$ itself. The main challenge with this approach is that the decomposition is not unique, and so depending on how we choose the decomposition to use for each doubly stochastic matrix $A$, the function $F$ may behave erratically. The core idea behind this paper is to introduce a total ordering over permutation matrices and to use that ordering to standardize the decomposition (by preferring to put weight on early permutation matrices in the ordering). The authors argue several nice properties hold with this modification: the extension $F$ is continuous and almost everywhere differentiable. Moreover, they define a collection of total orderings where the decomposition can also be computed in polynomial time (the orderings are those induced by inner products with scoring matrices, and as long as the scoring matrices have a bit of noise, this is guaranteed to be an ordering without ties). In this case, the decomposition of an $n \\times n$ doubly stochastic matrix can be computed on $O(n^5)$ time. This means that the extension can be evaluated in $O(n^5)$ time to obtain the decomposition, and $O(n^2)$ calls of the function $f$ over permutations.\n\nIn addition to the above results, the authors also provide a number of other useful results: a version of the Frank-Wolfe optimizer that can be used to optimize $F$ without needing an expensive projection back onto the set of doubly stochastic matrices, tricks for changing the scoring matrix used to define the ordering to prevent a current optimization iterate from being at a local minima of the extended function, and extensions for functions of rooted binary trees.\n\nFinally, the authors apply their extension to several combinatorial optimization problems (the quadratic assignment problem, traveling salesperson problem, and directed feedback arc set problem). Surprisingly, the extension finds better solutions than Gurobi and heuristic methods for QAP and DFASP under some regimes.", "questions": "Please see high level questions in strengths and weaknesses.", "limitations": "Yes.", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 4, "clarity": 4, "significance": 3, "originality": 3, "strengths_and_weaknesses": "# Strengths\n- The paper is very well written and despite being about a topic that I have relatively little background in, I was able to understand the main claims and arguments of the paper.\n- Deep learning has shown that gradient based optimization can often work well even when we cannot prove that it should (e.g., when run on non-convex functions with many local minima). This paper provides a framework allowing the application of that machinery to combinatorial optimization problems over permutations, and the experiments show some evidence that this approach can in some cases outperform state-of-the-art solvers.\n- The authors go well beyond establishing the basic idea and provide nice ideas for the practical application of their extension (e.g., updating the score matrix used to define the total ordering to escape local minima during optimization).\n\n# Weaknesses\n- The computational cost of computing the extension once is polynomial, but still relatively expensive ($O(n^5)$ time to compute the decomposition and $O(n^2)$ evaluations of the original function on permutations). In the experiments the authors only compute the first $k = 5$ coefficients in the decomposition. I am curious how much this cost can be reduced. I also wonder if there is any way to \"warm-start\" the decomposition so that if $A'$ is close to $A$ (i.e., maybe one step away in the gradient-based optimization), we can leverage the decomposition of $A$ to more efficiently decompose $A'$.\n- It seems like the choice of scoring matrix $S$ might affect some properties of the extension (e.g., is the Lipschitz constant of the extension dependent on the score $S$?). In addition to the paper's results on escaping local minima, it could be interesting to explore whether the choice of $S$ has other impacts that could be optimized for.", "was_revised": true, "final_justification": "The authors adequately addressed my comments and after reading the other reviews I do not have any concerns with the paper so I will maintain my positive score.", "rebuttal": "We are thankful for the especially thorough review, and glad the reviewer has a strong understanding of our work. We now address the questions raised.\n\nThe reviewer's suggestion of warm starting the decomposition for $X$ using the previous GD iterate $X'$ is very interesting. In an extreme case, e.g. if $X'$ is a permutation matrix, the previous iterate's decomposition may contain very little information about the current decomposition (as the previous decomposition only has a single element). Therefore it seems challenging to reduce the worst-case complexity from $O(n^3)$ using this approach. On the other hand, in most cases the decompositions of $X$ and $X'$ should be very similar, so using a warm-start could offer a substantial practical speedup. This is an encouraging direction for future improvements. We have explored this observation to some degree in our implementation, where we have made the engineering choice to use a hash table that saves the results of common maximum weight matching problems. Using this technique, any maximum weight matching problems solved for computing the decomposition in the previous iterate will not be recomputed if these problems arise in computing the new decomposition. \n\n It is compelling to consider, as the reviewer suggests, how the properties of the extension $F_S$ are dependent on $S$. We note a bound on the Lipschitz constant of the extension need not depend on $S$. The relevant quantity is instead $\\Delta_f = \\max_P f(P) - \\min_Pf(P)$. Roughly, this arises because if two doubly stochastic matrices $A$ and $A'$ are within an $L_1$ distance $\\epsilon$, then the $L_1$ difference in the weights  $\\alpha_i$ of their decompositions is $O(\\epsilon)$ (note that each coefficient is a 1-Lipschitz function of the doubly stochastic matrix). Thus, in the worst case, $|F(A) - F(A')|$ is $O(\\epsilon \\Delta_f)$, where we have suppressed the dependence on $n$. We will add a rigorous proof to the manuscript showing a Lipschitz bound independent of $S$. There may be other interesting ways in which the extension is dependent on $S$. These are compelling avenues for future work."}, {"review_id": "aBFYVJhe5G", "reviewer": "Reviewer_Qomh", "summary": "The paper proposes an approach for providing an differentiable extension to loss functions on permutation matrices using a novel decomposition of the Birkhoff polytope based on a total ordering of permutations created by the inner product to a particular 'score matrix'. The extension is differentiability almost everywhere since the decomposition is unique. Nonetheless, the extension may not be non-convex and the score matrix can be changed during operation to escape local minima. Additionally the extension has desirable rounding properties such that the permutation matrix obtained from the rounding operation to obtain the best permutation matrix within the decomposition is no worse than the extension's loss. Additionally if the score function is close to permutation than the rounding will only do better.", "questions": "The definition of $\\mathcal{P}(B_k)$ could be made clearer. I understand it as the set of all permutations that are non-negative components of $B_k$, but note that $B_k$ itself is not double stochastic. Line 219 discusses maximum weight matching problem for the Birkhoff, but the residual $B_k$ isn't in the Birkhoff polytope. This should be clarified. \n\nIn Algorithm 3 and 4,  there could be a discussion of the necessitate of starting from a *dense* random double stochastic matrix $A$ (which should really be denoted $A_1$) compared to either the dense barycenter (a uniform double stochastic matrix) or a sparse combination. Especially since only 5 terms are initially used. Clearly there will be issues if the initial $A_1$ is not dense. But does it need to be random? Especially since the score matrix is is used the barycenter seems a logic initialization. \n\nI assume that automatic gradient are necessary for the Frank-Wolfe algorithm as the relationship between extension and A through $\\alpha_i$ is non-trivial to express. \n\nMinor:\nLine 319 minor typo \"simple constraint\" -> \"simple constraints\"", "limitations": "Limitations are mentioned, but more discussion could be made.", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 4, "clarity": 3, "significance": 4, "originality": 3, "strengths_and_weaknesses": "Strength:\n\nThe paper contributes to combinatorial optimization through the use of differentiable extension. \n\n\nWeakness:\n\nThere a few times things could be made clearer as noted in Questions below. \n\nThe practical implications of the truncation to 5 terms in the Birkhoff decomposition are not addressed by any of the theoretical results. This creates a gap between the theory and practice that deserves more discussion.", "was_revised": true, "final_justification": "The paper provides a novel approach for an interesting combinatoric problems. The answers further clarified some minor points. There is some limitation for larger problems expressed by other reviewers. But this is a useful contribution. Strong accept.", "rebuttal": "We appreciate the thoughtful review and now respond to points raised. \n\nWe agree that our results could be strengthened by closing the gap between theory and practice created by truncation. To this end, we will introduce a section formally introducing the truncated extension $F_S^{(k)} := \\frac{1}{\\sum_{i=1}^{k} \\alpha_i} \\sum_{i=1}^{k}\\alpha_iP_i$, proving that all rounding properties  (Properties 3, 4-1, and 4-2) hold in the presence of truncation. \nThis result follows easily since these properties only utilize characteristics of the first term in the decomposition and that extension is formed from a convex combination of discrete objective values. We will also add a proof of Theorem 2.11 in the truncating case. Again this is straightforward as the proof of this theorem only requires evaluating properties of the first term in the decomposition (see duplicate comment to reviewer S7f4). \n\nThank you for suggesting the need for additional clarity with regards to $\\mathcal P$. We define (Definition 2.2) $\\mathcal P$ only as an operator on matrices $X$ with non-negative entries. For such matrices $\\mathcal P(X)$ is the set of all permutations $P$ such that $X(i,j) >0$ if $P(i,j) = 1.$ This set is different than the set of permutations appearing in the Birkhoff decomposition; it contains all permutations that could potentially appear in a Birkhoff decomposition of $X.$  Line 219 discusses finding the permutation in $\\mathcal P(X)$ that comes first in the total order. This operation does not require computing the decomposition, so there is not a remainder to address. \n\nRegarding the necessity of random initialization, since the extension is a.e. differentiable, random initialization ensures we initialize at a point that is differentiable (almost surely). This approach also has the advantage that optimizing from multiple initializations may yield superior results. However, it is very interesting to consider starting from a fixed doubly stochastic matrix. In this setting we can use multiple score matrix initializations to get multiple optimization trajectories. We performed experiments showing the effectiveness of initializing at the barycenter. We picked 14 problem instances from the QAPLIB with sizes from $n = 12$ to $n = 35$. We use identical hyperparameter (truncation depth $k = 5$, learning rate $\\eta = .001$,  time budget $t = 2n$ seconds) and, for random and barycenter initialization, we evaluate the final average assignment costs across all 14 problem instances. These experiments show that random initialization version has avg. gap to optimal of 4.14\\%, whereas the barycenter initialization has an avg. gap of 4.42\\%. We will include these experiments in the manuscript. \n\nWe do, indeed, use automatic gradients in practice and will include a comment stating this in the manuscript.  As the expression for the $k$th term in the decomposition includes $k$ nested min functions, the gradient is quite challenging to express."}, {"review_id": "K6fcyGQX8w", "reviewer": "Reviewer_s7f4", "summary": "The paper presents Birkhoff Extension (BE), an almost-everywhere differentiable continuous extension that lifts any function defined on permutation matrices to the Birkhoff polytope. The key idea is a score-induced Birkhoff decomposition, which selects permutations in a fixed order determined by a score matrix. The proposed method exhibits several appealing properties, such as minima correspondence and a rounding guarantee. The authors further develop a Frank-Wolfe optimizer with optional dynamic score updates to help escape local minima, and they apply this framework to unsupervised neural combinatorial optimization. Experimental results on QAP, TSP, and DFASP tasks demonstrate the effectiveness of the proposed approach.", "questions": "1. How does the number of decomposition terms $k$ affect solution quality? A theoretical bound or empirical curve would be helpful.\n\n2. Can the method scale to larger problems (e.g., $n=1000$) with manageable runtime?\n\n3. Can you provide sensitivity analysis for the key parameters like truncation depth and score update frequency?\n\n4. How does this method compare to direct optimization over the Birkhoff polytope (e.g., Sinkhorn, PGD)?", "limitations": "yes", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 2, "significance": 3, "originality": 3, "strengths_and_weaknesses": "## Strengths\n1. The score-induced Birkhoff decomposition is a novel construction that selects a permutation order based on a score matrix, effectively resolving the non-uniqueness of the decomposition. This design enables differentiable extensions of discrete objectives and facilitates gradient-based optimization over permutations.\n\n2. The theoretical analysis is solid, covering differentiability, minima correspondence, and a rounding guarantee that ensures relaxed solutions can be safely mapped back to valid permutations without degrading objective quality.\n\n3. The proposed Frank-Wolfe algorithm is general and versatile, applicable to both neural combinatorial optimization and classical matching-based problems. Its effectiveness is empirically validated through comprehensive experiments on QAP, TSP, and DFASP tasks.\n\n## Weaknesses\n1. The decomposition has $O(n^5)$ computational complexity, and while truncated in practice, the paper lacks an analysis of how such truncation affects solution accuracy or theoretical guarantees.\n\n2. The experimental insights are limited. The key hyperparameters, such as truncation depth and score update frequency, are fixed without accompanying ablation studies or sensitivity analysis.\n\n3. The advantage of using the Birkhoff decomposition over direct continuous optimization on the Birkhoff polytope (e.g., using Sinkhorn or projected gradient methods) is not clearly justified, and a more detailed discussion would be beneficial.", "was_revised": true, "final_justification": "This rebuttal address most of my concerns, and I decide to raise my rating.", "rebuttal": "We thank the reviewer for their constructive comments. We have performed a variety of ablation studies in response.\n\nIncluding evaluation of the dependence on truncation and score update frequency is a great suggestion. We have performed these experiments and will update the manuscript to include curves showing these dependencies. Our results show that, for $k \\in \\{3,5,10,15,20\\}$ the greatest difference in performance is about 0.88\\% in terms of avg. gap from optimal values. Furthermore, the difference from $k = 5$ to $k = 10$ (which achieves the best performance) is about 0.24\\%. Note that we do not include the plots here due to restrictions from Neurips on the response format.\n\n| **Truncation size (k)** | **Gap from optimal % ↓** |\n|-------------------------|--------------------------|\n| $k=3$                   | 4.78\\%                   |\n| $k=5$                   | 4.14\\%                   |\n| $k=10$                  | 3.90\\%                   |\n| $k=15$                  | 3.92\\%                   |\n| $k=20$                  | 4.34\\%                   |\n\nWe also performed ablation experiments on score update frequency that we will amend the paper to include. In these experiments `update_score` is set to `True` every $m\\in \\{1,5,20,10,50\\}$ epochs. These experiments show that the performance improves with higher frequency updates, that is, smaller $m$. Therefore, the best choice is always set `update_score = True`.\n\n| **Update frequency (m)** | **Gap from optimal % ↓** |\n|--------------------------|--------------------------|\n| $m=1$                    | 4.14\\%                   |\n| $m=5$                    | 4.56\\%                   |\n| $m=10$                   | 4.57\\%                   |\n| $m=20$                   | 5.42\\%                   |\n| $m=50$                   | 5.78\\%                   |\n\nDespite the $O(n^3)$ complexity (with truncation), our method does scale to large instances. For example, for QAP, we test on QAPLIB which has an instances of size $n =256$. For this problem we get within 0.8\\% of the best known solution using a runtime of 417 seconds.\n\nWe agree that additional analysis on the truncated case would be beneficial. To this end, we will introduce a section formally introducing the truncated extension $F_S^{(k)} := \\frac{1}{\\sum_{i=1}^{k} \\alpha_i} \\sum_{i=1}^{k}\\alpha_iP_i$. proving that all rounding properties  (Properties 3, 4-1, and 4-2) hold in the presence of truncation and renormalization of coefficients.  \nThis result follows easily since these properties only utilize characteristics of the first term in the decomposition and that extension is formed from a convex combination of discrete objective values. We will also add a proof of Theorem 2.11 in the truncating case. Again this is straightforward as the proof of this theorem only requires evaluating properties of the first term in the decomposition (ee duplicate comment to reviewer Qomh).\n\nThe reviewer brings up a valuable point that there are many ways to optimize over the Birkhoff polytope: Frank--Wolfe, gradient descent + Sinkhorn, or, potentially, projected gradient descent. We performed experiments with both Frank--Wolfe and gradient descent + Sinkhorn, and found that Frank--Wolfe offered superior optimization results. In particular, using Frank-Wolfe yields an avg. gap from optimal of 4.14\\%, whereas gradient descent + Sinkhorn gives a gap of 4.55\\%. We will edit the manuscript to include these in the appendix.\n\nFor each of these Birkhoff-polytope optimization methods, an objective over the polytope is required. Although there might be many options for this objective, we propose Birkhoff extension for its several advantageous properties, including differentiability and guaranteed rounding. For a comparison with an alternative objective, see the TSP experiments, where we evaluate optimization over a quadratic program relaxation. Here we demonstrate that superior results are achieved by initially optimizing the QP relaxation and then using this solution as the score matrix initialization. We invite further comments if there is a more direct objective we should consider."}, {"review_id": "f8C9jmuL33", "reviewer": "Reviewer_TNph", "summary": "The authors tackle optimization problems over the space of permutations.\n\nThey introduce **Birkhoff Extension (BE)**, a continuous, almost-everywhere differentiable relaxation of permutation functions defined on the Birkhoff polytope (the convex hull of permutation matrices).\n\nThey propose a scoring-matrix-based method for rounding.\n\nThey show good theoretical properties about rounding.\n\nThey also design optimization techniques, e.g., Frank-Wolfe for preserving doubly-stochasticness, and updating the scoring matrix for escaping from local minima.\n\nOn benchmark instances of QAP (from QAPLIB), random-Euclidean TSP, and Directed Feedback Arc Set (on random graphs), BE matches or outperforms Gurobi, classical heuristics, and existing approximation algorithms.", "questions": "- [Q1] Property 4-2 tells that the rounding result cannot be worse than the scoring matrix, but I am wondering whether the rounding can actually get better than the scoring matrix when S is “sufficiently close to some permutation P”.\n- [Q2] Theorem 2.11 only tells us the existence of such an S, but I am wondering whether we can have any results on, e.g., how easily we can find such an S.", "limitations": "Yes. Limitations are discussed in Section 4, “Concluding Remarks”.", "rating": 4, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 3, "significance": 2, "originality": 4, "strengths_and_weaknesses": "# Strengths\n\n- **[S1] Important Problem**: Optimization problems on permutations are very important with both theoretical and empirical values.\n- **[S2] Novel Differentiable Birkhoff Decomposition:** Although the classical Birkhoff Decomposition has been well known for a long time, the proposed differentiable Birkhoff decomposition based on the scoring matrix is novel and interesting.\n- **[S3] Experiments with Broad Applicability: The proposed framework works across a variety of problems.**\n\n# Weaknesses\n\n- **[W1] High Time Complexity**: O(n^5) time complexity can be prohibitive for large-scale problems. Also, Frank-Wolfe and dynamic score bring additional overheads.\n- **[W2] Probably Not-Very-Practically-Meaningful Theory:**\n    - See [Q1] and [Q2] below\n- **[W3] Missing References:**\n    - Min et al. \"Unsupervised learning for solving the travelling salesman problem.\" NeurIPS’23\n        - This work also discusses unsupervised learning for TSP\n    - Bu et al. \"Tackling prevalent conditions in unsupervised combinatorial optimization: Cardinality, minimum, covering, and more.\" ICML’24\n        - This work also discusses theoretical properties of objective extension and rounding in unsupervised combinatorial optimization, as well as continuous relaxation for problems on trees", "was_revised": true, "final_justification": "The authors tackle optimization problems over the space of permutations.\n\nThey introduce **Birkhoff Extension (BE)**, a continuous, almost-everywhere differentiable relaxation of permutation functions defined on the Birkhoff polytope (the convex hull of permutation matrices).\n\nThey propose a scoring-matrix-based method for rounding.\n\nThey show good theoretical properties about rounding.\n\nThey also design optimization techniques, e.g., Frank-Wolfe for preserving doubly-stochasticness, and updating the scoring matrix for escaping from local minima.\n\nOn benchmark instances of QAP (from QAPLIB), random-Euclidean TSP, and Directed Feedback Arc Set (on random graphs), BE matches or outperforms Gurobi, classical heuristics, and existing approximation algorithms.\n\n# Strengths\n\n- **[S1] Important Problem**: Optimization problems on permutations are very important with both theoretical and empirical values.\n- **[S2] Novel Differentiable Birkhoff Decomposition:** Although the classical Birkhoff Decomposition has been well known for a long time, the proposed differentiable Birkhoff decomposition based on the scoring matrix is novel and interesting.\n- **[S3] Experiments with Broad Applicability: The proposed framework works across a variety of problems.**", "rebuttal": "We greatly appreciate the reviewers insightful comments and suggestions. To begin, we are thankful for the suggested references in W3, and will add these to the manuscript.  \n\n\nRegarding Q1, if the score matrix is close to a permutation $P^\\*$, then the rounding scheme can indeed produce a solution better than $P^\\*$. All that is required for this to happen is that some permutation in the decomposition is better than $P^\\*$. While Property 4-2 gives that  at least one permutation will be $P^\\*$, there is nothing preventing the decomposition from containing better permutations. \n\nTheorem 2.11 shows the existence of an $S$ that allows for escaping a local minimum. However, the proof is not constructive, as, indeed, the score $S$ used in the proof is a permutation minimizing the permutation function $f$. Thus, the hardness of finding this $S$ is the same as optimizing $f$. However, as suggested in Q2, an interesting avenue for future research would be to investigate whether there are score matrices with desirable theoretical properties that are not computationally hard to find. For now, Theorem 2.11 only shows that the approach of changing $S$ always has the potential to escape local minima but does not provide any procedure for finding such an $S$."}]

Author rebuttal:

Meta-review: This paper presents a continuous and almost everywhere differentiable extension from functions over permutations to functions over doubly stochastic matrices, which is based on the Birkhoff decomposition. The proposed method exhibits several appealing properties, such as minima correspondence and a rounding guarantee. The authors further develop a Frank-Wolfe optimizer with optional dynamic score updates to help escape local minima, and they apply this framework to unsupervised neural combinatorial optimization. Experimental results on QAP, TSP, and DFASP tasks demonstrate the effectiveness of the proposed approach.

All reviewers agree that the paper makes an interesting and novel theoretical contribution on an important problem. Moreover, the experimental results are promising and clearly demonstrate the potential of the theoretical findings and the proposed method. All potential weaknesses identified by the reviewers (in particular, the high polynomial running time and the truncation in practice, which creates an inconsistency between the theoretical results and the experiments) have been addressed convincingly by the authors. Overall, after the discussion, no significant weaknesses remain.
Decision: Accept (poster)