\section{Advances in Combinatorial Approaches}
\label{sec:advances_in_combinatorial_approaches}

In this section, we review recent advances in combinatorial approaches, covering relaxations and heuristic approaches. Further, we report on methods ensuring geometric consistency and global optimality. The methods are summarized in \Cref{tab:combinatorial-comp}.

\begin{table*}[htb]
\centering
\renewcommand{\arraystretch}{1.3} % optional: increases row height
\begin{threeparttable}
\begin{tabular}{c||c|c|c|c}
\textbf{Method} & 
\makecell{\chg{Geometrically} \\ \chg{Consistent}} & 
\makecell{\chg{Globally} \\ \chg{Optimal}} & 
\makecell{\chg{Scalability}} &
\makecell{\chg{Shape}\\\chg{Collection}}\\
\hline
\hline
% RELAXATIONS
\makecell[l]{Bending Graphs \cite{saleh2022na-bendingGraphs}} & promoted via loss & \no  & - & \no \\ % gc: idk, probably not guaranteed but only encouraged, but couldnt find in paper, scalability: keine ahnung
\hline
\makecell[l]{HiPPI \cite{bernard2019-hippi}} & \no & \no & \makecell{few landmarks, \\ many shapes} & \yes \\ % gc should be loss like since they use energy that is maximized, hard to say whether outscaled since multi matching
\hline
\makecell[l]{\cite{benkner2020-adiabatic}} & \no & -\add{\tnote{\dag}} & logical qubits: \(n  = 3\) to \(4\)\add{\tnote{\(\ast\)}} & \no \\ % CHANGE, n=3
\hline % do not say what W is (is pairwise cost)

% HEURISTICS
\makecell[l]{\cite{roetzer2022in-scalableCombinatorial}} & \yes & \no & \makecell{order of magnitude: 
\(\sim 10^2\) faces} & \no \\ % sca: 30min=500 faces -> bad/no
\hline
\makecell[l]{Simulated Annealing \cite{holzschuh2020in-simulatedAnnealing}} & \no & \no & \makecell{\(\mathcal{O}\big(v \log v \sqrt{v}\big)\) \\ \(v\): number of vertices} & \no \\ % local extrema only less likely; no timings; use geo dist as pairwise dist
\hline
\makecell[l]{Q-Match \cite{benkner2021in-qmatch}} & \no & -\add{\tnote{\dag}} & logical qubits: \(n  = 16\)\add{\tnote{\(\ast\)}} & \no \\ % n=30-40?, only a probability to get globally optimal solution?; also geo dist as pairwise dist
\hline
\makecell[l]{CCuantuMM \cite{bhatia2023in-ccuantumm}} & \no & -\add{\tnote{\dag}} & logical qubits: \(n  = 40\)\add{\tnote{\(\ast\)}} & \yes \\ % n=?, do not mention global optimality at all -> is probably not guaranteed; W also based on geo dist
\hline
\makecell[l]{\cite{Sahillioglu2018-first-genetic}} & \no & \no & \makecell{\(\mathcal{O}\big(v n^2\big)
\) \\ \(v\): number of vertices \\ \(n\): number of landmarks} & \no \\
\hline
\makecell[l]{ENIGMA \cite{Edelstein2020-ENIGMA-elastic}} & promoted via loss & \no & \makecell{order of magnitude: \(\sim 10^3\) faces} & \no \\
\hline


% GEOMETRIC CONSISTENCY
\makecell[l]{DiscoMatch \cite{roetzer2024in-discomatch}} & \yes & \yes & \makecell{order of magnitude: \(\sim 10^3\) faces \\(multi-resolution)} & \no \\ % g opt in 90% (empirically) -> no; scalab: fastglobopt ist schneller, aber nicht super viel, aber discomatch braucht viel ram; 
\hline
\makecell[l]{\cite{ehm2024combi-geometricallyConsistentPartial}} & \yes & \yes & \makecell{order of magnitude: \(\sim 10^3\) faces \\(multi-resolution)} & \no \\ % go: low res map is go -> is then refined; scalab: handles hunderts of vertices, eval with 100 faces, time budget of 2 hours
\hline
\makecell[l]{SpiderMatch \cite{roetzer2024combi-spidermatch}} & (\yes)\add{\tnote{\ddag}} & \yes & \makecell{order of magnitude: \(\sim 10^3\) faces} & \no \\ % sc: ok, 1500f = 220s % gc: should be yes
\hline
\makecell[l]{\cite{roetzer2025combi-fastGlobOpti}} & \yes & \yes & \makecell{order of magnitude: \(\sim 10^3\) faces} & \no \\ 
\hline

% GLOBAL OPTIMALITY
\makecell[l]{MINA \cite{bernard2020ex-mina}} & promoted via loss & \yes & \makecell{order of magnitude: \(\sim 10^2\) control points} & \no \\ % go: yes, but not always if runtime limited?, sigma rather says no; also use gd?, can set max threshold to gd diff
\hline
\makecell[l]{SIGMA \cite{gao2023ex-sigma}} & promoted via loss & \yes & \makecell{order of magnitude: \(\sim 10^3\) faces} & \no \\ % go: not guaranteed, but only better than other methods (more likely go) % scalability: eig gut
% in spidermatch, mina and sigma are (yes) for global opt


%\hline
%\makecell{METHOD} & GEOM CONSISTENT & GLOBALLY OPTIMAL & SCALABLE & COLLECTION \\ 


\end{tabular}
\begin{tablenotes}
\item[\add{\(\dag\)}] \add{Measures optimal solution with certain probability (quantum method).}
\item[\add{\(\ast\)}] \add{Designed for quantum hardware.} % CHANGE
\item[\add{\ddag}] \add{Depends on choice of cycle.}
\end{tablenotes}
\end{threeparttable}
\caption{Summary of recent combinatorial methods. \add{We indicate whether a method enforces geometric consistency via matching constraints \yes, promotes it via a loss term (excluding losses based on geodesic distances), or does not address it explicitly \no. For global optimality, we mark \yes \,\!\! if a method can certify a globally optimal solution. Scalability indicates the problem size a method can handle (if available from the paper). The data provides only a rough estimate, as comparisons are difficult when different hardware is used across papers. Finally, we indicate if a method was designed for shape collections \yes or not \no.}}
\label{tab:combinatorial-comp}
\end{table*}

%
\subsection{Relaxations}
Bending Graphs \cite{saleh2022na-bendingGraphs} introduces a hierarchical learning framework that integrates both local patch-level and global shape-level information. The optimal transport problem is solved via the Sinkhorn algorithm, augmented with mesh connectivity information. Specifically, confidence values are used to propagate high-confidence features to neighboring low-confidence features along the mesh graph. This encourages nearby points to be matched to nearby locations on the target region. The entire pipeline is optimized in an end-to-end manner.

In HiPPI \cite{bernard2019-hippi}, a general framework for matching shape collections with theoretical convergence guarantees is proposed. At its core, the method parametrizes the pairwise matchings using a template guaranteeing cycle-consistency and enabling scalability for shape collections. This design choice leads to a fourth-order polynomial as objective function. The matching is found by using a higher-order projected power iteration method followed by a projection onto the valid set of matchings. 

\cite{benkner2020-adiabatic} reformulate the quadratic assignment problem with permutation matrix constraints as an unconstrained problem suitable for\del{efficient} execution on \add{adiabatic} quantum hardware. The key insight is that the permutation matrix constraints can be encoded as QUBO (Quadratic Unconstrained Binary Optimization) problem with soft permutation matrix constraints. It is worth noting that they do not project a relaxed solution to the space of permutation matrices, but that the soft constraints increase the probability to measure a valid solution, \ie, a permutation matrix. \add{They use a quantum annealer to solve the QUBO. As quantum annealers are limited, the authors can confirm the validity of their approach for \(n = 3\).}

\subsection{Heuristics} \
\cite{roetzer2022in-scalableCombinatorial} propose a combinatorial solver for non-rigid shape matching that can handle partial-to-partial correspondences. Using the geometric consistency formalism introduced by \cite{windheuser2011-geometricConsistency}, they relax the resulting linear program with a Lagrange dual formulation and use a new heuristic for primal rounding which respects the geometric consistency constraints. 

A common approach for finding a good matching is to start with a permutation matrix \(P\) and update the matrix iteratively\add{, for example} by applying an \emph{r-cycle}. An r-cycle permutes disjoint indices \(i_1, \dots, i_r\). Simulated Annealing \cite{holzschuh2020in-simulatedAnnealing} is used as an heuristic optimization scheme, where in each iteration, the current solution and a neighboring solution, obtained by applying one 2-cycle, are evaluated in terms of energy defined by the QAP. The neighbor is accepted with a certain \emph{probability} based on the energy change and current temperature (iteration). A generalization of this approach, in quantum computing context, is proposed in Q-Match \cite{benkner2021in-qmatch}. Instead of applying only one 2-cycle per iteration, Benkner \etal choose \(k\)-many disjoint candidate cycles \(c_i\) (Cyclic \(\alpha\)-Expansion) and use \add{adiabatic} quantum computing to solve the (still NP-hard, but smaller and unconstrained) subproblem of deciding whether to apply \(c_i\) or not. This allows them to explore larger solution spaces simultaneously and converge in fewer iterations. This setup was later extended to incorporate cycle-consistency \cite{bhatia2023in-ccuantumm} for the matching of three shapes allowing the computation of cycle-consistent shape correspondence for shape collections by choosing one template shape. \add{Similar to \cite{benkner2020-adiabatic}, these quantum methods rely on quantum annealers and are rather limited regarding problem sizes they can handle.}

\add{Alternatively, an appropriate permutation matrix \(P\) can be determined using a genetic algorithm. This paradigm was first used by \cite{Sahillioglu2018-first-genetic} in the context of isometric shape correspondence, aiming to preserve pairwise distances between sparse landmarks. \cite{Edelstein2020-ENIGMA-elastic} extend this idea to non-isometric matching of landmarks. To evaluate a candidate landmark correspondence, they extend the sparse landmark correspondences to a dense map and score it by the resulting local metric distortion, favoring dense maps with low distortion. This score guides the evolutionary update: In each generation, the best candidates are selected, and crossover and mutation operators are applied to produce the next generation of landmark matchings.}

\subsection{Geometric Consistency} \label{sec:geometric_consistency}
Combinatorial approaches offer the possibility to include geometric consistency constraints. In \cite{roetzer2024in-discomatch}, learning-based features are used together with geometric consistency constraints to leverage the state-of-the-art matching performance of learning-based approaches, while also ensuring structurally consistent matches.
The resulting integer linear program (ILP) is solved using a massively parallel Lagrange decomposition method.
\cite{ehm2024combi-geometricallyConsistentPartial} extend this idea to geometrically consistent matching between partial shapes together with learning-based features. To make the method applicable to high-resolution shapes, the authors propose a coarse-to-fine scheme. 
\cite{roetzer2022in-scalableCombinatorial} solve an ILP over the space of orientation-preserving diffeomorphisms, a formalism introduced in \cite{windheuser2011-geometricConsistency}. 

In \cite{roetzer2024combi-spidermatch}, a path-based formalism is used, where a surface is modeled as a long self-intersecting curve. The shape matching problem is then reformulated as finding the shortest path in the product graph of the curve and the graph of the second surface. Global geometric consistency is obtained by additional constraints such as self-intersection constraints.

\cite{roetzer2025combi-fastGlobOpti} represent the surface of the source shape as a collection of cyclic graphs, where each face is represented by a surface cycle. For each cycle, they build the product graph of that cycle with the target mesh. By coupling the collection of product graphs, they obtain a hyper product graph. Computing a matching in this hyper product graph together with injectivity constraints preserves neighborhoods between surface cycles and thus the matching is geometrically consistent.

\subsection{Global Optimality}
Instead of matching all elements of one surface to the other, it is possible to match only a sparse set of elements, \eg, vertices, and deduce an approximation of the correspondence for the other vertices via a deformation model. This formulation leads to a mixed integer problem (MIP) used in MINA \cite{bernard2020ex-mina} and makes it possible to solve for the global optimum, making the method independent of the initialization. Similarly, SIGMA \cite{gao2023ex-sigma} use a deformation model to formulate a MIP. Although they only match a sparse set of shapes, the deformation can be approximated by reconstructing the shape from the projected LBO and the sparse correspondences. Their formulation is provably invariant to rigid transformations and global scaling and capable of handling partial shapes.

\color{black}