# Optimize Planning Heuristics to Rank, not to Estimate Cost-to-Goal

Leah Chrestien
Czech Technical University in Prague
<EMAIL_ADDRESS>
&Tomás̆ Pevný
Czech Technical University in Prague,
and Gen Digital, Inc.
<EMAIL_ADDRESS>
Stefan Edelkamp
Czech Technical University in Prague
<EMAIL_ADDRESS>
&Antonín Komenda
Czech Technical University in Prague
<EMAIL_ADDRESS>

###### Abstract

In imitation learning for planning, parameters of heuristic functions are
optimized against a set of solved problem instances. This work revisits the
necessary and sufficient conditions of strictly optimally efficient heuristics
for forward search algorithms, mainly A* and greedy best-first search, which
expand only states on the returned optimal path. It then proposes a family of
loss functions based on ranking tailored for a given variant of the forward
search algorithm. Furthermore, from a learning theory point of view, it
discusses why optimizing cost-to-goal $h^{*}$ is unnecessarily difficult. The
experimental comparison on a diverse set of problems unequivocally supports
the derived theory.

##  1 Introduction

Automated planning finds a sequence of actions that will reach a goal in a
model of the environment provided by the user. It is considered to be one of
the core problems in Artificial Intelligence and it is behind some of its
successful applications [39, 27, 46]. Early analysis of planning tasks [31]
indicated that optimizing the heuristic function steering the search for a
given problem domain towards the goal can dramatically improve the performance
of the search. Automated optimization of the heuristic function therefore
becomes a central request in improving the performance of planners [6, 13].

Learning in planning means optimizing heuristic functions from plans of
already solved problems and their instances. This definition includes the
selection of proper heuristics in a set of pattern databases [17, 22, 33, 12],
a selection of a planner from a portfolio [25], learning planning operators
from instances [32, 56], learning for macro-operators and entanglements [9,
29], and learning heuristic functions by general function approximators (e.g.
neural networks) [45, 19, 16, 5].

The majority of research [45, 53, 19, 16] optimizes the heuristic function to
estimate the cost-to-goal, as it is well known that it is an optimal heuristic
for many best-first heuristic search algorithms including A* or IDA* [21, 28].
These works optimize the heuristic function to solve _regression problems_.
But even a true cost-to-goal $h^{*}$ does not guarantee that the _forward
search_ will find the optimal solution while expanding the minimal number of
states. This paper defines a stricter version of _optimally efficiency_ [10,
24] as follows. _Forward search_ (or its heuristic) is called _strictly
optimally efficient_ iff it expands _only_ states on one optimal solution path
returned by the search, i.e., if this optimal solution path has $l+1$ states,
the search will expand only $l$ states.

This work focuses exclusively on a _forward search_ with merit functions. It
presents theoretical arguments as to why the order (rank) of states provided
by the heuristic function is more important than the precise estimate of the
cost-to-goal and it formalizes the necessary and sufficient condition for
strictly optimally efficient search. It also argues why learning to rank is a
simpler problem than regression estimating the cost-to-goal from a statistical
learning theory point of view. The main motivation of the paper is similar to
that in [38], but the derived theory and implementation are much simpler. The
sufficiency of optimizing the rank has been already recognized in [59] for a
beam-search and in [18] for greedy best-first search (GBFS), though, unlike
this work, it cannot guarantee strict optimal efficiency.

We emphasize that this work neither deals with nor questions the existence of
a strictly optimally efficient heuristic function with the forward search. It
is assumed that the space of possible heuristic functions is sufficiently
large, such that there exist one that is sufficiently good. The theoretical
results are supported by experimental comparisons to the prior art on eight
domains.

The contributions of this paper are as follows.

  1. 1.

We state necessary and sufficient conditions for strictly optimally efficient
heuristic functions for forward search algorithms in deterministic planning.

  2. 2.

We show that when optimizing a heuristic function, it is better to solve the
ranking problem.

  3. 3.

We argue, why learning to rank is easier than regression used in learning the
cost-to-goal.

  4. 4.

We instantiate it for A* and GBFS searches leading to two slightly different
loss functions.

  5. 5.

We experimentally demonstrate on eight problems (three grid and five PDDL)
that optimizing rank is _always_ better.

##  2 Preliminaries

A search problem instance is defined by a directed weighted graph
$\Gamma=\langle\mathcal{S},\mathcal{E},w\rangle$, a distinct node
$s_{0}\in\mathcal{S}$ and a distinct set of nodes
$\mathcal{S}^{*}\subseteq\mathcal{S}$. The nodes $\mathcal{S}$ denote all
possible states $s\in\mathcal{S}$ of the underlying transition system
representing the graph. The set of edges $\mathcal{E}$ contains all possible
transitions $e\in\mathcal{E}$ between the states in the form
$e=(s,s^{\prime})$. $s_{0}\in\mathcal{S}$ is the initial state of the problem
instance and $\mathcal{S}^{*}\subseteq\mathcal{S}$ is a set of allowed goal
states. Problem instance graph weights (alias action costs) are mappings
$w:\mathcal{E}\rightarrow\mathbb{R}^{\geq 0}$.

A path (alias a plan) $\pi=({e}_{1},{e}_{2},\ldots,{e}_{l}),$
$e_{i}=(s_{i-1},s_{i})\in\mathcal{E},$ of length $l$ solves a task $\Gamma$
with $s_{0}$ and $\mathcal{S}^{*}$ iff
$\pi=((s_{0},s_{1}),(s_{1},s_{2}),\ldots,(s_{l-1},s_{l}))$ and
$s_{l}\in\mathcal{S}^{*}$. An optimal path $\pi^{*}$ is defined as a path
minimizing the cost of a problem instance $\Gamma,s_{0},\mathcal{S}^{*},$ its
value as $f^{*}=w(\pi^{*})=\sum_{i=1}^{l}w({e}_{i})$. When cost function
$w(e)=1$ optimal path corresponds to the shortest one.
$\pi_{:i}=((s_{0},s_{1}),(s_{1},s_{2}),\ldots,(s_{i-1},s_{i}))$ denotes
subplan of a plan $\pi$ containing its first $i$ actions.
$\mathcal{S}^{\pi}=\\{s_{i}\\}_{i=0}^{l}$ denotes the set of states from the
path $\pi=((s_{0},s_{1}),(s_{1},s_{2}),\ldots,(s_{l-1},s_{l})),$ and
correspondingly $\mathcal{S}^{\pi_{:i}}$ denotes the set of states in a
subplan $\pi_{:i}.$

###  2.1 Forward search algorithm

To align notation, we briefly review the forward search algorithm, of which A*
and GBFS are special cases. As we apply full duplicate detection in all of our
algorithms, there is a bijection of search tree nodes to planning states.
Forward Search for consistent heuristics, where $h(s)-h(s^{\prime})\leq
w(s,s^{\prime})$ for all edges $(s,s^{\prime})$ in the $w$-weighted state
space graph, mimics the working of Dijkstra’s shortest-path algorithm [11]. It
maintains two sets: the first called _Open list_ , $\mathcal{O},$ contains
generated but not expanded nodes; the second called _Closed list_ ,
$\mathcal{C},$ contains already expanded nodes. Parameters $\alpha$ and
$\beta$ of a merit function $f(s)=\alpha g(s)+\beta h(s)$ used to sort states
in Open list allow to conveniently describe different variants of the
algorithm as described in detail shortly below. The forward search works as
follows.

  1. 1.

Initialise Open list as $\mathcal{O}_{0}=\\{s_{0}\\}$.

  2. 2.

Set $g(s_{0})=0$

  3. 3.

Initiate the Closed list to empty, i.e. $\mathcal{C}_{0}=\emptyset.$

  4. 4.

For $i\in 1,\ldots$ until $\mathcal{O}_{i}=\emptyset$

    1. (a)

Select the state $s_{i}=\arg\min_{s\in\mathcal{O}_{i-1}}\alpha g(s)+\beta
h(s)$

    2. (b)

Remove $s_{i}$ from $\mathcal{O}_{i-1},$
$\mathcal{O}_{i}=\mathcal{O}_{i-1}\setminus\\{s_{i}\\}$

    3. (c)

If $s_{i}\in\mathcal{S}^{*}$, i.e. it is a goal state, go to 5.

    4. (d)

Insert the state $s_{i}$ to $\mathcal{C}_{i-1}$,
$\mathcal{C}_{i}=\mathcal{C}_{i-1}\cup\\{s_{i}\\}$

    5. (e)

Expand the state $s_{i}$ into states $s^{\prime}$ for which hold
$(s_{i},s^{\prime})\in\mathcal{E}$ and for each

      1. i.

set $g(s^{\prime})=g(s_{i})+w(s_{i},s^{\prime})$

      2. ii.

if $s^{\prime}$ is in the Closed list as $s_{c}$ and $g(s^{\prime})<g(s_{c})$
then $s_{c}$ is reopened (i.e., moved from the Closed to the Open list), else
continue with (e)

      3. iii.

if $s^{\prime}$ is in the Open list as $s_{o}$ and $g(s^{\prime})<g(s_{o})$
then $s_{o}$ is updated (i.e., removed from the Open list and re-added in the
next step with updated $g(\cdot)$), else continue with (e)

      4. iv.

add $s^{\prime}$ into the Open list

  5. 5.

Walk back to retrieve the solution path.

In the above algorithm, $g(s)$ denotes a function assigning an accumulated
cost $w$ for moving from the initial state ($s_{0}$) to a given state $s$.
During its execution, it always expands nodes with the lowest $f(s)=\alpha
g(s)+\beta h(s).$ Different settings of $\alpha$ and $\beta$ give rise to
different algorithms: for A*, $\alpha=\beta=1$; for GBFS $\alpha=0$ and
$\beta=1$.

_Consistent heuristics_ , which are of special interest, are called monotone
because the estimated cost of a partial solution $f(s)=g(s)+h(s)$ is
monotonically non-decreasing along the best path to the goal. More than this,
$f$ is monotone on all edges $(s,s^{\prime})$, if and only if $h$ is
consistent as we have $f(s^{\prime})=g(s^{\prime})+h(s^{\prime})\geq
g(s)+w(s,s^{\prime})+h(s)-w(s,s^{\prime})=f(s).$ For the case of consistent
heuristics, no reopening (moving back nodes from Closed to Open) in A* is
needed, as we essentially traverse a state-space graph with edge weights
$w(s,s^{\prime})+h(s^{\prime})-h(s)\geq 0$. For the trivial heuristic $h_{0}$,
we have $h_{0}(s)=0$ and for perfect heuristic $h^{*}$, we have
$f(s)=f^{*}=g(s)+h^{*}(s)$ for all nodes $s$. Both heuristics $h_{0}$ and
$h^{*}$ are consistent. For GBFS and other variants, reopening is usually
neglected. Even if the heuristic is not consistent, best-first algorithms
without the reopening, remain complete i.e., they find a plan if there is one.
Plans might not be provably optimal but are often good for planning practice
[14].

##  3 Conditions on strictly optimally efficient heuristic

Let $\mathcal{O}_{i}$ be an Open list in the $i^{\mathrm{th}}$ iteration of
the forward search expanding only states on the optimal path $\pi.$ Below
definition defines a _perfect ranking heuristic_ as a heuristic function where
state on some optimal path in the Open list,
$\mathcal{S}^{\pi_{:i}}\cap\mathcal{O}_{i},$ has _always_ strictly lower merit
value than other states in the Open list off the optimal path,
$\mathcal{O}_{i}\setminus\mathcal{S}^{\pi_{:i}}.$ The intersection
$\mathcal{S}^{\pi_{:i}}\cap\mathcal{O}_{i}$ contains always exactly one state.

######  Definition 1 (Perfect ranking heuristic).

A heuristic function $h(s)$ is a perfect ranking in forward search with a
merit function $f(s)=\alpha g(s)+\beta h(s)$ for a problem instance
$\Gamma=(\langle\mathcal{S},\mathcal{E},w\rangle,s_{0},\mathcal{S}^{*})$ if
and only if there exists an optimal plan
$\pi=((s_{0},s_{1}),(s_{1},s_{2}),\ldots,(s_{l-1},s_{l}))$ such that

  * •

$g(s)$ is the cost from $s_{0}$ to $s$ in a search-tree created by expanding
only states on the optimal path $\pi;$

  * •

$\forall i\in\\{1,\ldots,l\\}$ and $\forall
s_{j}\in\mathcal{O}_{i}\setminus\mathcal{S}^{\pi_{:i}}$ we have
$f(s_{j})>f(s_{i}).$

As an example consider a search-tree in Figure 1(a) with an optimal path
$((s_{0},s_{1}),(s_{1},s_{2}),(s_{2},s_{3}))$ and states
$\\{s_{4},s_{5},s_{6},s_{7}\\}$ off the optimal path. Then the Open list after
expanding $s_{0}$ is $\mathcal{O}_{1}=\\{s_{1},s_{4},s_{5}\\}$ and that after
expanding $s_{1}$ is $\mathcal{O}_{2}=\\{s_{2},s_{4},s_{5},s_{6},s_{7}\\}.$
The set of states on the optimal sub-path are
$\mathcal{S}^{\pi_{:1}}=\\{s_{0},s_{1}\\},$
$\mathcal{S}^{\pi_{:2}}=\\{s_{0},s_{1},s_{2}\\}.$ States on the optimal path
in the Open lists are $\mathcal{S}^{\pi_{:1}}\cap\mathcal{O}_{1}=\\{s_{1}\\}$
and $\mathcal{S}^{\pi_{:2}}\cap\mathcal{O}_{2}=\\{s_{2}\\}.$ States off the
optimal path in the Open lists are
$\mathcal{O}_{1}\setminus\mathcal{S}^{\pi_{:1}}=\\{s_{4},s_{5}\\}$ and
$\mathcal{O}_{2}\setminus\mathcal{S}^{\pi_{:2}}=\\{s_{4},s_{5},s_{6},s_{7}\\}.$
Function $f$ is perfectly ranking problem instance in Figure 1(a) iff
following inequalities holds:

$\displaystyle f(s_{1})$ $\displaystyle<f(s_{4}),\quad\quad$ $\displaystyle
f(s_{2})$ $\displaystyle<f(s_{4}),\quad\quad$ $\displaystyle f(s_{2})$
$\displaystyle<f(s_{6}),$ $\displaystyle f(s_{1})$
$\displaystyle<f(s_{5}),\quad\quad$ $\displaystyle f(s_{2})$
$\displaystyle<f(s_{5}),\quad\quad$ $\displaystyle f(s_{2})$
$\displaystyle<f(s_{7}).$

Notice that the heuristic values of states on the optimal path are never
compared. This is because if the forward search expands only states on the
optimal path, then its Open list always contains only one state from the
optimal path. The definition anticipates the presence of multiple optimal
solutions and even multiple goals, which requires the heuristic to break ties
and have a perfect rank with respect to one of them.

######  Theorem 1.

The forward search with a merit function $f(s)=\alpha g(s)+\beta h(s)$ and a
heuristic $h$ is strictly optimally efficient on a problem instance
$(\langle\mathcal{S},\mathcal{E},w\rangle,s_{0},\mathcal{S}^{*})$ if and only
if $h$ is a perfect ranking on it.

###### Proof.

_Sufficiency:_ If the conditions of a perfect ranking heuristic hold, then
there exists an optimal path such that a state on it that is in the Open list
always has the strictly lowest heuristic value. Therefore, forward search will
never move off the optimal path and is, therefore, strictly optimally
efficient.

_Necessity:_ If $h$ is a strictly optimally efficient heuristic, then the
forward search always selects the state on the optimal path, which means that
the state has the lowest heuristic value of all the states in the Open list,
which is precisely the condition of a perfect ranking heuristic. ∎

Theorem 1, despite being trivial, allows certifying that the forward search
with a given heuristic in a given problem instance is strictly optimally
efficient. This property has been discussed in [59] for Beam search, but the
conditions stated there are different because Beam search prunes states in the
Open list. Recall that the complexity class of computing a perfect ranking
heuristic function is the same as solving the problem because one implies the
other. We now remind the reader that the popular cost-to-goal does not always
lead to a strictly optimally efficient forward search.

$s_{4}$$s_{5}$$s_{6}$$s_{7}$$s_{0}$$s_{1}$$s_{2}$$s_{3}$ (a) An example of a
search tree created by expanding only states on the optimal path
$(s_{0},s_{1},s_{2},s_{3}).$ Grey nodes denote states on the optimal path and
pink nodes denote states off the optimal path.

A10B3C8D40E28434 (b) Problem instance where perfect heuristic is not strictly
optimally efficient with GBFS. Numbers on the edges denote the cost of action
and red numbers next to nodes denote the minimal cost-to-goal.

######  Example 1.

While cost-to-goal $h^{*}$ is the best possible heuristic for algorithms like
A* (up to tie-breaking) in terms of nodes being expanded, for GBFS, $h^{*}$
does not necessarily yield optimal solutions.

See Figure 2 for a counterexample. The complete proof is in the Appendix.
Appendix 8.4 also gives an example of a problem instance, where an optimally
efficient heuristic for GBFS returning the optimal solution does not exist.

##  4 Loss functions

###  4.1 Ranking loss function

Let’s now design a loss function minimizing the number of violated conditions
in Definition 1 for a problem instance $(\Gamma,s_{0},\mathcal{S}^{*})$ and
its optimal plan $\pi.$ Assuming a heuristic function $h(s,\theta)$ with
parameters $\theta\in\mathbf{\Theta},$ the number of violated conditions can
be counted as

$\mathrm{L}_{01}(h,\Gamma,\pi)=\sum_{s_{i}\in\mathcal{S}^{\pi}}\sum_{s_{j}\in\mathcal{O}_{i}\setminus\mathcal{S}^{\pi_{:i}}}\llbracket
r(s_{i},s_{j},\theta)>0\rrbracket,$ (1)

where

$r(s_{i},s_{j},\theta)=\alpha(g(s_{i})-g(s_{j}))+\beta(h(s_{i},\theta)-h(s_{j},\theta)),$
(2)

$\mathcal{O}_{i},$ $\mathcal{S}^{\pi_{:i}}$ as used in Definition 1,
$\llbracket\cdot\rrbracket$ denotes Iverson bracket being one if the argument
is true and zero otherwise, and $\alpha$ and $\beta$ are parameters
controlling the type of search, as introduced above.

In imitation learning, we assume that we have a training set
$\mathcal{T}^{\mathrm{trn}}$ consisting of tuples of problem instances and
optimal plans
$\mathcal{T}^{\mathrm{trn}}=\\{\left(\langle\mathcal{S}_{i},\mathcal{E}_{i},w_{i}\rangle,s_{0,i},\mathcal{S}_{i}^{*},\pi_{i}\right)\\}_{i=1}^{n}.$
To optimize parameters $\theta$ of a heuristic function $h(s,\theta),$ we
propose to minimize the number of wrongly ranked states over all the problem
instances in the training set,

$\arg\min_{\theta\in\mathbf{\Theta}}\sum_{(\Gamma,\pi)\in\mathcal{T}^{\mathrm{trn}}}\mathrm{L}_{01}(h,\Gamma,\pi).$
(3)

In practice, one often wants to solve problem (3) by efficient gradient
optimization methods. To do so, the Iverson bracket
$\llbracket\cdot\rrbracket$ (also called 0-1 loss) is usually replaced by a
convex surrogate such as the hinge-loss or the logistic loss used in the
Experimental section below, in which case Equation (3) becomes

$\arg\min_{\theta\in\mathbf{\Theta}}\sum_{(\Gamma,\pi)\in\mathcal{T}^{\mathrm{trn}}}\sum_{s_{i}\in\mathcal{S}^{\pi}}\sum_{s_{j}\in\bar{\mathcal{S}}^{i}}\log(1+\exp(r(s_{i},s_{j},\theta))).$
(4)

While optimization of the surrogate loss might look as if one optimizes a
different problem, it has been shown that for certain classes of functions,
optimization of surrogate losses solves the original problem (3) as the number
of training samples approach infinity [50, 35] .

###  4.2 Regression loss function

For the sake of comparison, we review a loss function used in the majority of
works optimizing the cost-to-goal in imitation learning for planning (some
notable exceptions are in the related work section). The optimization problem
with $\mathrm{L}_{2}$ loss function is defined as

$\arg\min_{\theta\in\mathbf{\Theta}}\sum_{(\Gamma,\pi)\in\mathcal{T}^{\mathrm{trn}}}\sum_{s_{i}\in\mathcal{S}^{\pi}}(h(s_{i},\theta)-h^{*}(s_{i})^{2},$
(5)

where $h^{*}(s_{i})$ is the true cost-to-goal. The optimal solution of the
optimization is the heuristic function $h^{*}.$ What is important for the
discussion below is that the $\mathrm{L}_{2}$ loss function uses _only states
on the optimal path_ and the heuristic function is optimized to solve the
_regression_ problem. Other variants of regression loss, such as asymmetric
$\mathrm{L}_{1}$ in [51], do not solve issues of regression losses discussed
in this paper.

###  4.3 Advantages and disadvantages of loss functions

#### Cost-to-goal provides a false sense of optimality

It has been already mentioned above that A* with $h^{*}$ is strictly optimally
efficient up to resolving states with the same value of $g(s)+h^{*}(s).$ This
means that for problems with a large (even exponential) number of optimal
solutions of equal cost (see [23] for examples), A* will be very inefficient
unless some mechanism for breaking ties is used. Moreover, since learning is
inherently imprecise, a heuristic close but not equal to $h^{*}$ will likely
be less efficient than the one close to a perfect-ranking heuristic.

#### Size of the solution set

The set of all functions satisfying Equation (5) is likely smaller than that
satisfying Equation (3), because the solution set of (3) is invariant to
transformations by _stricly monotone_ functions, unlike the solution set of
(5). From that, we conjecture that finding the solution of (5) is more
difficult. The precise theorem is difficult to prove since the solution of (5)
is not a solution of (3) due to tie resolution in Definition 1.

#### Optimizing cost-to-goal does not utilize states off the solution path

Assume a problem instance
$(\langle\mathcal{S},\mathcal{E},w\rangle,s_{0},\mathcal{S}^{*})$ and an
optimal solution path
$\pi=((s_{0},s_{1}),(s_{1},s_{2}),\ldots,(s_{l-1},s_{l})).$ Loss function
optimized in (5) uses only states $\\{s_{i}\\}_{i=0}^{l}$ on the optimal path.
This means that even when the loss is zero, the search can be arbitrarily
inefficient if states $\mathcal{S}\backslash\mathcal{S}^{\pi}$ off the optimal
path will have heuristic values smaller than
$\min_{s\in\mathcal{S}^{\pi}}h(s).$ The proposed ranking loss (3) uses states
on and off the optimal path and therefore, does not suffer from this problem.
This issue can be fixed if the training set is extended to contain heuristic
values for all states off the optimal path (states
$\\{s_{4},s_{5},s_{6},s_{7}\\}$ in Figure 1(a)), as has been suggested in [5,
49]. With respect to the above definitions, this is equal to adding more
problem instances to the training set sharing the same graph $\Gamma$ and set
of goal states $\mathcal{S}^{*},$ but differing in initial states. Solving all
of them to obtain true cost-to-goal _greatly increases the cost of creating
the training set_.

#### Heuristic value for dead-end states

should be sufficiently high to ensure they are never selected. An infinity in
the $\mathrm{L}_{2}$ loss would always lead to an infinite gradient in
optimization. In practice, the $\infty$ is replaced by a sufficiently large
value, but too large values can cause large gradients of $\mathrm{L}_{2},$
which might negatively impact the stability of convergence of gradient-based
optimization methods. Proposed ranking losses do not suffer this problem.

#### Goal-awareness

The drawback of the heuristic function $h(s,\theta)$ obtained by optimizing
ranking losses is that they are not goal-aware unlike the heuristic optimizing
the estimate of cost-to-goal.

#### Speed of convergence against the size of training set

From classical bounds on true error [54] and [7], we can derive (see Appendix
for details) that the excess error of ranking loss converges to zero at a rate
$\frac{\sqrt{\ln n_{p}}}{\sqrt{n_{p}}},$ which is slightly faster than that of
regression $\frac{\sqrt{\ln n_{s}}}{\sqrt{n_{s}}-\sqrt{\ln n_{s}}},$ where
$n_{s}$ and $n_{p}$ denotes the number of states and state pairs respectively
in the training set. But, the number of state pairs in the training set grows
quadratically with the number of states; therefore, the excess error of
ranking loss for the number of states $n_{s}$ can be expressed as
$\frac{\sqrt{2\ln n_{s}}}{n_{s}},$ which would be by at least
$\frac{1}{\sqrt{n_{s}}}$ factor faster than that of regression. Note though
that these rates are illustrative since the independency of samples (and
sample pairs) is in practice, violated since even samples from the same
problem instance are not independent.

#### Conflicts in training set

Unlike regression loss, the proposed family of ranking losses is potentially
sensitive to conflicts in the training set, which occur when the training set
contains examples of two (or more) different solution paths of _the same_
problem instance. This might prevent achieving a zero loss. The effect on the
heuristic depends on the composition of the training set. Either solution
paths more frequent in the training set will be preferred, or there will be
ties, but their resolution do not affect optimally efficiency. Appendix 8.5
provides an illustrative example.

##  5 Related Work

The closest work to ours is [36, 38], which adjusts heuristic value by a
policy estimated from the neural network. The theory therein has the same
motivation, to minimize the number of expanded states, though it is more
complicated than the theory presented here. For efficiency, [38] needs to
modify the GBFS such that search nodes store policy and probability of
reaching the node from the start. Moreover, in some variants (used in the
experimental comparison below) the neural network has to estimate both
heuristic value and policy.

Inspired by statistical surveys of heuristic functions in [57, 58], [18]
proposes to optimize the rank of states preserving order defined by the true
cost-to-goal. The formulation is therefore valid only for GBFS and not for A*
search. Ties are not resolved, which means the strict optimal efficiency
cannot be guaranteed. Similarly, heuristic values for domains where all
actions have zero cost cannot be optimized. Lastly and importantly, the need
to know the true cost-to-goal for all states in the training set greatly
increases the cost of its construction as discussed above. From Bellman’s
equation [49] arrives at the requirement that the smallest heuristic value of
child nodes has to be smaller than that of the parent node. The resulting loss
is regularized with $\mathrm{L}_{1}$ forcing the value of the heuristic to be
close to cost-to-goal.

In [5], A* is viewed as a Markov Decision Process with the Q-function equal to
the number of steps of A* reaching the solution. This detaches the heuristic
values from the cost-to-goal cost, but it does not solve the problem with
dead-ends and ties, since the optimization is done by minimizing
$\mathrm{L}_{2}$ loss.

Symmetry of $\mathrm{L}_{2}$ (and of $\mathrm{L}_{1}$) loss are discussed in
[51], as it does not promote the admissibility of the heuristic. It suggests
asymmetric $\mathrm{L}_{1}$ loss with different weights on the left and right
parts, but this does not completely solve the problem of estimating cost-to-
goal (for example with ties, dead-ends).

Parameters of potential heuristics [41] are optimized by linear programming
for each problem instance to meet admissibility constraints while helping the
search. On contrary, this work focuses on the efficiency of the search for
many problem instances but cannot guarantee admissibility in general.

Combining neural networks with discrete search algorithms, such that they
become an inseparable part of the architecture is proposed in [55] and [60].
The gradient with respect to parameters has to be propagated through the
search algorithm which usually requires both approximation and the search
algorithm to be executed during every inference.

A large body of literature [46, 20, 15, 2] improves the Monte Carlo Tree
Search (MCTS), which has the advantage to suppress the noise in estimates of
the solution costs. These methods proposed there are exclusive to MCTS and do
not apply to variants of the forward search algorithm. Similarly, a lot of
works [45, 53, 61, 19, 8, 16, 5] investigate architectures of neural networks
implementing the heuristic function $h(s,\theta)$, ideally for arbitrary
planning problems. We use [45, 8] to implement $h(s,\theta),$ but remind our
readers that our interest is the loss function minimized when optimizing
$\theta.$

Finally, [1] proposes to construct a training set similarly to pattern
databases [12] by random backward movements from the goal state and to adapt
A* for multi-core processors by expanding $k>1$ states in each iteration.
These improvements are independent of those proposed here and can be combined.

##  6 Experiments

The experimental protocol was designed to compare the loss functions while
removing all sources of variability. It therefore uses imitation learning in
planning [18, 19] and leaves more interesting bootstrap protocol for future
work [3]. The goal was not to deliver state-of-the-art results, therefore the
classical solvers were omitted, but their results are in the Appendix for
completeness.

The proposed ranking losses are instantiated with the logistic loss function
as in Equation (4) for A* search by setting $\alpha=\beta=1$ in (2), called
$\mathrm{L}^{*}$, and for GBFS, by setting $\alpha=0$, $\beta=1$ in (2),
called $\mathrm{L}_{\mathrm{gbfs}}.$ $\mathrm{L}_{2}$ is used as defined in
Equation (5). Since the heuristics used in this paper are not admissible, A*
search will not provide optimal solutions. It was included in the comparison
to demonstrate that heuristics optimized with respect to
$\mathrm{L}_{\mathrm{gbfs}}$ might be inferior in A* search. To observe the
advantage of ranking over regression, we defined a $\mathrm{L}_{\mathrm{rt}}$
loss function

$\mathrm{L}_{\mathrm{rt}}(h,\Gamma,\pi)=\sum_{(s_{i-1},s_{i})\in\pi}\log(1+\exp(h(s_{i},\theta)-h(s_{i-1},\theta)))$
(6)

comparing only states on the optimal trajectory. We also compare to the loss
function proposed in [49] defined as
$\mathrm{L}_{\mathrm{be}}(h,\Gamma,\pi)=\sum_{s\in\mathcal{S}^{\pi}}\max\\{1+\min_{s^{\prime}\in\mathcal{N}(s)}h(s^{\prime},\theta)-h(s,\theta)\
,0\\}+\max\\{0,h^{*}(s)-h(s,\theta)\\}+\max.\\{0,h(s,\theta)-2h^{*}(s)\\},$
where $\mathcal{N}(s)$ denotes child-nodes of state $s.$ Furthermore, we
compare to the policy-guided heuristic search [38] (denoted as
$\mathrm{L}_{\mathrm{le}}$) with GBFS modified for efficiency as described in
the paper while simultaneously adapting the neural network to provide the
policy and heuristic value.

The heuristic function $h(s,\theta)$ for a given problem domain was optimized
on the training set containing problem instances and their solutions, obtained
by SymBA* solver [52] for Sokoban and Maze with teleports and [42] 111runs the
”LAMA 2011 configuration” of the planner for Sliding Tile. Each minibatch in
SGD contained states from exactly a single problem instance needed to
calculate the loss. This means that for $\mathrm{L}_{2},$
$\mathrm{L}_{\mathrm{rt}}$, $\mathrm{L}_{\mathrm{le}}$, only states on the
solution trajectory were included; for $\mathrm{L}^{*}$,
$\mathrm{L}_{\mathrm{gbfs}}$, and $\mathrm{L}_{\mathrm{be}}$ there were
additional states of distance one (measured by the number of actions) from
states on the solution trajectory. Optimized heuristic functions were always
evaluated within A* and GBFS forward searches, such that we can analyze how
the heuristic functions specialize for a given search. To demonstrate the
generality, experiments contained eight problem domains of two types (grid and
general PDDL) using two different architectures of neural networks described
below.

#### Neural network for grid domains

States in _Sokoban_ , _Maze_ with teleports, and _Sliding-puzzle_ can be
easily represented in a grid (tensor). The neural network implementing the
heuristic function $h(s,\theta)$ was copied from [8]. Seven convolution
layers, $P_{1}\ldots P_{7}$ of the same shape $3\times 3$ with 64 filters are
followed by four CoAT blocks (CoAT block contains a convolution layer followed
by a multi-head attention operation with 2 heads and a positional encoding
layer). Convolution layers in CoAT blocks have size $3\times 3$ and have 180
filters. The network to this stage preserves the dimension of the input. The
output from the last block is flattened by mean pooling along $x$ and $y$
coordinates before being passed onto a fully connected layer (FC) predicting
the heuristic. All blocks have skip connections.

#### Neural network for PDDL domains

For _Blocksworld_ , _Ferry_ , _Spanner_ , _N-Puzzle_ , and _Elevators_ , where
the state cannot be easily represented as a tensor, we have used Strips-HGN
[44] and we refer the reader therein for details. Strips-HGN was implemented
in Julia [4] using PDDL.jl [62] extending [30]. Importantly, we have used the
hyper-graph reduplication proposed in [48] to improve computational
efficiency. Since we did not know the precise architecture (unlike in the grid
domain), we performed a grid-search over the number of hyper-graph
convolutions in $\\{1,2,3\\},$ dimensions of filters in $\\{4,8,16\\},$ and
the use of residual connections between features of vertices. The
representation of vertices was reduced by mean and maximum pooling and passed
to a two-layer feed-forward layer with hidden dimension equal to the number of
filters in hyper-graph convolutions. All non-linearities were of relu type.
The best configuration was selected according to the number of solved problems
in the validation set. Forward search algorithms were given 5s to solve each
problem instance chosen in a manner to emphasize differences between loss
functions. All experiments were repeated 3 times. All experiments on the PDDL
domains including development took 17724 CPU hours. Optimizing heuristic for a
single problem instance took approximately 0.5 CPU/h except Gripper domain,
where the training took about 8h due to large branch factors. The rest of this
subsection gives further details about domains.

#### Sokoban

The training set for the Sokoban domain contained 20000 mazes of grid size
$10\times 10$ with a single agent and 3 boxes. The testing set contained 200
mazes of the same size $10\times 10,$ but with 3, 4, 5, 6, 7 boxes. Assuming
the complexity of to be correlated with the number of boxes, we can study
generalization to more difficult problems. The problem instances were
generated by [40].

#### Maze with teleports

The training set contained 20000 mazes of size $15\times 15$ with the agent in
the upper-left corner and the goal in the lower-right corner. The testing set
contained 200 mazes of size $50\times 50,$ $55\times 55,$ and $60\times 60.$
Therefore, as in the case of Sokoban, the testing set was constructed such
that it evaluates the generalization of the heuristic to bigger mazes. The
mazes were generated using an open-source maze generator [47], where walls
were randomly broken and 4 pairs of teleports were randomly added to the maze
structure.

#### Sliding puzzle

The training set contained 20000 puzzles of size $5\times 5$. The testing set
contained 200 mazes of size $5\times 5$. These puzzles were taken from [37].
Furthermore, we tested our approach on puzzles of higher dimensions such
$6\times 6$ and $7\times 7$, all of which were generated with [37]. Thus, in
total, the testing set contained 200 sliding tile puzzles for each of the
three dimensions.

#### Blocksworld, N-Puzzle, Ferry, Elevators, and Spanner

The problem instances were copied from [44] (Blocksworld, N-Puzzle, and Ferry)
and from [34] (elevators-00-strips). Problem instances of Spanner were
generated by [43]. Blocksworld contained 732 problem instances of size 3–15,
N-Puzzle contained 80, Ferry contained 48, Elevators contained 80, and Spanner
contained 440. Table 5 in Appendix shows the fraction of -olved mazes by
breadth-first search for calibration of the difficulty. The problem instances
were randomly divided into training, validation, and testing sets, each
containing 50% / 25% / 25% of the whole set of problems. The source code of
all experiments is available at https://github.com/aicenter/Optimize-Planning-
Heuristics-to-Rank with MIT license.

|  | A* |  | GBFS
---|---|---|---|---
problem | complx. | $\mathrm{L}^{*}$ | $\mathrm{L}_{\mathrm{gbfs}}$ | $\mathrm{L}_{\mathrm{rt}}$ | $\mathrm{L}_{2}$ | $\mathrm{L}_{\mathrm{be}}$ |  | $\mathrm{L}^{*}$ | $\mathrm{L}_{\mathrm{gbfs}}$ | $\mathrm{L}_{\mathrm{rt}}$ | $\mathrm{L}_{2}$ | $\mathrm{L}_{\mathrm{be}}$ | $\mathrm{L}_{\mathrm{le}}$
Blocks |  | 100 | 100 | 100 | 99 | 100 |  | 100 | 100 | 100 | 100 | 100 | 99
Ferry |  | 98 | 98 | 100 | 92 | 100 |  | 98 | 100 | 100 | 100 | 98 | 98
N-Puzzle |  | 89 | 87 | 88 | 83 | 89 |  | 92 | 89 | 89 | 89 | 92 | 88
Spanner |  | 100 | 89 | 100 | 84 | 92 |  | 100 | 100 | 100 | 100 | 100 | 100
Elevators |  | 91 | 85 | 75 | 36 | 66 |  | 92 | 85 | 79 | 76 | 67 | 58
Sokoban | 3 boxes | 99 | 98 | 96 | 97 | 92 |  | 98 | 100 | 94 | 95 | 92 | 98
| 4 boxes | 89 | 89 | 85 | 81 | 82 |  | 87 | 91 | 84 | 83 | 84 | 84
| 5 boxes | 80 | 75 | 72 | 72 | 73 |  | 78 | 77 | 74 | 72 | 72 | 73
| 6 boxes | 76 | 69 | 59 | 51 | 53 |  | 73 | 71 | 56 | 51 | 54 | 64
| 7 boxes | 55 | 49 | 47 | 42 | 45 |  | 51 | 49 | 48 | 43 | 45 | 49
Maze w. t. | $50\times 50$ | 92 | 91 | 88 | 87 | 87 |  | 89 | 90 | 89 | 84 | 85 | 89
| $55\times 55$ | 78 | 75 | 73 | 72 | 74 |  | 74 | 75 | 74 | 72 | 75 | 74
| $60\times 60$ | 49 | 37 | 35 | 32 | 31 |  | 42 | 48 | 36 | 34 | 32 | 42
Sliding puzzle | $5\times 5$ | 88 | 83 | 84 | 80 | 82 |  | 86 | 87 | 84 | 84 | 84 | 85
| $6\times 6$ | 51 | 48 | 49 | 45 | 46 |  | 47 | 49 | 45 | 43 | 46 | 48
| $7\times 7$ | 39 | 35 | 36 | 32 | 34 |  | 35 | 36 | 35 | 32 | 34 | 35

Table 1: Fraction of solved problem instances in percent from the testing set.
The top row denotes the type of the search and the second top row denotes the
loss function against which the heuristic function was optimized. The
complexity is explicitly shown for grid domains. Standard deviations are most
of the times smaller than one percent and are provided in Table in the
Appendix.

###  6.1 Experimental results

Table 1 shows the fraction of solved mazes in percent for all combinations of
search algorithms (A* and GBFS) and heuristic functions optimized against
compared loss functions. The results match the above theory, as the proposed
ranking losses are _always better_ than the $\mathrm{L}_{2}$ regression loss.

For A* search, the heuristic function optimized against the $\mathrm{L}^{*}$
is clearly the best, as it is only once inferior to $\mathrm{L}_{\mathrm{rt}}$
and $\mathrm{L}_{\mathrm{be}}$ on the Ferry domain. As expected, heuristic
functions optimized with respect to $\mathrm{L}_{\mathrm{gbfs}}$ behave
erratically in A* search, because they were optimized for a very different
search.

For GBFS, heuristic functions optimized against the $\mathrm{L}^{*}$ or
against $\mathrm{L}_{\mathrm{gbfs}}$ behave the best. It can be shown (the
proof is provided in Apendix) that since each action has a constant and
positive cost, the heuristic function optimizing $\mathrm{L}^{*}$ optimizes
$\mathrm{L}_{\mathrm{gbfs}}$ as well, but the opposite is not true.

Heuristic functions optimized against $\mathrm{L}_{2}$ estimating cost-to-goal
deliver most of the time the worst results. Heuristic functions optimized
against other ranking losses $\mathrm{L}_{\mathrm{rt}}$ and
$\mathrm{L}_{\mathrm{be}}$ frequently perform well, most of the time better
than $\mathrm{L}_{2}$, but sometimes (Elevators, Sokoban, Sliding puzzle) are
much worse than proposed losses. These experiments further confirm the [57,
18, 59] stating that optimizing the exact estimate of cost-to-goal is
unnecessary, but they also show the advantage of including states off the
optimal path in the loss function. Finally, the method combining policy and
heuristic in GBFS, $\mathrm{L}_{\mathrm{le}},$ [38] performs better that
$\mathrm{L}_{2},$ but worse than the proposed ranking losses. Since both
approaches are derived from the same motivation, we attribute this to
difficulties in training neural networks with two heads. Ranking losses
optimize just the heuristic which is directly used in the search, which seems
to us to be simpler. Additional statistics, namely the average number of
expanded states (Table 6) and the average length of plans (Table 8) are
provided in the Appendix.

##  7 Conclusion

The motivation of this paper stemmed from the observation that even the cost-
to-goal, considered to be an optimal heuristic, fails to guarantee a strictly
optimally efficient search. Since a large body of existing research optimizes
this quantity, we are effectively lost with respect to what should be
optimized. To fill this gap, we have stated the necessary and sufficient
conditions guaranteeing the forward search to be strictly optimally efficient.
These conditions show that the absolute value of the heuristic is not
important, but that the ranking of states in the Open list is what controls
the efficiency. Ranking can be optimized by minimizing the ranking loss
functions, but its concrete implementation needs to correspond to a variant of
the forward search. In case of mismatch, the resulting heuristic can perform
poorly, which has been demonstrated when the heuristic optimized for BGFS
search was used with A* search. The other benefit of ranking losses is that
from the point of view of statistical learning theory, they solve a simpler
problem than ubiquitous regression in estimating the cost-to-goal.

The experimental comparisons on eight problem domains convincingly support the
derived theory. Heuristic functions optimized with respect to ranking loss
$\mathrm{L}^{*}$ instantiated for A* search perform almost always the best
with A* search (except for one case where it was the second best). In the case
of GBFS, heuristic functions optimizing ranking losses
$\mathrm{L}_{\mathrm{gbfs}}$ and $\mathrm{L}^{*}$ instantiated for GBFS and A*
worked the best. This is caused by the fact that the heuristic optimizing
$\mathrm{L}^{*}$ for A* optimizes $\mathrm{L}_{\mathrm{gbfs}}$ for GBFS as
well.

We do not question the existence of a strictly optimally efficient heuristic.
Given our experiments, we believe that if the space of heuristic functions
over which the loss is optimized is sufficiently rich, the result will be
sufficiently close to the optimal for the needs of the search.

###  7.1 Acknowledgements

This work has been supported by project numbers 22-32620S and 22-30043S from
Czech Science Foundation and OP VVV project CZ.02.1.01/0.0/0.0/16_019/0000765
"Research Center for Informatics". This work has also received funding from
the European Union’s Horizon Europe Research and Innovation program under the
grant agreement TUPLES No 101070149.

## References

  * [1] Forest Agostinelli, Stephen McAleer, Alexander Shmakov, and Pierre Baldi.  Solving the rubik’s cube with deep reinforcement learning and search.  Nature Machine Intelligence, 1(8):356–363, 2019.
  * [2] Thomas Anthony, Zheng Tian, and David Barber.  Thinking fast and slow with deep learning and tree search.  In Isabelle Guyon, Ulrike von Luxburg, Samy Bengio, Hanna M. Wallach, Rob Fergus, S. V. N. Vishwanathan, and Roman Garnett, editors, Advances in Neural Information Processing Systems 30: Annual Conference on Neural Information Processing Systems 2017, December 4-9, 2017, Long Beach, CA, USA, pages 5360–5370, 2017.
  * [3] Shahab Jabbari Arfaee, Sandra Zilles, and Robert C Holte.  Learning heuristic functions for large state spaces.  Artificial Intelligence, 175(16-17):2075–2098, 2011.
  * [4] Jeff Bezanson, Alan Edelman, Stefan Karpinski, and Viral B Shah.  Julia: A fresh approach to numerical computing.  SIAM review, 59(1):65–98, 2017.
  * [5] Mohak Bhardwaj, Sanjiban Choudhury, and Sebastian Scherer.  Learning heuristic search via imitation.  In Conference on Robot Learning, pages 271–280. PMLR, 2017.
  * [6] Blai Bonet and Héctor Geffner.  Planning as heuristic search.  Artificial Intelligence. 2001 Jun; 129 (1-2): 5-33., 2001.
  * [7] Vladimir Cherkassky, Xuhui Shao, Filip M Mulier, and Vladimir N Vapnik.  Model complexity control for regression using vc generalization bounds.  IEEE transactions on Neural Networks, 10(5):1075–1089, 1999.
  * [8] Leah Chrestien, Tomáš Pevný, Antonín Komenda, and Stefan Edelkamp.  Heuristic search planning with deep neural network using imitation, attention and curriculum learning.  In PRL-IJCAI workshop, 2021.
  * [9] Lukás Chrpa.  Combining learning techniques for classical planning: Macro-operators and entanglements.  In ICTAI, pages 79–86. IEEE Computer Society, 2010.
  * [10] Rina Dechter and Judea Pearl.  The optimality of A* revisited.  In Michael R. Genesereth, editor, AAAI, pages 95–99. AAAI Press, 1983.
  * [11] Edsger W. Dijkstra.  A note on two problems in connexion with graphs.  Numerische Mathematik, 1:269–271, 1959.
  * [12] Stefan Edelkamp.  Automated creation of pattern database search heuristics.  In Model Checking and Artificial Intelligence, pages 35–50, 2006.
  * [13] Stefan Edelkamp and Stefan Schrödl.  Heuristic Search - Theory and Applications.  Academic Press, 2012.
  * [14] Ariel Felner, Uzi Zahavi, Robert Holte, Jonathan Schaeffer, Nathan Sturtevant, and Zhifu Zhang.  Inconsistent heuristics in theory and practice.  Artificial Intelligence, 175(9):1570–1603, 2011.
  * [15] Dieqiao Feng, Carla P. Gomes, and Bart Selman.  Solving hard AI planning instances using curriculum-driven deep reinforcement learning.  In Christian Bessiere, editor, Proceedings of the Twenty-Ninth International Joint Conference on Artificial Intelligence, IJCAI 2020, pages 2198–2205. ijcai.org, 2020.
  * [16] Patrick Ferber, Malte Helmert, and Jörg Hoffmann.  Neural network heuristics for classical planning: A study of hyperparameter space.  In ECAI, pages 2346–2353. IOS Press, 2020.
  * [17] Santiago Franco, Álvaro Torralba, Levi H. S. Lelis, and Mike Barley.  On creating complementary pattern databases.  In IJCAI, pages 4302–4309, 2017.
  * [18] Caelan Reed Garrett, Leslie Pack Kaelbling, and Tomás Lozano-Pérez.  Learning to rank for synthesizing planning heuristics.  page 3089–3095, 2016.
  * [19] Edward Groshev, Aviv Tamar, Maxwell Goldstein, Siddharth Srivastava, and Pieter Abbeel.  Learning generalized reactive policies using deep neural networks.  2018\.
  * [20] Arthur Guez, Théophane Weber, Ioannis Antonoglou, Karen Simonyan, Oriol Vinyals, Daan Wierstra, Rémi Munos, and David Silver.  Learning to search with mctsnets.  In International conference on machine learning, pages 1822–1831. PMLR, 2018.
  * [21] Peter E Hart, Nils J Nilsson, and Bertram Raphael.  A formal basis for the heuristic determination of minimum cost paths.  IEEE transactions on Systems Science and Cybernetics, 4(2):100–107, 1968.
  * [22] Patrik Haslum, Adi Botea, Malte Helmert, Blai Bonet, and Sven Koenig.  Domain-independent construction of pattern database heuristics for cost-optimal planning.  In AAAI, pages 1007–1012, 2007.
  * [23] Malte Helmert and Gabriele Röger.  How good is almost perfect?  In AAAI, pages 944–949. AAAI Press, 2008.
  * [24] Robert C. Holte and Sandra Zilles.  On the optimal efficiency of cost-algebraic A.  In The Thirty-Third AAAI Conference on Artificial Intelligence, AAAI 2019, pages 2288–2295. AAAI Press, 2019.
  * [25] Michael Katz, Shirin Sohrabi, Horst Samulowitz, and Silvan Sievers.  Delfi: Online planner selection for cost-optimal planning.  IPC-9 planner abstracts, pages 57–64, 2018.
  * [26] Diederik P Kingma and Jimmy Ba.  Adam: A method for stochastic optimization.  arXiv:1412.6980, 2014.
  * [27] Donald E Knuth and Ronald W Moore.  An analysis of alpha-beta pruning.  Artificial intelligence, 6(4):293–326, 1975.
  * [28] Richard E. Korf.  Iterative-deepening-A*: An optimal admissible tree search.  In IJCAI, pages 1034–1036. Morgan Kaufmann, 1985.
  * [29] Richard E Korf.  Macro-operators: A weak method for learning.  Artificial intelligence, 26(1):35–77, 1985.
  * [30] Šimon Mandlík, Matěj Račinskỳ, Viliam Lisỳ, and Tomáš Pevnỳ.  Jsongrinder. jl: automated differentiable neural architecture for embedding arbitrary json data.  Journal of Machine Learning Research, 23(298):1–5, 2022.
  * [31] Drew V. McDermott.  A heuristic estimator for means-ends analysis in planning.  In AIPS, pages 142–149. AAAI, 1996.
  * [32] David Ménager, Dongkyu Choi, Mark Roberts, and David W. Aha.  Learning planning operators from episodic traces.  In 2018 AAAI Spring Symposia. AAAI Press, 2018.
  * [33] Ionut Moraru, Stefan Edelkamp, Santiago Franco, and Moisés Martínez.  Simplifying automated pattern selection for planning with symbolic pattern databases.  In KI, pages 249–263, 2019.
  * [34] Christian Muise.  Collection of pddl files of classical problems.  https://github.com/AI-Planning/classical-domains, 2023.
  * [35] XuanLong Nguyen, Martin J Wainwright, and Michael I Jordan.  On surrogate loss functions and f-divergences.  The Annals of Statistics, pages 876–904, 2009.
  * [36] Laurent Orseau, Levi Lelis, Tor Lattimore, and Théophane Weber.  Single-agent policy tree search with guarantees.  Advances in Neural Information Processing Systems, 31, 2018.
  * [37] Laurent Orseau and Levi HS Lelis.  Implementation of the algorithms described in "policy-guided heuristic search with guarantees" by l. orseau and l. lelis, published at aaai’21.  https://github.com/levilelis/h-levin, 2021.
  * [38] Laurent Orseau and Levi HS Lelis.  Policy-guided heuristic search with guarantees.  In Proceedings of the AAAI Conference on Artificial Intelligence, volume 35, pages 12382–12390, 2021.
  * [39] Arthur L Samuel.  Some studies in machine learning using the game of checkers. ii—recent progress.  IBM Journal of research and development, 11(6):601–617, 1967.
  * [40] Max-Philipp B. Schrader.  gym-sokoban.  https://github.com/mpSchrader/gym-sokoban, 2018.
  * [41] Jendrik Seipp, Florian Pommerening, and Malte Helmert.  New optimization functions for potential heuristics.  In ICAPS, 2015.
  * [42] Jendrik Seipp and Silvan Sievers.  Fast downward stone soup 2023.
  * [43] Jendrik Seipp, Álvaro Torralba, and Jörg Hoffmann.  PDDL generators.  https://doi.org/10.5281/zenodo.6382173, 2022.
  * [44] William Shen.  Strips-hgn.  https://github.com/williamshen-nz/STRIPS-HGN, 2021.
  * [45] William Shen, Felipe Trevizan, and Sylvie Thiébaux.  Learning domain-independent planning heuristics with hypergraph networks.  In ICAPS, volume 30, pages 574–584, 2020.
  * [46] David Silver, Julian Schrittwieser, Karen Simonyan, Ioannis Antonoglou, Aja Huang, Arthur Guez, Thomas Hubert, Lucas Baker, Matthew Lai, Adrian Bolton, et al.  Mastering the game of go without human knowledge.  nature, 550(7676):354–359, 2017.
  * [47] Kristian Smith.  Maze creator / solver.  https://github.com/ravenkls/Maze-Generator-and-Solver, 2022.
  * [48] Gustav Sourek, Vojtech Aschenbrenner, Filip Zelezny, Steven Schockaert, and Ondrej Kuzelka.  Lifted relational neural networks: Efficient learning of latent relational structures.  Journal of Artificial Intelligence Research, 62:69–100, 2018.
  * [49] Simon Ståhlberg, Blai Bonet, and Hector Geffner.  Learning general optimal policies with graph neural networks: Expressive power, transparency, and limits.  In Akshat Kumar, Sylvie Thiébaux, Pradeep Varakantham, and William Yeoh, editors, Proceedings of the Thirty-Second International Conference on Automated Planning and Scheduling, ICAPS 2022, Singapore (virtual), June 13-24, 2022, pages 629–637. AAAI Press, 2022.
  * [50] Ingo Steinwart.  Consistency of support vector machines and other regularized kernel classifiers.  IEEE transactions on information theory, 51(1):128–142, 2005.
  * [51] Takeshi Takahashi, He Sun, Dong Tian, and Yebin Wang.  Learning heuristic functions for mobile robot path planning using deep neural networks.  In ICAPS, volume 29, pages 764–772, 2019.
  * [52] Alvaro Torralba, Vidal Alcázar, Daniel Borrajo, Peter Kissmann, and Stefan Edelkamp.  SymBA*: A symbolic bidirectional A* planner.  In International Planning Competition, pages 105–108, 2014.
  * [53] Sam Toyer, Sylvie Thiébaux, Felipe Trevizan, and Lexing Xie.  Asnets: Deep learning for generalised planning.  Journal of Artificial Intelligence Research, 68:1–68, 2020.
  * [54] Vladimir Vapnik.  Principles of risk minimization for learning theory.  Advances in neural information processing systems, 4, 1991.
  * [55] Marin Vlastelica, Michal Rolinek, and Georg Martius.  Neuro-algorithmic policies enable fast combinatorial generalization.  In International Conference on Machine Learning, pages 10575–10585. PMLR, 2021.
  * [56] Xuemei Wang.  Learning planning operators by observation and practice.  In Kristian J. Hammond, editor, AIPS, pages 335–340. AAAI, 1994.
  * [57] Christopher Wilt and Wheeler Ruml.  When does weighted A* fail?  In International Symposium on Combinatorial Search, volume 3, 2012.
  * [58] Christopher Makoto Wilt and Wheeler Ruml.  Effective heuristics for suboptimal best-first search.  J. Artif. Intell. Res., 57:273–306, 2016.
  * [59] Yuehua Xu, Alan Fern, and Sungwook Yoon.  Learning linear ranking functions for beam search with application to planning.  Journal of Machine Learning Research, 10(7), 2009.
  * [60] Ryo Yonetani, Tatsunori Taniai, Mohammadamin Barekatain, Mai Nishimura, and Asako Kanezaki.  Path planning using neural A* search.  In International Conference on Machine Learning, pages 12029–12039. PMLR, 2021.
  * [61] Ziqi Zhang and Florian Geißer.  Extending graph neural networks for generalized stochastic planning.  In Bridging the Gap Between AI Planning and Reinforcement Learning. 2021.
  * [62] Tan Zhi-Xuan.  PDDL. jl: An Extensible Interpreter and Compiler Interface for Fast and Flexible AI Planning.  PhD thesis, Massachusetts Institute of Technology, 2022.
  * [63] Juntang Zhuang, Tommy Tang, Yifan Ding, Sekhar C Tatikonda, Nicha Dvornek, Xenophon Papademetris, and James Duncan.  Adabelief optimizer: Adapting stepsizes by the belief in observed gradients.  Advances in neural information processing systems, 33:18795–18806, 2020.

##  8 Appendix

###  8.1 Speed of convergence against the size of training set:

We use the classical results of statistical learning theory to show that
estimating the cost-to-go converges more slowly with respect to the size of
the training set than estimating the rank.

From Equation (1) in the main text it should be obvious that ranking is in its
essence a classification problem if one considers a pair of states
$(s,s^{\prime})$ as a single sample. For classification/ranking problems, a
following bound on true error rate [54] holds with probability $1-\eta$

${R}_{c}\leq\hat{R}_{c}+\sqrt{\frac{1}{n_{p}}\left[\left(1+\kappa\ln\frac{2n_{p}}{\kappa}\right)-\ln{\eta}\right]},$

where ${R}_{c}$ denotes the true loss (Equation (1) in the main text) and
$\hat{R}_{c}$, its estimate from $n_{p}$ state pairs $(s_{i},s_{j}),$ and
$\kappa$, the Vapnik-Chervonenkis dimension [54] of the hypothesis space.

Optimizing $h(s,\theta)$ with respect to cost-to-goal (Equation (5) in the
main text) is a regression problem for which a different generalization bound
on prediction error [7] holds with probability $1-\eta$

$R_{r}\leq\hat{R}_{r}\left[1-\sqrt{\frac{1}{n_{s}}\left[\kappa\left(1+\ln\frac{n_{s}}{\kappa}\right)-\ln\eta\right]}\right]^{-1}_{+}$

where again $R_{r}$ is the error of the estimator of cost-to-go and
$\hat{R}_{r}$ is its estimate from $n_{s}$ states (Equation (5) in the main
text),222The formulation in [7] contains constants $c$ and $a,$ but authors
argue they can be set to $a=c=1$ and hasubve been therefore omitted here. and
$[x]_{+}$ is a shorthand for $\max\\{0,x\\}.$

From the above, we can see that excess error in the ranking case converges to
zero at a rate $\frac{\sqrt{\ln n_{p}}}{\sqrt{n_{p}}},$ which is slightly
faster than that of regression $\frac{\sqrt{\ln n_{s}}}{\sqrt{n_{s}}-\sqrt{\ln
n_{s}}}.$ But, the number of state pairs in the training set grows
quadratically with the number of states; therefore, the convergence rate of
the ranking problem for the number of states $n_{s}$ can be expressed as
$\frac{\sqrt{2\ln n_{s}}}{n_{s}},$ which would be by at least
$\frac{1}{\sqrt{n_{s}}}$ factor faster than that of regression. We note that
bounds are illustrative, since the independency of samples is in practice
violated, since samples from the same problem-instance are not independent.

###  8.2 Suboptimality of perfect heuristic in GBFS

A10B3C8D40E28434 Figure 2: Problem instance where perfect heuristic is not
strictly optimally efficient with GBFS. Numbers on the edges denote the cost
of action and red numbers next to nodes denote the minimal cost-to-go.

######  Example 2.

While cost-to-goal $h^{*}$ is the best possible heuristic for algorithms like
A* (up to tie-breaking) in terms of nodes being expanded, for GBFS, $h^{*}$
does not necessarily yield optimal solutions.

###### Proof.

A* explores all the nodes with $f(s)<f^{*}(s)=g(s)+h^{*}(s)$ and some with
$f(s)=f^{*}(s)=g(s)+h^{*}(s)$, so nodes can only be saved with optimal
$f^{*}(s)=g(s)+h^{*}(s)$. In fact, when given $h^{*}$as the heuristic, _only_
nodes $s$ with $f(s)=f^{*}(s)$ are expanded. Depending on the strategy of tie-
breaking, the solution path can be found in the minimal number of node
expansions or take significantly longer (e.g., in lower $g$-value first
exploration of the search frontier). Any heuristic other than $h^{*}$ is
either overestimating, and, therefore, may lead to either non-optimal
solutions in A*, or weaker than $h^{*}$, leading to more nodes being expanded.

Even if $h^{*}$ is given, GBFS is not guaranteed to be optimal. Consider the
following graph with five nodes $A,B,C,D,E$, and $w(A,B)=8$, $w(B,E)=3$,
$w(A,C)=2$, $w(C,D)=4$, $w(D,E)=4$, and $h^{*}(A)=10$, $h^{*}(B)=3$,
$h^{*}(C)=8$, $h^{*}(D)=4$, $h^{*}(E)=0$ (see Figure 2), initial node
$s_{0}=A$, goal node $E\in{\cal S^{*}}$. The numbers are the actual costs and
the red numbers are the exact heuristic function. For finding a path from node
$A$ to node $E$, GBFS would return $(A,B,E)$ following the heuristic function.
However, the path $(A,C,D,E)$ has cost $10$ instead of $11$. ∎

###  8.3 Theory

Below, we prove the claim in made in Experimental section stating that if a
heuristic $h$ is strictly optimally efficient for A* search, then it is also
strictly optimally efficient for GBFS.

######  Theorem 2.

Let a heuristic $h$ is a perfect ranking for A* search on a problem instance
$\Gamma=(\langle\mathcal{S},\mathcal{E},w\rangle,s_{0},\mathcal{S}^{*})$ with
a constant non-negative cost of actions $\left(\left(\exists c\geq
0\right)\left(\forall e\in\mathcal{E}\right)\left(w(e)=c\right)\right)$. Then
$h$ is a perfect ranking for GBFS on $\Gamma.$

###### Proof.

Let $\pi=((s_{0},s_{1}),(s_{1},s_{2}),\ldots,(s_{l-1},s_{l}))$ be an optimal
plan such that $\forall i\in\\{1,\ldots,l\\}$ and $\forall
s_{j}\in\left\\{s_{j}\mid\exists(s_{k},s_{j})\in\mathcal{E}\wedge
s_{k}\in\mathcal{S}^{\pi_{:i-1}}\wedge
s_{j}\notin\mathcal{S}^{\pi_{:i}}\right\\}$ we have
$g(s_{j})+h(s_{j})>g(s_{i})+h(s_{i}),$ where $g(s)$ is the distance from
$s_{0}$ to $s$ in a search-tree created by expanding only states on the
optimal path $\pi.$ We want to proof that if all actions have the same
positive costs, then $h(s_{j})>h(s_{i})$ as well.

We carry the proof by induction with respect to the number of expanded states.

At the initialization (step 0) the claim trivially holds as the Open list
contains just a root node and the set of inequalities is empty.

Let’s now make the induction step and assume the theorem holds for the first
$i-1$ step. We divide the proof to two parts. At first, we prove the claim for
$\mathcal{(}O)_{i}\setminus\mathcal{N}(s_{i-1})$ and then we proof the claim
for $\mathcal{N}(s_{i-1}),$ where $\mathcal{N}(s)$ denotes child states of the
state $s.$

1) Assume $s_{j}\in\mathcal{(}O)_{i}\setminus\mathcal{N}(s_{i-1}).$ Since $h$
is strictly optimally efficient for A*, it holds that

$\displaystyle g(s_{j})+h(s_{j})>$ $\displaystyle g(s_{i})+h(s_{i})$
$\displaystyle h(s_{j})>$ $\displaystyle(g(s_{i})-g(s_{j}))+h(s_{i})$
$\displaystyle h(s_{j})>$ $\displaystyle h(s_{i}),$

where the last inequality is true because $g(s_{i})-g(s_{j})\geq 0.$

Assume $(s_{j}\in\mathcal{N}(s_{i-1}))(s_{j}\neq s_{i}).$ Since $h$ is
strictly optimall efficient for A*, it holds that

$g(s_{j})+h(s_{j})>g(s_{i})+h(s_{i}).$ (7)

Since $g(s_{j})=w((s_{i-1},s_{j}))=w((s_{i-1},s_{i}))=g(s_{i}),$ it holds

$h(s_{j})>h(s_{i}),$ (8)

which finishes the proof of the theorem. ∎

###  8.4 Optimally efficient heuristic might not exists for GBFS

A0B1C1D2A11991 Figure 3: Problem instance where optimally efficient heuristic
does not exists for GBFS.

Consider the following graph in Figure 3 with four nodes $A,B,C,D,$, and
$w(A,B)=1$, $w(B,D)=1$, $w(A,C)=1$, $w(A,D)=9$, $w(B,C)=9$, and $h^{*}(A)=0$,
$h^{*}(B)=1$, $h^{*}(C)=1$, $h^{*}(D)=2$ where A is the goal state and D is
the initial state.

We can see that for GBFS, the perfect heuristic does not exist for D. On
expansion, the GBFS algorithm will put A and B to the open list with heuristic
values h(A) = 0 and h(B) = 1. GBFS takes A from the open list and checks if it
is a goal state. Since A is goal state, GBFS terminates returning path (D,A)
as a solution. However, this is not the optimal path as a better (optimal)
solution exists (D, B, A). Since the definition of optimal ranking requires
the inequalities to hold for this optimal path, in GBFS, the perfect heuristic
does not exist for all initial states.

The problem can be fixed if the definition of optimal ranking is changed to
consider two cases in the merit function $f(s)=\alpha g(s)+\beta h(s)$:
$\alpha>0$ and $\beta>0$ and $\alpha=0$ and $\beta>0$ . In the first case, the
optimal ranking should be defined with respect to the "optimal path" (this is
the A*); in the latter case, it should be the path with minimal number of
expansions. With GBFS, the user simply wants to find a solution but not care
about its optimality. With this change, the heuristic function will exist for
Figure 3.

###  8.5 Training set with multiple solution paths with the same length of the
plan

The behavior of the learned heuristic function depends on the composition of
the training set, which is illustrated below on a simple grid problem. The
agent starts at position (4,4) and has to move to the goal position (0,0).
There are no obstacles and the agent can only move one tile down or one tile
left (the effects on his positions are either (-1,0) or (0,-1)), the cost of
the action is one. The number of different solutions path grows exponentially
with the size of the grid and all solution path has the same cost (for problem
of size 5x5, there are 70 solutions). This problem is interesting, since merit
function in A* with optimal heuristic function (cost to goal) is constant,
equal to 8, for all states.

For the size of the grid (4,4), the heuristic is determined by a table with 25
values. Below, these values are determined by minimizing the proposed loss for
A* algorithm with logistic loss surrogate by BFGS (with all zeros as initial
solution). Since the loss function is convex, BFGS finds the solution quickly.

In the first case, the training set contains one solution path ([4, 4], [3,
4], [2, 4], [1, 4], [0, 4], [0, 3], [0, 2], [0, 1], [0, 0]), where the agent
first goes left and then down. The learnt heuristic values h is shown in Table
2.

y/x | 0 | 1 | 2 | 3 | 4
---|---|---|---|---|---
4 | -191.53 | -131.99 | -76.94 | -31.25 | 0.0
3 | 0.0 | 31.25 | 76.94 | 131.99 | 191.53
2 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0
1 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0
0 | 0.0 | 0.0 | 0.0 | 0.0 | 0.0
Table 2: Table showing one solution path.

An A* search using these heuristic values will first always take action moving
the agent to left and then down. When moving down, the agent does not have
another choice. Notice that the heuristic values of many states are not
affected by the optimization, because they are not needed to effectively solve
the problem.

In the second case, the training set contains two solution paths: the first is
in Table 2 and in the second table 3, the agent goes first down and then left.
The learnt heuristic values are shown in 3.

y/x | 0 | 1 | 2 | 3 | 4
---|---|---|---|---|---
4 | -95.76 | -66.0 | -38.47 | 80.14 | 0.0
3 | 0.0 | 15.62 | 38.47 | 131.99 | 80.14
2 | 0.0 | 0.0 | 0.0 | 38.47 | -38.47
1 | 0.0 | 0.0 | 0.0 | 15.62 | -66.0
0 | 0.0 | 0.0 | 0.0 | 0.0 | -95.76
Table 3: Table showing two solution paths.

An A* search will now face tie at state (4,4), since states (3,4) and (4,3)
have the same heuristic value. But the presence of the tie does not affect the
optimal efficiency of the search, as A* will always expand states on one of
the optimal path from the training set.

Finally let’s consider a case, where all 70 possible solutions are in the
training set. The learned heuristic values are shown in Table 4.

y/x | 0 | 1 | 2 | 3 | 4
---|---|---|---|---|---
4 | 3.74 | 92.62 | 134.66 | 163.17 | 0.0
3 | -46.61 | 2.35 | 91.93 | 134.38 | 163.17
2 | -167.96 | -47.7 | 1.94 | 91.93 | 134.66
1 | -210.03 | -168.66 | -47.7 | 2.35 | 92.62
0 | 0.0 | -210.03 | -167.96 | -46.61 | 3.74
Table 4: Table showing multiple solution paths.

Rank of heuristic values "roughly" corresponds to cost to goal. The reason,
why some states are preferred over the others despite their true cost-to-goal
being the same is that they appear in more solution paths. As shown in Table
4, the A* search with learnt heuristic values is strictly optimally efficient.

###  8.6 Baseline Comparison to breadth-first search

The fraction of solved problems for breadth-first search (5s time limit as
used by solvers in the paper) is shown in Table 5.

domain | fraction
---|---
blocks | 0.35
ferry | 0.31
npuzzle | 0.14
spanner | 0.63
elevators | 0.32
Table 5:  Fraction of solved mazes by breadth-first search.

###  8.7 Training Details

For the grid domains, ADAM [26] training algorithm was run for $100\times
20000$ steps for the grid domains.The experiments were conducted in the
Keras-2.4.3 framework with Tensorflow-2.3.1 as the backend. While all solvers
were always executed on the CPU, the training used an NVIDIA Tesla GPU model
V100-SXM2-32GB. Forward search algorithms were given 10 mins to solve each
problem instance.

For the PDDL domains, the training consumed approximately $100$ GPU hours and
evaluation consumed 1000 CPU hours.All training and evaluation were done on
single-core Intel Xeon Silver 4110 CPU 2.10GHz with a memory limit of 128GB.
The training algorithm AdaBelief [63] was allowed to do $10000$ steps on the
CPU. We emphasize though, that the training time does not include the cost of
creating the training set.

|  | A* |  | GBFS
---|---|---|---|---
problem | complx. | $\mathrm{L}^{*}$ | $\mathrm{L}_{\mathrm{gbfs}}$ | $\mathrm{L}_{\mathrm{rt}}$ | $\mathrm{L}_{2}$ | $\mathrm{L}_{\mathrm{be}}$ |  | $\mathrm{L}^{*}$ | $\mathrm{L}_{\mathrm{gbfs}}$ | $\mathrm{L}_{\mathrm{rt}}$ | $\mathrm{L}_{2}$ | $\mathrm{L}_{\mathrm{be}}$ | $\mathrm{L}_{\mathrm{le}}$
blocks |  | 37 | 54 | 27 | 137 | 54 |  | 33 | 28 | 29 | 32 | 32 | 127
ferry |  | 53 | 43 | 36 | 339 | 20 |  | 48 | 34 | 31 | 33 | 20 | 51
npuzzle |  | 294 | 660 | 843 | 1936 | 641 |  | 297 | 311 | 333 | 591 | 272 | 418
spanner |  | 55 | 546 | 53 | 807 | 416 |  | 61 | 56 | 53 | 117 | 65 | 148
elevators |  | 33 | 35 | 52 | 657 | 310 |  | 33 | 33 | 43 | 115 | 198 | 73
Sokoban | 3 boxes | 14 | 14 | 14 | 17 | 14 |  | 14 | 14 | 14 | 17 | 14 | 14
| 4 boxes | 32 | 35 | 37 | 44 | 35 |  | 34 | 32 | 36 | 43 | 35 | 33
| 5 boxes | 61 | 67 | 68 | 72 | 63 |  | 65 | 62 | 66 | 77 | 64 | 61
| 6 boxes | 171 | 179 | 180 | 210 | 177 |  | 175 | 179 | 181 | 214 | 180 | 174
| 7 boxes | 643 | 651 | 653 | 755 | 654 |  | 645 | 641 | 640 | 754 | 655 | 643
Maze w. t. | $50\times 50$ | 34 | 37 | 34 | 41 | 40 |  | 34 | 33 | 35 | 43 | 40 | 35
| $55\times 55$ | 51 | 59 | 52 | 63 | 60 |  | 54 | 52 | 55 | 65 | 61 | 58
| $60\times 60$ | 72 | 78 | 75 | 83 | 78 |  | 73 | 71 | 77 | 89 | 79 | 84
Sliding puzzle | $5\times 5$ | 1521 | 1558 | 1534 | 1559 | 1545 |  | 1524 | 1539 | 1533 | 1556 | 1544 | 1531
| $6\times 6$ | 2322 | 2353 | 2334 | 2439 | 2388 |  | 2321 | 2326 | 2329 | 2334 | 2329 | 2329
| $7\times 7$ | 3343 | 3375 | 3347 | 3431 | 3411 |  | 3421 | 3356 | 3449 | 3512 | 3448 | 3379

Table 6:  Average number of expanded states.

|  |  |  |  | A* |  | GBFS
---|---|---|---|---|---|---|---
problem | complx. | SBA* | FDSS |  | $\mathrm{L}^{*}$ | $\mathrm{L}_{\mathrm{gbfs}}$ | $\mathrm{L}_{\mathrm{rt}}$ | $\mathrm{L}_{2}$ | $\mathrm{L}_{\mathrm{be}}$ |  | $\mathrm{L}^{*}$ | $\mathrm{L}_{\mathrm{gbfs}}$ | $\mathrm{L}_{\mathrm{rt}}$ | $\mathrm{L}_{2}$ | $\mathrm{L}_{\mathrm{be}}$ | $\mathrm{L}_{\mathrm{le}}$
blocks |  | 100 | 100 |  | $100\pm 0$ | $100\pm 1$ | $100\pm 0$ | $99\pm 1$ | $100\pm 0$ |  | $100\pm 0$ | $100\pm 0$ | $100\pm 0$ | $100\pm 0$ | $100\pm 0$ | $99\pm 1$
ferry |  | 100 | 95 |  | $98\pm 3$ | $98\pm 3$ | $100\pm 0$ | $92\pm 8$ | $100\pm 0$ |  | $98\pm 3$ | $100\pm 0$ | $100\pm 0$ | $100\pm 0$ | $98\pm 4$ | $98\pm 3$
npuzzle |  | 100 | 89 |  | $89\pm 1$ | $87\pm 2$ | $88\pm 0$ | $83\pm 4$ | $89\pm 1$ |  | $92\pm 5$ | $89\pm 1$ | $89\pm 1$ | $89\pm 1$ | $92\pm 7$ | $88\pm 3$
spanner |  | 100 | 90 |  | $100\pm 0$ | $89\pm 2$ | $100\pm 0$ | $84\pm 6$ | $92\pm 6$ |  | $100\pm 0$ | $100\pm 0$ | $100\pm 0$ | $100\pm 0$ | $100\pm 0$ | $100\pm 0$
elevators |  | 100 | 37 |  | $91\pm 2$ | $85\pm 10$ | $75\pm 8$ | $36\pm 6$ | $66\pm 3$ |  | $92\pm 3$ | $85\pm 11$ | $79\pm 8$ | $76\pm 3$ | $67\pm 4$ | $58\pm 11$
Sokoban | 3 boxes | 100 | 100 |  | $99\pm 0$ | $98\pm 0$ | $96\pm 0$ | $97\pm 0$ | $92\pm 0$ |  | $98\pm 0$ | $100\pm 0$ | $94\pm 0$ | $95\pm 0$ | $92\pm 0$ | $98\pm 0$
| 4 boxes | 100 | 81 |  | $89\pm 2$ | $89\pm 1$ | $85\pm 2$ | $81\pm 0$ | $82\pm 3$ |  | $87\pm 2$ | $91\pm 1$ | $84\pm 4$ | $83\pm 3$ | $84\pm 3$ | $84\pm 3$
| 5 boxes | 97 | 67 |  | $80\pm 1$ | $75\pm 2$ | $72\pm 2$ | $72\pm 1$ | $73\pm 0$ |  | $78\pm 1$ | $77\pm 1$ | $74\pm 3$ | $72\pm 3$ | $72\pm 2$ | $73\pm 1$
| 6 boxes | 55 | 49 |  | $76\pm 2$ | $69\pm 1$ | $59\pm 3$ | $51\pm 0$ | $53\pm 3$ |  | $73\pm 1$ | $71\pm 1$ | $56\pm 0$ | $512$ | $54\pm 0$ | $64\pm 1$
| 7 boxes | 46 | 31 |  | $55\pm 4$ | $49\pm 3$ | $47\pm 4$ | $42\pm 3$ | $45\pm 5$ |  | $51\pm 3$ | $49\pm 3$ | $48\pm 3$ | $43\pm 2$ | $45\pm 2$ | $49\pm 3$
Maze w. t. | $50\times 50$ | 92 | 75 |  | $92\pm 0$ | $91\pm 0$ | $88\pm 1$ | $87\pm 1$ | $87\pm 0$ |  | $89\pm 0$ | $90\pm 1$ | $89\pm 1$ | $84\pm 0$ | $85\pm 1$ | $89\pm 0$
| $55\times 55$ | 52 | 50 |  | $78\pm 1$ | $75\pm 4$ | $73\pm 1$ | $72\pm 3$ | $74\pm 3$ |  | $74\pm 0$ | $75\pm 3$ | $74\pm 2$ | $72\pm 3$ | $75\pm 3$ | $74\pm 2$
| $60\times 60$ | 0 | 0 |  | $49\pm 1$ | $37\pm 0$ | $35\pm 3$ | $32\pm 2$ | $31\pm 1$ |  | $42\pm 3$ | $48\pm 1$ | $36\pm 3$ | $34\pm 3$ | $32\pm 2$ | $42\pm 0$
Sliding puzzle | $5\times 5$ | - | 1 |  | $88\pm 1$ | $83\pm 3$ | $84\pm 3$ | $80\pm 2$ | $82\pm 0$ |  | $86\pm 3$ | $87\pm 1$ | $84\pm 0$ | $84\pm 1$ | $84\pm 1$ | $85\pm 1$
| $6\times 6$ | - | - |  | $51\pm 2$ | $48\pm 3$ | $49\pm 2$ | $45\pm 4$ | $46\pm 3$ |  | $47\pm 3$ | $49\pm 3$ | $45\pm 2$ | $43\pm 3$ | $66\pm 2$ | $48\pm 3$
| $7\times 7$ | - | - |  | $39\pm 3$ | $35\pm 3$ | $36\pm 3$ | $32\pm 3$ | $34\pm 3$ |  | $35\pm 3$ | $36\pm 3$ | $35\pm 3$ | $32\pm 3$ | $34\pm 3$ | $35\pm 3$
Table 7: Average fraction of solved problem instances in percent with standard
deviation. SBA* and FDSS denotes Fast Downward Stone Soup, They are domain
independent planners.

|  | A* |  | GBFS
---|---|---|---|---
problem | complx. | $\mathrm{L}^{*}$ | $\mathrm{L}_{\mathrm{gbfs}}$ | $\mathrm{L}_{\mathrm{rt}}$ | $\mathrm{L}_{2}$ | $\mathrm{L}_{\mathrm{be}}$ |  | $\mathrm{L}^{*}$ | $\mathrm{L}_{\mathrm{gbfs}}$ | $\mathrm{L}_{\mathrm{rt}}$ | $\mathrm{L}_{2}$ | $\mathrm{L}_{\mathrm{be}}$ | $\mathrm{L}_{\mathrm{le}}$
blocks |  | 21.6 | 22.7 | 21.4 | 22.5 | 22.9 |  | 21.8 | 22.7 | 22.1 | 22.5 | 23.1 | 34.1
ferry |  | 16.8 | 16.9 | 16.7 | 16.8 | 16.8 |  | 16.8 | 16.9 | 16.8 | 16.8 | 16.8 | 16.9
npuzzle |  | 19.2 | 23.5 | 18.2 | 17.5 | 17.3 |  | 37.6 | 36.4 | 39.6 | 35.1 | 34.1 | 124.6
spanner |  | 49.1 | 49.2 | 49.2 | 49.2 | 49.1 |  | 49.0 | 49.2 | 49.1 | 49.2 | 49.3 | 49.1
elevators |  | 16.1 | 16.4 | 14.7 | 17.9 | 15.8 |  | 19.4 | 16.6 | 15.9 | 18.0 | 16.1 | 21.4
Sokoban | 3 boxes | 13.1 | 13.3 | 13.8 | 13.2 | 13.4 |  | 13.1 | 13.2 | 13.2 | 13.5 | 13.2 | 13.2
| 4 boxes | 15.4 | 15.2 | 16.1 | 15.7 | 16.8 |  | 16.7 | 15.1 | 16.1 | 15.6 | 16.3 | 15.8
| 5 boxes | 20.1 | 19.2 | 21.3 | 22.1 | 20.9 |  | 19.8 | 20.4 | 21.3 | 22.1 | 19.9 | 20.1
| 6 boxes | 29.6 | 29.3 | 27.7 | 28.2 | 26.8 |  | 28.3 | 28.1 | 29.6 | 29.4 | 29.9 | 30.1
| 7 boxes | 31.9 | 31.4 | 35.4 | 33.1 | 34.1 |  | 30.1 | 33.3 | 35.2 | 32.7 | 34.0 | 35.5
Maze w. t. | $50\times 50$ | 24.1 | 25.3 | 25.1 | 24.3 | 24.3 |  | 24.3 | 24.5 | 25.4 | 24.3 | 25.4 | 24.7
| $55\times 55$ | 34.1 | 33.2 | 35.0 | 34.2 | 33.9 |  | 33.2 | 33.1 | 34.6 | 34.6 | 36.5 | 36.3
| $60\times 60$ | 41.2 | 42.9 | 41.4 | 43.2 | 42.8 |  | 42.1 | 43.6 | 44.2 | 45.2 | 45.3 | 45.1
Sliding puzzle | $5\times 5$ | 150.1 | 153.7 | 155.2 | 154.6 | 154.5 |  | 154.5 | 153.5 | 153.9 | 155.0 | 151.1 | 152.1
| $6\times 6$ | 252.3 | 254.2 | 253.8 | 254.8 | 254.0 |  | 255.9 | 256.4 | 255.3 | 256.2 | 254.9 | 256.3
| $7\times 7$ | 321.1 | 324.1 | 322.2 | 324.3 | 320.4 |  | 322.9 | 324.1 | 323.4 | 327.1 | 324.6 | 323.7

Table 8:  Average number of length of the solution. The average is computed
only over problem instances solved by all $11$ variants of forward search.


# Simulating gauge theories on Lefschetz thimbles

Jan M. Pawlowski
Institut für Theoretische Physik, Universität Heidelberg,
Philosophenweg 16, D-69120 Heidelberg, Germany.
E-mail<EMAIL_ADDRESS>Manuel Scherzer
Institut für Theoretische Physik, Universität Heidelberg,
Philosophenweg 16, D-69120 Heidelberg, Germany.
E-mail<EMAIL_ADDRESS>Christian Schmidt
Fakultät für Physik, Universität Bielefeld, Postfach 100131, D-33501
Bielefeld, Germany.
E-mail<EMAIL_ADDRESS>Felix P.G. Ziegler
Institut for Matematik og Datalogi, Syddansk Universitet,
Campusvej 55, 5230 Odense, Danmark.
E-mail<EMAIL_ADDRESS>
Fakultät für Physik, Universität Bielefeld, Postfach 100131, D-33501
Bielefeld, Germany.
E-mail The authors C. Schmidt and F. Ziesché acknowledge support by the
Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) through the
CRC-TR 211 ’Strong-interaction matter under extreme conditions’ project number
315477589 TRR 211. This work is further supported by EMMI, the BMBF grant
05P18VHFCA and is part of and supported by the DFG Collaborative Research
Centre SFB 1225 (ISOQUANT) as well as by DFG under Germany’s Excellence
Strategy EXC-2181/1-390900948 (the Heidelberg Excellence Cluster STRUCTURES).
M. Scherzer acknowledges support from DFG under grant STA 283/16-2. F. P. G.
Ziegler is supported by Heidelberg University. Additionally we thank Andrei
Alexandru, I.O. Stamatescu and Alexander Lindemeier for helpful discussions.

###### Abstract:

Lefschetz thimbles have been proposed recently as a possible solution to the
complex action problem (sign problem) in Monte Carlo simulations. Here we
discuss pure abelian gauge theory with a complex coupling $\beta$ and apply
the concept of Generalized Lefschetz thimbles. We propose to simulate the
theory on the union of the tangential manifolds to the thimbles. We construct
a local Metropolis-type algorithm, that is constrained to a specific
tangential manifold. We also discuss how, starting from this result,
successive subleading tangential manifolds can be taken into account via a
reweighting approach. We demonstrate the algorithm on $U(1)$ gauge theory in
1+1 dimensions and investigate the residual sign problem.

##  1 Generalized Lefschetz thimbles

The numerical sign problem plagues many theories from being simulated at
certain parameters with conventional Monte-Carlo techniques [1]. Since this
problem is representation dependent [2], there are different ways to alleviate
or even eliminate it such as dual representations, specialised Monte-Carlo
techniques (density of states) or methods based on extending the configuration
space (Complex Langevin [3], Lefschetz thimbles [4]). If we want to calculate
a multi-dimensional integral

$\int\mathrm{d}\phi~{}e^{-S[\phi]}\mathcal{O}(\phi)$

with compact real fields $\phi$, where the integrand is holomorphic, we can
complexify the fields $\phi$ and according to Cauchy’s theorem choose a
submanifold in complexified space homotopic to real subspace to get the same
result. Dealing with the sign problem means in this case choosing a
submanifold, where the fluctuations of the phase of the integrand are reduced
(see e.g.[5]). Since the phase itself is non-holomorphic, it depends on the
integration manifold. In other words, the sign problem is representation
dependent.
Lefschetz thimbles are originally a basis of homology classes for complex
varieties. In our case, their representatives can be chosen to keep the phase
of our integrand (or just $e^{-iS_{I}}$) constant. Being a basis, one can
build a submanifold homotopic to the original integration space from these.
A Lefschetz thimble is generally defined to be the union of flowlines
generated by the steepest descent equation of the action $S$

$\frac{\mathrm{d}\phi}{\mathrm{d}t}=-\left(\frac{\delta
S}{\delta\phi}\right)^{*}$ (1)

which end in a non-degenerate critical point $\phi_{\sigma}$. Since we are
looking at a gauge theory, every critical point is naturally degenerate and
the classical Picard-Lefschetz theory does not apply. But still we have the
concept of Generalized Lefschetz thimbles [6], which was also outlined for QCD
[4]: Instead of critical points, we have seperate critical manifolds spanned
by the zero modes of the action. Complementary to the zero modes on the
critical manifold are the Takagi and Anti-Takagi modes, which classically span
the tangent spaces of the thimble and the anti-thimble, if there are no zero
modes. In the case of compact gauge groups, one can choose a compact
submanifold of the critical manifold, whose dimension plus the number of the
Takagis gives the real dimension of the original integration space. If the
degeneracy comes from the gauge degrees of freedom, this can be spanned by the
real gauge transformations and is called _gauge orbit_. In the following, we
will use this freedom to construct a local update algorithm on the tangent
space.

##  2 (1+1)d-U(1) Lattice Yang-Mills theory

We discretize the Yang-Mills action on a two dimensional Euclidean space time
lattice $\Lambda$, i.e. we consider Wilson’s plaquette action [7]

$S=\beta\sum_{x}\left\\{1-\frac{1}{2}\left(P_{01}(x)+P^{-1}_{01}(x)\right)\right\\},$
(2)

where $P_{01}(x)=U_{0}(x)U_{1}(x+\hat{0})U^{-1}_{0}(x+\hat{1})U^{-1}_{1}(x)$
denotes the elementary plaquette in the ($0,1$)-plane at site $x$. The link
variables $U_{\mu}(x)$ are elements of the gauge group, which we consider to
be $\text{U}(1)$. Luckily, we have a formal solution for the partition sum

$Z=\int\mathrm{d}U\exp\left({-S[U]}\right)=\sum_{n=-\infty}^{+\infty}\left[I_{n}(\beta)\right]^{V}$
(3)

being a series in modified Bessel functions $I_{n}(\beta)$, where $V$ is the
number of plaquettes [8, 9]. The sign problem is introduced by generalising to
complex couplings $\beta$. This corresponds physically to interpolating
between imaginary and real time: Imaginary $\beta$ corresponds to the real-
time case, while real $\beta$ is the imaginary time case. In principal this
allows to study thermal physics using paths in the complex time plane which
for example approximate the Schwinger-Keldysh contour [10, 11]. The critical
manifolds obtained by setting the gradient of the action to zero can be
described by the following relations

$P_{01}(x)=P_{01}(x-\hat{\mu})\>\>\text{or}\>\>P_{01}(x)=-P^{-1}_{01}(x-\hat{\mu}),\;\;\mu\in\\{0,1\\}$
(4)

between neighboring plaquettes variables. Additionally, we get the constraint

$\prod_{x}P_{01}(x)=1$ (5)

by periodic boundary conditions. This still leaves us a great amount of
critical manifolds, since e.g. the constraint gives us the possibility for
every selection of $V$-th roots of unity. We have to take a good selection of
critical manifolds, in that sense that the overall manifold, we are going to
create is homotopic to $[\mathrm{U}(1)]^{2V}$ modulo copies. For inspiration,
we look at the action for one plaquette

$S=-\beta/2(P+P^{-1}),$ (6)

where we omitted a volume factor. We have naturally two critical points $P=\pm
1$ with the respective imaginary parts of the action being
$S_{I}=\mp\beta_{I}$ and their attached thimbles $\mathcal{J}_{0}$,
$\mathcal{J}_{1}$ , see Figure 1.

Figure 1: This is the thimble manifold for $\beta=(1+3i)$ and their
corresponding plot of $e^{-S_{R}}$, which is proportional to the probability
density.

We can take these possibilities for values of $P_{01}(x)$. Then we are
restricted by (5) to configurations, where an even number of plaquettes can be
$-1$. These configurations clearly fulfill (4).

###  2.1 Critical manifolds, their tangent spaces and local updates

These are critical manifolds, since our degrees of freedom are still links. We
have to compute their Takagi and Anti-Takagis, i.e. vectors, which are
solution to the equation

$H^{*}\Delta z^{*}=\lambda\ \Delta z\>\>\text{with}\>\>\lambda\in\mathbb{R},$
(7)

where $H$ is the Hessian of the action. Positive $\lambda$ refer to thimble
directions, negative $\lambda$ to anti-thimble directions and $\lambda=0$
refers to zero-modes, which come e.g. from the gauge degrees of freedom (see
e.g. [4, 12]). For our selection of critical manifolds $P_{i}=\pm 1$ the
Hessian splits into a real matrix with a complex prefactor $H=\beta M$, whose
eigenvectors $v$ and eigenvalues $\alpha$ can be computed. For $\alpha\neq 0$
we get for the Takagis

$\Delta z=\sqrt{\frac{\pm\mathrm{sign}(\alpha)\beta^{*}}{|\beta|}}v,$ (8)

where the sign denotes, if it is a Takagi or Anti-Takagi vector.
We observe, that our Hessian is independent from the actual configuration in
the critical manifold. Therefore, we can deduce that the projection of the
subspace spanned by its zero modes in the Lie algebra is the critical manifold
itself.

$\left\\{U_{\mu}(x)^{\mathrm{crit}}=U_{\mu}(x)^{\mathrm{crit},0}\exp\left({i\sum_{k=1}^{\\#(\alpha=0)}c^{k}v^{x,\mu}_{k}(\alpha=0)}\right)\;\;|c^{k}\in\mathbb{C},\;\mu\in\\{0,1\\},\;x\in\Lambda\right\\}$

Normally, we have to choose a gauge orbit of real dimension by allowing only
$c^{k}\in\mathbb{R}$, where we can span the thimble using the Takagi
directions. For our main critical manifold $\\{P_{i}=1\;\forall i\\}$, the
complex prefactor in (8) is the same for all Takagi modes. By tilting the real
zero modes with the same factor, we differ from a gauge orbit, but get a
manifold, which still gives the same expectation values for observables
invariant under the zero modes. The set of tilted zero modes and Takagis is
equivalent to the tilted unit basis. This gives us the possibility to have a
local update algorithm on the tangent space. This is naturally computationally
far less demanding than sampling on the thimbles themselves. Restricting
ourselves to the main critical manifold, the Jacobian is constant and drops
out for expectation values. Another reason, why this could be feasible is the
closeness of the local tangent space to the thimble in the one-plaquette model
(see Figure 2) in their important regions.

Figure 2: Left: Comparison of the thimbles and their tangent spaces for
$\beta=3+3i$. Middle: Thimble hierarchy depending on $\beta_{R}$ in the
approximation. Right: Configurations on the tangent space of two plaquettes
related to one link.

###  2.2 Hierarchy of critical manifolds

Our selection of critical manifolds has a natural hierarchy, which is
reflected by the values of the action
$S=4k\beta,\>\>\>k=0,\ldots,\lfloor{\frac{V}{2}}\rfloor,$ where $k$ is the
number of turned plaquette pairs, which are $-1$. With increasing $\beta_{R}$
their weights in the partition sum decrease exponentially. Note that on the
thimble the real part of the action becomes minimal at the critical manifold.
Consequently for pure imaginary $\beta$, every thimble contributes equally.
Otherwise we get a close result by taking into account only a few thimbles,
since the others are exponentially suppressed. To get a hint on how strong
this is the case, we approximate our model by just taking the leading order
contribution of our formal expansion in modified Bessel functions (3)

$Z=\left[\int_{\mathrm{U}(1)}\mathrm{d}P\
e^{\beta/2(P+P^{-1})}\right]^{V}=\left[I_{0}(\beta)\right]^{V}.$ (9)

Physically this corresponds to removing periodic boundary conditions. Since
this is the One-plaquette model to the power of the volume, we can expand this
in term of its thimbles $\mathcal{J}_{0}$, $\mathcal{J}_{1}$

$Z=\left[\int_{\mathcal{J}_{0}}\mathrm{d}P\
e^{\beta/2(P+P^{-1})}+\int_{\mathcal{J}_{1}}\mathrm{d}P\
e^{\beta/2(P+P^{-1})}\right]^{V}\\\
=:\left[Z_{0}+Z_{1}\right]^{V}=\sum_{k=0}^{V}\begin{pmatrix}V\\\
k\end{pmatrix}Z_{0}^{V-k}\ Z_{1}^{k}.$ (10)

Using these, we can calculate approximate values for our observables and their
dependence on the thimble hierarchy (see Figure 2).

##  3 Simulation and comparison

###  3.1 Algorithm

Since our critical manifolds and tangent spaces depend on the values of the
plaquettes, we use them to confine the regions where we sample on them. As one
can observe in Figure 2 the two tangent spaces for one plaquette intersect and
we can glue them together. This union is homotopic to the orginal
$\mathrm{U}(1)$ group. We will use these intersections to limit the region of
tangent space we explore. All in all we do a local Metropolis update on the
links, where we control the values of the associated plaquettes, preventing
them from wandering off the designated main tangent space shown on the right
hand side of Figure 2. This is guaranteed by defining that candidates which
would land beyond the edges have probability zero and will be rejected by the
Metropolis.

###  3.2 Results for the main tangent space

Figure 3:  First row: plaquette for $\beta_{R}=2$ with reweighting compared to
simulation on the main tangential manifold on a 4 x 4 lattice. Second row:
volume dependence of the formal solution and the average sign of the main
tangential manifold. Third row: plaquette for $\beta_{R}=2$ with reweighting
compared to simulation on the main tangential manifold on a 16 x 16 lattice.

We calculate the expectation value of the average plaquette $<1/2(P+P^{-1})>$.
Since $\beta$ is complex, this has a real and an imaginary part. We first note
that simulating on the main tangent space alleviates the sign problem in
comparison with normal reweighting. Since we only take the main tangent space
into account, the results can be considered ’right’ only for large enough beta
and sufficiently large volumes. Another thing we note is that for constant
$\beta_{R}$, the average sign $|<e^{-iS_{I}}>|$ has a minimum on the
$\beta_{I}$ range for every observed volume. So it gets better again for
higher $\beta_{I}$.
Looking at this dip at $\beta=2+1.4i$, we look at the volume dependence of
$|<e^{-iS_{I}}>|$ using 4-by-4 to 16-by-16 lattices. This seems strictly
exponential like it is predicted for the full theory by considering the free
energy. But even for a 16-by-16 lattice, where reweighting fails very early
for small $\beta_{I}$, we can get a large range where the sign problem is
quite mild on the main tangent space. Since higher orders in the modified
Bessel function expansion (3) are suppressed by the volume, the leading order
approximation coincides stronger with the simulation result.

###  3.3 Reweighting onto the other tangent spaces

For taking another tangential manifold $\tau_{1}$ into account, e.g. the one
related to the critical manifold, where two plaquettes are $-1$, which would
be the next one in the hierarchy, we need the ratio $Z_{1}/Z_{0}$ of the
partition sums, since we have

$<\mathcal{O}>_{\tau_{1}\cup\tau_{2}}=\frac{\int_{\tau_{1}}\mathrm{d}U~{}\mathcal{O}[U]e^{-S[U]}+\int_{\tau_{2}}\mathrm{d}U~{}\mathcal{O}[U]e^{-S[U]}}{\int_{\tau_{1}}\mathrm{d}U~{}e^{-S[U]}+\int_{\tau_{2}}\mathrm{d}U~{}e^{-S[U]}}=\frac{<\mathcal{O}>_{\tau_{0}}+(Z_{1}/Z_{0})<\mathcal{O}>_{\tau_{2}}}{1+(Z_{1}/Z_{0})}.$
(11)

Calculating this, we can follow the method proposed in [13]: Suppose, we have
a mapping

$f\colon\tau_{0}\longrightarrow\;\;\tau_{1}.$

Then we can write

$\frac{Z_{1}}{Z_{0}}=\frac{\int_{\tau_{0}}\mathrm{d}Ue^{-S[f(U)]+S[U]}\det[\mathrm{d}f]e^{-S[U]}}{\int_{\tau_{0}}\mathrm{d}Ue^{-S[U]}}=<e^{-S\circ
f+S}\det[\mathrm{d}f]>_{0}$ (12)

The problem here is to find a suitable $f$, which we discuss in our upcoming
paper [14]. But we can already say, that since we consider only tangent
spaces, $f$ is linear and therefore the Jacobian $\det[\mathrm{d}f]$ is be a
constant factor.

##  4 Summary and Outlook

We have simulated a two dimensional U(1) lattice gauge theory on the tangent
space of its main thimble. Hereby, the sign problem is drastically reduced,
while the computational complexity has stayed the same as for a local
Metropolis update. We proposed taking into account subleading thimbles by a
reweighting approach, which is under construction and discussion. We will
pursue the technique in the future by looking at other gauge groups including
also fermionic determinants and a chemical potential.

## References

  * [1] P. de Forcrand,  _Proceedings, 27th International Symposium on Lattice field theory (Lattice 2009): Beijing, P.R. China, July 26-31, 2009_ , PoS LAT2009, 010 (2009), arXiv:1005.0539 [hep-lat]
  * [2] C. Gattringer and K. Langfeld, Int. J. Mod. Phys. A31, 1643007 (2016), arXiv:1603.09517 [hep-lat]
  * [3] G. Aarts and I.-O. Stamatescu, JHEP 09, 018 (2008), arXiv:0807.1597 [hep-lat]
  * [4] M. Cristoforetti, et al. (AuroraScience), Phys. Rev. D86, 074506 (2012), arXiv:1205.3996 [hep-lat]
  * [5] Y. Mori, et al., Phys. Rev. D96, 111501 (2017), arXiv:1705.05605 [hep-lat]
  * [6] E. Witten, AMS/IP Stud. Adv. Math. 50, 347 (2011), arXiv:1001.2933 [hep-th]
  * [7] K. G. Wilson, Phys. Rev. D10, 2445 (1974)
  * [8] R. Balian, et al., Phys. Rev. D10, 3376 (1974), [,74(1974)]
  * [9] A. A. Migdal, Sov. Phys. JETP 42, 413 (1975), [,114(1975)]
  * [10] J. Berges, et al., Phys. Rev. D75, 045007 (2007), arXiv:hep-lat/0609058 [hep-lat]
  * [11] H. Makino, et al., Phys. Rev. D92, 085020 (2015), arXiv:1503.00417 [hep-lat]
  * [12] A. Alexandru, et al., Phys. Rev. D93, 014504 (2016), arXiv:1510.03258 [hep-lat]
  * [13] S. Bluecher, et al., SciPost Phys. 5, 044 (2018), arXiv:1803.08418 [hep-lat]
  * [14] J. Pawlowski, et al., “A local update for pure Yang-Mills theories with a complex coupling $\beta$,” In preparation