Latent Heuristic Search: Continuous Optimization for Automated Algorithm Design

# Latent Heuristic Search: Continuous Optimization for Automated Algorithm Design
Source: [https://arxiv.org/html/2605.17137](https://arxiv.org/html/2605.17137)
11institutetext:Huawei Technologies Canada, Burnaby, Canada###### Abstract

The integration of Large Language Models \(LLMs\) into evolutionary frameworks has established a new paradigm for automated heuristic discovery\. Despite their promise, these methods typically search in the discrete space of program syntax, relying on stochastic sampling to navigate a highly non\-convex optimization landscape\. This work proposes a continuous heuristic discovery framework that shifts optimization to a learned latent manifold\. We employ an encoder to map discrete programs into continuous embeddings and train a differentiable surrogate model to predict performance, enabling gradient\-based search\. To regularize the optimization trajectory, an invertible normalizing flow maps these embeddings to a structured Gaussian prior, where we perform gradient ascent\. The resulting optimized latent vectors are projected through a learned mapper into soft prompts, which condition a frozen LLM to synthesize novel executable heuristics\. We evaluate the proposed method on the Traveling Salesman Problem \(TSP\), the Capacitated Vehicle Routing Problem \(CVRP\), the Knapsack Problem \(KSP\), and Online Bin Packing \(OBP\)\. Empirical results demonstrate that continuous latent\-space optimization achieves performance competitive with state\-of\-the\-art discrete evolutionary baselines while offering a complementary methodological alternative for automated algorithm design\. The implementation code is available at[https://github\.com/cheikh025/LHS](https://github.com/cheikh025/LHS)\.

## 1Introduction

Heuristic algorithms are central to combinatorial optimization, enabling high\-quality solutions for problems such as the Traveling Salesman Problem \(TSP\)\[[10](https://arxiv.org/html/2605.17137#bib.bib24)\], the Capacitated Vehicle Routing Problem \(CVRP\)\[[21](https://arxiv.org/html/2605.17137#bib.bib25)\], and online bin packing \(OBP\)\[[7](https://arxiv.org/html/2605.17137#bib.bib26)\]under practical time constraints\. Despite their importance, designing strong heuristics remains largely expert\-driven, requiring iterative refinement, deep domain knowledge, and extensive empirical testing\.

Recent work has shown that large language models \(LLMs\) can support algorithm discovery when embedded in a rigorous*search\-and\-evaluate*loop\. FunSearch couples LLM\-based code generation with evolutionary sampling and execution\-based evaluation, demonstrating improvements on mathematical construction tasks and on classic heuristic settings such as online bin packing\[[19](https://arxiv.org/html/2605.17137#bib.bib1)\]\. Building on this idea, Evolution of Heuristics \(EoH\) introduces a framework that co\-evolves high\-level natural\-language “thoughts” describing heuristic strategies alongside executable code, with the aim of improving search efficiency and diversity\[[12](https://arxiv.org/html/2605.17137#bib.bib2)\]\. ReEvo further frames LLMs as hyper\-heuristics and uses reflective feedback to guide evolutionary search, providing an explicit mechanism for steering generation toward more effective heuristic behaviors\[[22](https://arxiv.org/html/2605.17137#bib.bib3)\]\. Complementary lines of work explore alternative search operators and exploration strategies, including MCTS\-based exploration\[[23](https://arxiv.org/html/2605.17137#bib.bib10)\], diversity\-driven evolutionary search with Harmony Search\[[2](https://arxiv.org/html/2605.17137#bib.bib15)\], and more general evolutionary coding agents for scientific and algorithmic discovery\[[16](https://arxiv.org/html/2605.17137#bib.bib16)\]\. Collectively, these works highlight a key lesson: LLMs become substantially more reliable for heuristic design when paired with systematic exploration and objective evaluation, rather than one\-shot generation\.

However, many LLM\-based automated heuristic design methods still operate primarily in the*discrete space of programs*, relying on mutations, recombination, and repeated LLM sampling over code tokens\. This often yields a challenging optimization landscape and can require many expensive evaluations to discover consistent improvements\. In parallel, research on*latent program representations*suggests that program search can be made more structured by learning continuous embeddings of code and performing optimization directly in that learned space\. Approaches such as LEAPS learn embeddings and decoders that enable search over latent variables rather than explicit syntactic edits\[[20](https://arxiv.org/html/2605.17137#bib.bib4)\], and subsequent work explores latent program space search as a mechanism for generalization and adaptation\[[14](https://arxiv.org/html/2605.17137#bib.bib5)\]\. Closely related,*program synthesis*has been explicitly framed as continuous optimization under specifications given by input/output examples or test\-based error signals: NPO\[[11](https://arxiv.org/html/2605.17137#bib.bib21)\]trains a neural program autoencoder to embed programs into a continuous space and then applies derivative\-free optimization \(e\.g\., CMA\-ES\) over latents, decoding candidates back to programs for execution\-based scoring\. In contrast, GENESYS\[[15](https://arxiv.org/html/2605.17137#bib.bib22)\]avoids a learned autoencoder by directly parameterizing token choices with continuous variables and similarly optimizes these parameters with CMA\-ES augmented by restart policies\. These results motivate the hypothesis that automated heuristic design may benefit from shifting the primary search process from code token space to a continuous representation space that captures semantic variation among heuristics; however, our setting differs from I/O consistency synthesis in that we search for*high\-performing heuristic programs*under execution\-based benchmarking, rather than programs that satisfy a given set of I/O examples\.

In this work, we propose a framework for discovering improved heuristic algorithms by*optimizing in a learned latent space of heuristic programs*\. Our approach represents candidate programs as continuous embeddings and learns a surrogate ranking model that predicts heuristic performance from these representations\. To regularize the optimization procedure and preserve decodability, we fit an invertible normalizing flow that maps program embeddings to a structured prior space, enabling gradient\-based search in a well\-conditioned space while remaining close to the manifold of plausible heuristics\.

To translate optimized latent codes into executable algorithms, we introduce a decoding pathway in which candidate embeddings are mapped to soft prompts that condition an LLM\-based code generator\. The resulting code is validated for syntactic correctness and executed under the benchmark protocol to obtain its performance score\. In contrast to token\-level evolutionary search, our approach leverages a*continuous, learnable search geometry*over heuristic designs while retaining strict execution\-based evaluation\.

## 2Problem Formulation

Let𝒫\\mathcal\{P\}denote the set of syntactically valid programs, written in a fixed programming language, that implement heuristics for a target combinatorial optimization problem \(e\.g\., TSP, OBP\)\. For anyp∈𝒫p\\in\\mathcal\{P\}, we evaluate performance by executingppon a benchmark set of instancesDbench=\{xj\}j=1MD\_\{\\mathrm\{bench\}\}=\\\{x\_\{j\}\\\}\_\{j=1\}^\{M\}under a standardized protocol \(e\.g\., fixed time limits, feasibility checks, and scoring rules\), which yields instance\-level objective valuesCost​\(p;xj\)∈ℝ\\mathrm\{Cost\}\(p;x\_\{j\}\)\\in\\mathbb\{R\}\.

To unify notation across minimization and maximization tasks, we define

y​\(p\):=1M​∑j=1MCost​\(p;xj\),s​\(p\):=α​y​\(p\),α∈\{\+1,−1\},y\(p\)\\;:=\\;\\frac\{1\}\{M\}\\sum\_\{j=1\}^\{M\}\\mathrm\{Cost\}\(p;x\_\{j\}\),\\qquad s\(p\)\\;:=\\;\\alpha\\,y\(p\),\\quad\\alpha\\in\\\{\+1,\-1\\\},\(1\)whereα=\+1\\alpha=\+1for maximization tasks andα=−1\\alpha=\-1for minimization tasks\. The mappingp↦s​\(p\)p\\mapsto s\(p\)is treated as a black box and is non\-differentiable with respect to program tokens\.

Using this evaluator, we maintain a dataset of scored programs𝒟=\{\(pi,si\)\}i=1N\\mathcal\{D\}=\\\{\(p\_\{i\},s\_\{i\}\)\\\}\_\{i=1\}^\{N\}, wheresi:=s​\(pi\)s\_\{i\}:=s\(p\_\{i\}\)\. The goal is to discover an improved heuristic program by solving

p⋆∈arg⁡maxp∈𝒫⁡s​\(p\),p^\{\\star\}\\in\\arg\\max\_\{p\\in\\mathcal\{P\}\}\\;s\(p\),\(2\)using search in a learned continuous representation space rather than discrete, token\-level edits in program space\.

## 3Method

We propose a framework for heuristic discovery via continuous optimization in a learned latent program space\. Figure[3](https://arxiv.org/html/2605.17137#S3)provides an overview of the proposed pipeline\. Specifically, we \(i\) encode programs into continuous representations, \(ii\) train a differentiable surrogate that predicts program performance from these representations, and \(iii\) perform gradient\-based optimization to identify improved latent candidates, which are subsequently decoded into executable programs by a large language model \(LLM\) and evaluated under the benchmark protocol\.

![[Uncaptioned image]](https://arxiv.org/html/2605.17137v1/LHS_v1.png)### 3\.1Latent Heuristic Search

To discover high\-performing heuristics, we shift the search process from the discrete, non\-differentiable space of code to a continuous latent manifold\. This transition enables the use of continuous optimization techniques, such as gradient\-based updates, in the learned representation space, The complete discovery procedure is summarized in Algorithm[1](https://arxiv.org/html/2605.17137#alg1)\. Each heuristic programp∈𝒫p\\in\\mathcal\{P\}is mapped to a latent vector using a code encoder:

z=E​\(p\),z∈ℝd\.z\\;=\\;E\(p\),\\qquad z\\in\\mathbb\{R\}^\{d\}\.\(3\)Our objective is to find programs with maximal benchmark performance scores​\(p\)s\(p\)\(Section[2](https://arxiv.org/html/2605.17137#S2)\)\. Sinces​\(p\)s\(p\)is obtained by executingppunder the benchmark protocol, it is treated as a black\-box function and does not provide gradients with respect tozz\.

To obtain a differentiable search objective, we learn a surrogate predictorfθ:ℝd→ℝf\_\{\\theta\}:\\mathbb\{R\}^\{d\}\\to\\mathbb\{R\}\. Using the dataset of evaluated programs𝒟=\{\(pi,si\)\}i=1N\\mathcal\{D\}=\\\{\(p\_\{i\},s\_\{i\}\)\\\}\_\{i=1\}^\{N\}, we encode each program aszi=E​\(pi\)z\_\{i\}=E\(p\_\{i\}\)and trainfθf\_\{\\theta\}such that:

fθ​\(zi\)≈si\.f\_\{\\theta\}\(z\_\{i\}\)\\approx s\_\{i\}\.\(4\)
Given a seed latentz\(0\)z^\{\(0\)\}\(e\.g\., the encoding of a reference program\), we generate improved latent candidates by gradient ascent on the surrogate:

z\(t\+1\)=z\(t\)\+η​∇zfθ​\(z\(t\)\),t=0,…,T−1\.z^\{\(t\+1\)\}\\;=\\;z^\{\(t\)\}\+\\eta\\,\\nabla\_\{z\}f\_\{\\theta\}\\\!\\left\(z^\{\(t\)\}\\right\),\\quad t=0,\\ldots,T\-1\.\(5\)
After optimization, the final latent vectorz\(T\)z^\{\(T\)\}must be converted back into executable code\. We projectz\(T\)z^\{\(T\)\}into the embedding space of a code\-generating LLM via a learnable projection functionϕω\\phi\_\{\\omega\}:

h=ϕω​\(z\(T\)\)∈ℝK×e,h\\;=\\;\\phi\_\{\\omega\}\(z^\{\(T\)\}\)\\in\\mathbb\{R\}^\{K\\times e\},\(6\)which producesKKcontinuous conditioning embeddings\. The sequencehhis used as a soft prompt to condition a frozen decoder\-only LLM𝒢\\mathcal\{G\}for code generation using prefix/prompt tuning\. Given a context promptII\(problem description, specification, constraints\), we generate a candidate program:

p^∼𝒢​\(I⊕h\)\.\\hat\{p\}\\sim\\mathcal\{G\}\\\!\\left\(I\\oplus h\\right\)\.\(7\)where⊕\\oplusdenotes the concatenation of the task promptIIand the soft prompthh\. We then validate and executep^\\hat\{p\}to obtain its benchmark scores​\(p^\)s\(\\hat\{p\}\)\.

### 3\.2Training the Latent\-to\-Prompt Mapper

The mapper networkϕω\\phi\_\{\\omega\}serves as a semantic bridge between the compact latent space and the high\-dimensional input space of the decoder\. It projects a program’s latent vectorzzinto a sequence ofKKcontinuous embedding vectors, essentially a soft prefix, that conditions the frozen Large Language Model \(LLM\)\. This approach adopts a prefix\-tuning strategy\[[9](https://arxiv.org/html/2605.17137#bib.bib6)\], allowing us to steer the generation process without modifying the parameters of the LLM\.

For each programpip\_\{i\}in our dataset, we first compute its latent representationzi=E​\(pi\)z\_\{i\}=E\(p\_\{i\}\)using the encoder\. The mapper then transforms this vector into a sequence of soft tokens:

hi=ϕω​\(zi\)∈ℝK×e,h\_\{i\}\\;=\\;\\phi\_\{\\omega\}\(z\_\{i\}\)\\in\\mathbb\{R\}^\{K\\times e\},\(8\)whereeedenotes the embedding dimension of the decoder𝒢\\mathcal\{G\}\.

To ensure thathih\_\{i\}retains the semantic information required to reconstruct the original program, we trainψ\\psiusing a supervised reconstruction objective\. Let the programpip\_\{i\}be tokenized as a sequence\(ti,1,…,ti,Ti\)\(t\_\{i,1\},\\dots,t\_\{i,T\_\{i\}\}\), and letIiI\_\{i\}denote the fixed textual task prompt \(specifying the problem definition, interface, and constraints\)\. The decoder𝒢\\mathcal\{G\}defines an autoregressive probability distribution over the tokens:

P𝒢​\(ti,1:Ti∣Ii,hi\)=∏k=1TiP𝒢​\(ti,k\|ti,<k,Ii⊕hi\)\.P\_\{\\mathcal\{G\}\}\(t\_\{i,1:T\_\{i\}\}\\mid I\_\{i\},h\_\{i\}\)\\;=\\;\\prod\_\{k=1\}^\{T\_\{i\}\}P\_\{\\mathcal\{G\}\}\\\!\\left\(t\_\{i,k\}\\;\\big\|\\;t\_\{i,<k\},\\;I\_\{i\}\\oplus h\_\{i\}\\right\)\.We optimize the parametersω\\omegaby minimizing the negative log\-likelihood of the ground\-truth program tokens under teacher forcing:

ℒmap​\(ω\)=−1N​∑i=1N1Ti​∑k=1Tilog⁡P𝒢​\(ti,k\|ti,<k,Ii⊕ϕω​\(E​\(pi\)\)\)\.\\mathcal\{L\}\_\{\\mathrm\{map\}\}\(\\omega\)\\;=\\;\-\\frac\{1\}\{N\}\\sum\_\{i=1\}^\{N\}\\frac\{1\}\{T\_\{i\}\}\\sum\_\{k=1\}^\{T\_\{i\}\}\\log P\_\{\\mathcal\{G\}\}\\\!\\left\(t\_\{i,k\}\\;\\big\|\\;t\_\{i,<k\},\\;I\_\{i\}\\oplus\\phi\_\{\\omega\}\(E\(p\_\{i\}\)\)\\right\)\.\(9\)By minimizingℒmap\\mathcal\{L\}\_\{\\mathrm\{map\}\}, the mapper learns to translate the geometric coordinates of the latent space into a conditioning signal that reliably elicits the corresponding syntax and logic from the frozen decoder𝒢\\mathcal\{G\}\.

### 3\.3Latent optimization with a normalizing flow

While the encoderEEmaps discrete heuristics to continuous vectors, performing gradient\-based search directly in the raw embedding space𝒵\\mathcal\{Z\}is suboptimal\. Pre\-trained encoder representations are known to suffer from the anisotropy problem: their empirical distribution is typically highly non\-isotropic and may contain low\-density regions where the decoder becomes unreliable\. Direct gradient ascent in𝒵\\mathcal\{Z\}can therefore drift toward out\-of\-distribution latents, yielding invalid or low\-quality generated heuristics\. Prior work in representation learning and controllable generation reports analogous issues, such as anisotropy, asymmetry, and non\-convex support in latent spaces, and shows that learning an invertible transformation to a simple Gaussian prior space substantially stabilizes downstream optimization and control\[[8](https://arxiv.org/html/2605.17137#bib.bib12),[5](https://arxiv.org/html/2605.17137#bib.bib13)\]\.

#### 3\.3\.1Normalizing flows\.

A*normalizing flow*is an invertible mappingFφ:ℝd→ℝdF\_\{\\varphi\}:\\mathbb\{R\}^\{d\}\\to\\mathbb\{R\}^\{d\}that transforms samples between a complex data distribution and a simple base distribution \(typically a standard Gaussian\), while providing an exact likelihood via the change\-of\-variables rule\. Concretely, flows are commonly implemented as a composition of invertible layersFφ=fL∘⋯∘f1F\_\{\\varphi\}=f\_\{L\}\\circ\\cdots\\circ f\_\{1\}with efficient inversion and Jacobian log\-determinant computation\[[17](https://arxiv.org/html/2605.17137#bib.bib18),[3](https://arxiv.org/html/2605.17137#bib.bib19),[8](https://arxiv.org/html/2605.17137#bib.bib12)\]\. Givenu=Fφ​\(z\)u=F\_\{\\varphi\}\(z\)and a base densityp𝒰​\(u\)p\_\{\\mathcal\{U\}\}\(u\), the induced density onzzis:

log⁡p𝒵​\(z\)=log⁡p𝒰​\(u\)\+log⁡\|det∂Fφ​\(z\)∂z\|,u=Fφ​\(z\)\.\\log p\_\{\\mathcal\{Z\}\}\(z\)\\;=\\;\\log p\_\{\\mathcal\{U\}\}\(u\)\\;\+\\;\\log\\left\|\\det\\frac\{\\partial F\_\{\\varphi\}\(z\)\}\{\\partial z\}\\right\|,\\qquad u=F\_\{\\varphi\}\(z\)\.\(10\)Invertibility is crucial in our setting because it allows us to optimize in the structured prior space and then map the optimized point back throughFφ−1F\_\{\\varphi\}^\{\-1\}for decoding\.

#### 3\.3\.2Density Modeling

We introduce a flow networkFφ:ℝd→ℝdF\_\{\\varphi\}:\\mathbb\{R\}^\{d\}\\to\\mathbb\{R\}^\{d\}that maps a code embeddingzzto a latent variableuuin a structured search space:

u=Fφ​\(z\),where​u∼𝒩​\(0,I\)\.u\\;=\\;F\_\{\\varphi\}\(z\),\\quad\\text\{where \}u\\sim\\mathcal\{N\}\(0,I\)\.\(11\)By enforcing a standard Gaussian prior onuu, the flow transform the irregular manifold of valid programs into a smooth, convex space\. Optimization steps taken in𝒰\\mathcal\{U\}correspond to traversing the data manifold in𝒵\\mathcal\{Z\}, preventing the search from drifting into low\-density regions where the decoder behavior is undefined\.

The flow is trained via maximum likelihood estimation on the dataset𝒟\\mathcal\{D\}\. We minimize the negative log\-likelihood of the observed program embeddings under the change\-of\-variables formula:

ℒflow​\(φ\)=−1N​∑i=1N\[log⁡p𝒩​\(Fφ​\(zi\)\)\+log⁡\|det∂Fφ​\(zi\)∂zi\|\],\\mathcal\{L\}\_\{\\mathrm\{flow\}\}\(\\varphi\)\\;=\\;\-\\frac\{1\}\{N\}\\sum\_\{i=1\}^\{N\}\\Big\[\\log p\_\{\\mathcal\{N\}\}\(F\_\{\\varphi\}\(z\_\{i\}\)\)\+\\log\\left\|\\det\\frac\{\\partial F\_\{\\varphi\}\(z\_\{i\}\)\}\{\\partial z\_\{i\}\}\\right\|\\Big\],\(12\)wherep𝒩p\_\{\\mathcal\{N\}\}is the standard normal density\. We implementFφF\_\{\\varphi\}using a sequence of invertible coupling layers \(e\.g\., RealNVP\[[4](https://arxiv.org/html/2605.17137#bib.bib14)\]\), ensuring that the Jacobian determinant is computationally efficient to compute\.

#### 3\.3\.3Search in the Prior Space

With the flow trained, we refine the search procedure described in Section[3\.1](https://arxiv.org/html/2605.17137#S3.SS1)\. Instead of optimizingzzdirectly, we perform gradient ascent in the prior space𝒰\\mathcal\{U\}\. Given a seedz\(0\)z^\{\(0\)\}, we map it tou\(0\)=Fφ​\(z\(0\)\)u^\{\(0\)\}=F\_\{\\varphi\}\(z^\{\(0\)\}\)\. The update rule becomes:

u\(t\+1\)=u\(t\)\+η​∇ufθ​\(u\(t\)\)\.u^\{\(t\+1\)\}\\;=\\;u^\{\(t\)\}\+\\eta\\,\\nabla\_\{u\}f\_\{\\theta\}\(u^\{\(t\)\}\)\.\(13\)This formulation naturally constrains the generated codes to remain "program\-like," as the Gaussian prior ensures we stay within the support of the training data\. After optimization, the final candidate is obtained by inverting the flow:z∗=Fφ−1​\(u\(T\)\)z^\{\*\}=F\_\{\\varphi\}^\{\-1\}\(u^\{\(T\)\}\)\.

### 3\.4Surrogate Training

To guide the gradient search, we require a differentiable surrogatefθf\_\{\\theta\}that approximates the black\-box performance scores​\(p\)s\(p\)\. We define this surrogate on the prior space𝒰\\mathcal\{U\}, such thatfθ:𝒰→ℝf\_\{\\theta\}:\\mathcal\{U\}\\to\\mathbb\{R\}\.

A straightforward approach would be to trainfθf\_\{\\theta\}to regress the absolute scores​\(p\)s\(p\)via Mean Squared Error \(MSE\)\. However, we reject this approach for three critical reasons\. Pointwise regression with an MSE loss attempts to match absolute scores, which can be brittle when evaluation is noisy and when the scale ofs​\(p\)s\(p\)varies across instance sets\. In contrast, the search procedure only requires a consistent ordering of candidates\. Pairwise ranking losses focus precisely on this relative information, yielding a stable learning signal even when absolute score calibration is unreliable\. A further advantage of ranking is that it provides a natural way to expand the supervised training signal from a small set of evaluated programs\. FromNNscored programs\{\(ui,si\)\}i=1N\\\{\(u\_\{i\},s\_\{i\}\)\\\}\_\{i=1\}^\{N\}, we construct preference pairs

𝒟r=\{\(ui,uj\)\|si\>sj\}\.\\mathcal\{D\}\_\{r\}=\\Big\\\{\(u\_\{i\},u\_\{j\}\)\\;\\big\|\\;s\_\{i\}\>s\_\{j\}\\Big\\\}\.This yields up toO​\(N2\)O\(N^\{2\}\)training comparisons derived from onlyNNevaluations, substantially increasing the number of supervised constraints available to fitfθf\_\{\\theta\}\. We emphasize that these pairwise examples are not independent; accordingly, we evaluate generalization using program\-level splits\.

We trainfθf\_\{\\theta\}to minimize the RankNet loss\[[1](https://arxiv.org/html/2605.17137#bib.bib8)\]\. The parametersω\\omegaare optimized to maximize the likelihood of correctly predicting the superior candidate:

ℒrank​\(ω\)=−𝔼\(ui,uj\)∼𝒟r​\[log⁡σ​\(fθ​\(ui\)−fθ​\(uj\)\)\],\\mathcal\{L\}\_\{\\mathrm\{rank\}\}\(\\omega\)=\-\\mathbb\{E\}\_\{\(u\_\{i\},u\_\{j\}\)\\sim\\mathcal\{D\}\_\{r\}\}\\Big\[\\log\\sigma\\Big\(f\_\{\\theta\}\(u\_\{i\}\)\-f\_\{\\theta\}\(u\_\{j\}\)\\Big\)\\Big\],\(14\)whereσ\\sigmais the sigmoid function\.

Algorithm 1Latent Heuristic Search \(LHS\)1:Pretrained encoder

EE, flow

FφF\_\{\\varphi\}, surrogate

fθf\_\{\\theta\}, mapper

ϕω\\phi\_\{\\omega\}, and frozen code LLM

𝒢\\mathcal\{G\}
2:Initial seed corpus

𝒫0\\mathcal\{P\}\_\{0\}; task prompt

II; benchmark evaluator

s​\(⋅\)s\(\\cdot\)
3:Budgets: outer iterations

RR, candidates per round

BB, ascent steps

TT, step size

η\\eta
4:

𝒟←∅\\mathcal\{D\}\\leftarrow\\emptyset
5:for all

p∈𝒫0p\\in\\mathcal\{P\}\_\{0\}do

6:Execute

ppto obtain score

s​\(p\)s\(p\)
7:

𝒟←𝒟∪\{\(p,s​\(p\)\)\}\\mathcal\{D\}\\leftarrow\\mathcal\{D\}\\cup\\\{\(p,s\(p\)\)\\\}
8:endfor

9:for

r=1r=1to

RRdo

10:Select seed programs

\{p\(b\)\}b=1B\\\{p^\{\(b\)\}\\\}\_\{b=1\}^\{B\}from

𝒟\\mathcal\{D\}\(e\.g\., top\-

kk\)

11:for

b=1b=1to

BBdo

12:

z\(0\)←E​\(p\(b\)\)z^\{\(0\)\}\\leftarrow E\(p^\{\(b\)\}\);

u\(0\)←Fφ​\(z\(0\)\)u^\{\(0\)\}\\leftarrow F\_\{\\varphi\}\(z^\{\(0\)\}\)
13:for

t=0t=0to

T−1T\-1do

14:

u\(t\+1\)←u\(t\)\+η​∇ufθ​\(u\(t\)\)u^\{\(t\+1\)\}\\leftarrow u^\{\(t\)\}\+\\eta\\,\\nabla\_\{u\}f\_\{\\theta\}\\big\(u^\{\(t\)\}\\big\)
15:endfor

16:

z∗←Fφ−1​\(u\(T\)\)z^\{\*\}\\leftarrow F\_\{\\varphi\}^\{\-1\}\(u^\{\(T\)\}\);

h←ϕω​\(z∗\)h\\leftarrow\\phi\_\{\\omega\}\(z^\{\*\}\)
17:Sample candidate code

p^∼𝒢​\(Concat​\(I,h\)\)\\hat\{p\}\\sim\\mathcal\{G\}\(\\mathrm\{Concat\}\(I,h\)\)
18:if

p^\\hat\{p\}is validthen

19:Execute

p^\\hat\{p\}to obtain

s​\(p^\)s\(\\hat\{p\}\)
20:

𝒟←𝒟∪\{\(p^,s​\(p^\)\)\}\\mathcal\{D\}\\leftarrow\\mathcal\{D\}\\cup\\\{\(\\hat\{p\},s\(\\hat\{p\}\)\)\\\}
21:endif

22:endfor

23:endfor

24:returnbest program

p∗=arg⁡max\(p,s\)∈𝒟⁡sp^\{\*\}=\\arg\\max\_\{\(p,s\)\\in\\mathcal\{D\}\}s

## 4Experiments

We empirically validate the proposed latent\-space heuristic discovery framework across a diverse suite of combinatorial optimization tasks\.

### 4\.1Data Generation

To train the LHS components, we constructed a dataset of heuristics across ten combinatorial optimization tasks\. We first generated approximately 100 seed heuristics per task using GPT\-5\.2\. We then expanded this set using a lightweight model \(GPT\-OSS\-20B\) via three augmentation strategies: syntactic rewriting, parametric tuning, and behavioral variation\. After validating the generated heuristics for correctness and uniqueness, we obtained a final corpus of approximately 4,500 heuristics\. This pooled dataset was used to train the global components \(normalizing flow and mapper\), while task\-specific subsets were used to train the individual performance surrogates\.

### 4\.2Baselines

We compare against representative LLM\-based heuristic search methods, including*FunSearch*\[[19](https://arxiv.org/html/2605.17137#bib.bib1)\],*Evolution of Heuristics \(EoH\)*\[[12](https://arxiv.org/html/2605.17137#bib.bib2)\],*ReEvo*\[[22](https://arxiv.org/html/2605.17137#bib.bib3)\], and*MCTS\-AHD*\[[23](https://arxiv.org/html/2605.17137#bib.bib10)\]\. We implement and evaluate all baselines using theLLM4ADplatform\[[13](https://arxiv.org/html/2605.17137#bib.bib11)\]\. Our experimental protocol is designed to ensure a fair comparison by isolating the effect of the search strategy from other sources of variation\. Across all tasks, every method \(including ours\) uses the same code\-generation backbone,Qwen3\-4B\-Instruct\-2507, and is allocated an identical LLM budget of 100 inference calls\. In addition, all methods are initialized from the same single reference program, ensuring a common starting point and comparable prior program distribution\. Since several baseline methods rely on an LLM\-querying procedure to generate an initial population, we adopt the same initialization strategy for our method to ensure methodological parity\. For all methods, the population size is set to 10\. To account for stochasticity in LLM sampling and search, each task\-method experiment is repeated with five independent random seeds; objective values are reported as the mean±\\pmstandard deviation across these five seeds, and runtimes are reported as averages over the same runs\. Additional implementation details are provided in the appendix\. For our approach, programs are embedded using aQwen3\-Embedding\-0\.6Bembedding model, which serves as the encoderEE\. All candidate heuristics are evaluated using a shared execution\-based protocol in the same sandbox\[[13](https://arxiv.org/html/2605.17137#bib.bib11)\], with identical instance distributions, time limits, and feasibility checks\.

### 4\.3Tasks

We formulate each task as a sequential constructive process executed by a programp∈𝒫p\\in\\mathcal\{P\}\. Given an instancexxand a partial solution stateσ\\sigma, the program selects an action \(e\.g\., the next city in a tour or the bin assignment for an incoming item\), which deterministically updates the state toσ′\\sigma^\{\\prime\}\. This rollout continues until a complete feasible solution is produced or a preset resource limit is reached\. The resulting solution quality is measured by the same execution\-based evaluator used throughout the paper, yielding an instance\-level objectiveCost​\(p;x\)\\mathrm\{Cost\}\(p;x\)\. Overall heuristic performance is then summarized by aggregatingCost​\(p;x\)\\mathrm\{Cost\}\(p;x\)over a fixed validation set, consistent with the definition ofy​\(p\)y\(p\)in Section[2](https://arxiv.org/html/2605.17137#S2)\. The five\-seed protocol described above is applied uniformly to the final evaluations for TSP, CVRP, OBP, and Knapsack; the task\-specific settings below therefore specify only the instance distributions and time limits\. For the routing tasks, search is performed on size\-50 validation instances, while final evaluation also includes larger held\-out instance sizes \(TSP\-100, TSP\-200, and CVRP\-100\), allowing us to assess cross\-size generalization beyond the search setting\.

##### Traveling Salesperson Problem\.

The Traveling Salesperson Problem \(TSP\) requires constructing a single continuous route that visits every city in a given set exactly once before returning to the starting point\. The constructive heuristic incrementally selects the next city to visit, with the sole objective of minimizing the total distance of the completed tour\.Settings:During the search phase, we utilize a validation set of 16 randomly generated instances of size 50 with a 30\-second timeout per instance\. For final evaluation \(Table[1](https://arxiv.org/html/2605.17137#S4.T1)\), we report performance on 100 instances of size 50, 100 instances of size 100, and 100 instances of size 200, all with a timeout of 60 seconds\.

Table 1:TSP Performance\.Comparison on TSP\-50, TSP\-100, and TSP\-200\. Objective values are mean±\\pmstandard deviation over seeds; time values are seed averages\. \(↓\\downarrow\) indicates lower is better\.
##### Capacitated Vehicle Routing Problem\.

The Capacitated Vehicle Routing Problem \(CVRP\) extends the routing challenge to a fleet of vehicles based at a central depot\. The goal is to design a set of routes that serve all customers without exceeding the maximum load capacity of any single vehicle\. The heuristic must ensure all demands are met while minimizing the total accumulated travel cost across the entire fleet\.Settings:For the search phase, we use 16 instances of size 50 with a capacity of 40 and a 30\-second timeout\. For final evaluation \(Table[2](https://arxiv.org/html/2605.17137#S4.T2)\), we use 100 instances of size 50 and 100 instances of size 100, keeping the capacity fixed at 40 and increasing the timeout to 60 seconds\.

Table 2:CVRP Performance\.Comparison on CVRP\-50 and CVRP\-100\. Objective values are mean±\\pmstandard deviation over five seeds; time values are seed averages\. \(↓\\downarrow\) indicates lower is better\.
##### Online Bin Packing\.

The Online Bin Packing Problem involves the sequential arrival of items, each requiring immediate and irrevocable assignment to a bin of fixed capacity\. The constructive heuristic must determine whether to place each incoming item into an existing bin or open a new one, with the primary objective of minimizing the total number of bins used to accommodate the stream of items\.Settings:Item sizes are drawn from a Weibull distribution\. During search, we validate on 5 instances of length 5k with a capacity of 100 \(timeout 30s\)\. For evaluation \(Table[3](https://arxiv.org/html/2605.17137#S4.T3)\), we report on 10 instances of length 1k and 10 instances of length 5k \(capacity 100\), with a 30\-second timeout\.

Table 3:Bin Packing Performance\.Comparison on OBP\-1k and OBP\-5k\. Objective values are mean±\\pmstandard deviation over five seeds; time values are seed averages\. \(↓\\downarrow\) indicates lower is better\.
##### Knapsack\.

The Knapsack Problem requires selecting a subset of items, each characterized by a specific weight and value, to maximize the total accumulated value without exceeding a strictly defined weight capacity\. The heuristic operates constructively by evaluating and selecting items for inclusion, ensuring that the aggregate weight remains within the feasible limit\.Settings:We follow the hard instance generation protocols of Pisinger\[[18](https://arxiv.org/html/2605.17137#bib.bib17)\]\. During search, we use a mixed validation set of 32 instances \(uncorrelated, weakly correlated, and strongly correlated\) with capacity 100 and 50 items, utilizing a 20\-second timeout\. For evaluation \(Table[4](https://arxiv.org/html/2605.17137#S4.T4)\), we test on 100 instances of each family \(capacity 100, 50 items\) with a 20\-second timeout\.

Table 4:Knapsack Performance\.Comparison by correlation type\. Objective values are mean±\\pmstandard deviation over five seeds; time values are seed averages\. \(↑\\uparrow\) indicates higher is better\.
##### Results\.

Tables[1](https://arxiv.org/html/2605.17137#S4.T1),[2](https://arxiv.org/html/2605.17137#S4.T2),[4](https://arxiv.org/html/2605.17137#S4.T4), and[3](https://arxiv.org/html/2605.17137#S4.T3)summarize the quantitative comparisons across four representative domains\. On the routing tasks \(Tables[1](https://arxiv.org/html/2605.17137#S4.T1)and[2](https://arxiv.org/html/2605.17137#S4.T2)\), LHS achieves the best average objective on both TSP and CVRP, consistently outperforming token\-space evolutionary baselines\[[19](https://arxiv.org/html/2605.17137#bib.bib1),[12](https://arxiv.org/html/2605.17137#bib.bib2),[22](https://arxiv.org/html/2605.17137#bib.bib3)\]and MCTS\-AHD\[[23](https://arxiv.org/html/2605.17137#bib.bib10)\]without increasing the execution time of the discovered heuristics\. For Knapsack \(Table[4](https://arxiv.org/html/2605.17137#S4.T4)\), all methods are tightly clustered; ReEvo attains the best objective, while LHS remains competitive\. For Online Bin Packing \(Table[3](https://arxiv.org/html/2605.17137#S4.T3)\), MCTS\-AHD achieves the best objective but with substantially slower runtime, whereas LHS matches the evolutionary baselines’ objective with moderate heuristic runtime\. Overall, these results demonstrate that continuous optimization over learned program representations yields state\-of\-the\-art heuristic quality on core routing domains while maintaining practical efficiency\.

### 4\.4Ablation Study

To investigate the contributions of the normalizing flow and the gradient\-based optimization strategy, we conduct an ablation study comparing our proposed framework \(LHS\) against two variants:

- •No\-Flow \(Optimization in𝒵\\mathcal\{Z\}\):We remove the normalizing flow and perform gradient ascent directly on the raw encoder embeddingszz\. This tests the hypothesis that the flow provides a necessary regularization of the search space\.
- •No\-Grad \(Interpolation\):We remove the differentiable surrogate model to test if simple geometric exploration is sufficient\. Instead of gradient ascent, we generate candidates by computing the linear midpoint between pairs of top\-performing programs:unew←0\.5​ua\+0\.5​ubu\_\{\\text\{new\}\}\\leftarrow 0\.5u\_\{a\}\+0\.5u\_\{b\}\. This represents a conservative baseline that explores the center of the convex hull defined by the parent programs\.

We evaluate these variants on the TSP task using the same computational budget\. We report theSuccess Rate\(percentage of generated heuristics that compile and execute within time limits\) and theBenchmark Score\(average objective cost, lower is better\)\.

Table 5:Ablation Analysis\.Comparison of search strategies on TSP\.LHS\(Ours\) achieves the best \(lowest\) objective cost\.No\-Flowsuffers from low validity due to latent space anisotropy\.No\-Gradshows high stability but poor optimization capability\.##### Impact of the Normalizing Flow\.

Comparing LHS to theNo\-Flowvariant reveals the critical role of the prior space𝒰\\mathcal\{U\}in stabilizing optimization\. Searching directly in the raw encoder space𝒵\\mathcal\{Z\}results in a significantly lower success rate \(61%61\\%vs74%74\\%\)\. Raw encoder embeddings often suffer from anisotropy and "holes" in the manifold; gradient steps in𝒵\\mathcal\{Z\}are prone to traversing low\-density regions where the decoder behavior is undefined, leading to invalid syntax or runtime errors\. The normalizing flow maps these embeddings to a smooth, convex Gaussian prior, ensuring that optimization trajectories remain closer to the manifold of valid programs\.

##### Gradient Search vs\. Interpolation\.

The comparison with theNo\-Gradbaseline highlights a fundamental trade\-off between optimization power and stability\. Interpolation achieves the highest success rate \(86%86\\%\) because interpolating between valid points generally keeps the candidate within the data hull\. However, it yields the worst performance \(6\.79\), as it lacks the directional signal required to discover novel regions of the landscape\. Conversely, LHS leverages the surrogate gradient to actively minimize cost, achieving the best score \(6\.61\), but with a slightly lower success rate than interpolation \(74%74\\%\)\. This indicates that while the gradient is essential for quality, it introduces a risk of pushing the latent vector "off\-manifold" if the optimization trajectory is too aggressive\. Consequently, the learning rateη\\etaand the number of optimization stepsTTare critical hyperparameters that must be tuned carefully to balance the aggressiveness of the search with the validity of the generated programs\.

### 4\.5Visualization of Heuristics

To analyze the topology of the program space, we project the TSP heuristic embeddings \(learned by the encoderEE\) into two dimensions using t\-SNE\. Figure[1](https://arxiv.org/html/2605.17137#S4.F1)The t\-SNE plot displays distinct clusters of heuristics\. Notably, areas with concentrated yellow points represent heuristics producing shorter routes \(better solutions\), while darker clusters indicate less effective ones\. This pattern suggests that similar heuristics tend to group together, both in terms of their characteristics and their relative performance on the TSP This behavior is consistent with the hypothesis thatEEclusters semantically similar heuristics i\.e\., programs implementing related construction strategies into nearby latent neighborhoods\. Consequently, continuous latent\-space optimization is able to exploit this structure by moving candidates toward high\-quality neighborhoods, increasing the likelihood that decoding yields strong executable heuristics\.

![Refer to caption](https://arxiv.org/html/2605.17137v1/tsp_construct_tsne.png)Figure 1:Latent Space Visualization\.A 2D t\-SNE projection of program embeddings for the TSP\. The color gradient represents the heuristic score, showing that high\-performing programs \(lighter colors\) cluster together rather than being randomly dispersed\.

## 5Discussion

This work introduces a continuous alternative to token\-level evolutionary search by shifting heuristic discovery into a learned latent manifold\. Our results suggest that continuous optimization is a powerful complement to discrete search; it enables systematic navigation via gradient information from a differentiable surrogate while maintaining the rigor of execution\-based evaluation\.

Empirical results demonstrate that LHS is competitive with state\-of\-the\-art baselines, often achieving superior objectives within a fixed evaluation budget\. A key finding is the necessity of optimizing within the flow prior space𝒰\\mathcal\{U\}\. While searching in𝒰\\mathcal\{U\}significantly improves stability and validity compared to the raw embedding space𝒵\\mathcal\{Z\}, it also highlights a critical trade\-off: overly aggressive ascent steps can push candidates into regions that, while favored by the surrogate, remain brittle during decoding\. This suggests that latent optimization benefits from a "trust\-region" style approach to ensure candidates stay on the manifold of executable code\.

The current implementation of LHS relies on a relatively small initial corpus to train the flow, mapper, and surrogate\. While we expect that moving to a larger, more versatile dataset would significantly improve the stabilization of the latent manifold and the robustness of the search, scaling the data comes with a practical trade\-off in terms of upfront generation and evaluation costs\. This initial investment is most justifiable in scenarios involving repeated searches across related tasks or long\-term algorithmic refinement\. However, because the normalizing flow and mapper learn general properties of program syntax rather than task\-specific logic, these components can be amortized across different optimization domains\. In this multi\-task regime, only the lightweight surrogate needs to be retrained for new problems, making the overall framework increasingly efficient as it is applied to broader classes of heuristics\. This efficiency is critical because the efficacy of LHS is fundamentally tied to the surrogate modelfθf\_\{\\theta\}\. As optimization moves away from the initial training distribution, the surrogate may suffer from distribution shift, leading to “false ascents” where predicted improvements do not manifest in execution\. While our ranking\-based supervision provides a good signal, we expect that scaling the initial program corpus and incorporating active learning where the surrogate is iteratively retrained on newly discovered heuristics will further stabilize the search\.

In conclusion, by treating programs as continuous latent vectors, we move toward a more differentiable pipeline for automated algorithm design\. Future work will explore hybridizing LHS with discrete mutations to combine the global exploration of evolution with the local refinement of gradient\-based search, aditionaly othrs evolution search in latent space ideas can be explored

In conclusion, by treating programs as continuous latent vectors, we move toward a more differentiable pipeline for automated algorithm design\. Future work will explore hybridizing LHS with discrete evolutionary operators \(e\.g\., mutation and recombination\) to combine the global exploration of evolutionary search with the local refinement of gradient\-based optimization\[[19](https://arxiv.org/html/2605.17137#bib.bib1),[12](https://arxiv.org/html/2605.17137#bib.bib2),[22](https://arxiv.org/html/2605.17137#bib.bib3)\]\. In addition, broader families of latent\-space evolutionary strategies and population\-based search methods could be incorporated\[[6](https://arxiv.org/html/2605.17137#bib.bib20)\]\.

## References

- \[1\]C\. Burges, T\. Shaked, E\. Renshaw, A\. Lazier, M\. Deeds, N\. Hamilton, and G\. Hullender\(2005\)Learning to rank using gradient descent\.InProceedings of the 22nd International Conference on Machine Learning,ICML ’05,New York, NY, USA,pp\. 89–96\.External Links:ISBN 1595931805,[Link](https://doi.org/10.1145/1102351.1102363),[Document](https://dx.doi.org/10.1145/1102351.1102363)Cited by:[§3\.4](https://arxiv.org/html/2605.17137#S3.SS4.p3.2)\.
- \[2\]P\. V\. T\. Dat, L\. Doan, and H\. T\. T\. Binh\(2025\)Hsevo: elevating automatic heuristic design with diversity\-driven harmony search and genetic algorithm using llms\.InProceedings of the AAAI Conference on Artificial Intelligence,Vol\.39,pp\. 26931–26938\.Note:[https://github\.com/datphamvn/HSEvo](https://github.com/datphamvn/HSEvo)Cited by:[§1](https://arxiv.org/html/2605.17137#S1.p2.1)\.
- \[3\]L\. Dinh, D\. Krueger, and Y\. Bengio\(2015\)NICE: non\-linear independent components estimation\.External Links:1410\.8516,[Link](https://arxiv.org/abs/1410.8516)Cited by:[§3\.3\.1](https://arxiv.org/html/2605.17137#S3.SS3.SSS1.p1.5)\.
- \[4\]L\. Dinh, J\. Sohl\-Dickstein, and S\. Bengio\(2017\)Density estimation using real NVP\.InInternational Conference on Learning Representations,External Links:[Link](https://openreview.net/forum?id=HkpbnH9lx)Cited by:[§3\.3\.2](https://arxiv.org/html/2605.17137#S3.SS3.SSS2.p2.3)\.
- \[5\]Y\. Gu, X\. Feng, S\. Ma, L\. Zhang, H\. Gong, W\. Zhong, and B\. Qin\(2023\-07\)Controllable text generation via probability density estimation in the latent space\.InProceedings of the 61st Annual Meeting of the Association for Computational Linguistics \(Volume 1: Long Papers\),A\. Rogers, J\. Boyd\-Graber, and N\. Okazaki \(Eds\.\),Toronto, Canada,pp\. 12590–12616\.External Links:[Link](https://aclanthology.org/2023.acl-long.704/),[Document](https://dx.doi.org/10.18653/v1/2023.acl-long.704)Cited by:[§3\.3](https://arxiv.org/html/2605.17137#S3.SS3.p1.3)\.
- \[6\]N\. Hansen\(2023\)The cma evolution strategy: a tutorial\.External Links:1604\.00772,[Link](https://arxiv.org/abs/1604.00772)Cited by:[§5](https://arxiv.org/html/2605.17137#S5.p5.1)\.
- \[7\]C\. C\. Lee and D\. T\. Lee\(1985\-07\)A simple on\-line bin\-packing algorithm\.J\. ACM32\(3\),pp\. 562–572\.External Links:ISSN 0004\-5411,[Link](https://doi.org/10.1145/3828.3833),[Document](https://dx.doi.org/10.1145/3828.3833)Cited by:[§1](https://arxiv.org/html/2605.17137#S1.p1.1)\.
- \[8\]B\. Li, H\. Zhou, J\. He, M\. Wang, Y\. Yang, and L\. Li\(2020\-11\)On the sentence embeddings from pre\-trained language models\.InProceedings of the 2020 Conference on Empirical Methods in Natural Language Processing \(EMNLP\),B\. Webber, T\. Cohn, Y\. He, and Y\. Liu \(Eds\.\),Online,pp\. 9119–9130\.External Links:[Link](https://aclanthology.org/2020.emnlp-main.733/),[Document](https://dx.doi.org/10.18653/v1/2020.emnlp-main.733)Cited by:[§3\.3\.1](https://arxiv.org/html/2605.17137#S3.SS3.SSS1.p1.5),[§3\.3](https://arxiv.org/html/2605.17137#S3.SS3.p1.3)\.
- \[9\]X\. L\. Li and P\. Liang\(2021\)Prefix\-tuning: optimizing continuous prompts for generation\.External Links:2101\.00190Cited by:[§3\.2](https://arxiv.org/html/2605.17137#S3.SS2.p1.3)\.
- \[10\]S\. Lin and B\. W\. Kernighan\(1973\-04\)An effective heuristic algorithm for the traveling\-salesman problem\.Oper\. Res\.21\(2\),pp\. 498–516\.External Links:ISSN 0030\-364X,[Link](https://doi.org/10.1287/opre.21.2.498),[Document](https://dx.doi.org/10.1287/opre.21.2.498)Cited by:[§1](https://arxiv.org/html/2605.17137#S1.p1.1)\.
- \[11\]P\. Liskowski, K\. Krawiec, N\. E\. Toklu, and J\. Swan\(2020\)Program synthesis as latent continuous optimization: evolutionary search in neural embeddings\.Proceedings of the 2020 Genetic and Evolutionary Computation Conference\.External Links:[Link](https://api.semanticscholar.org/CorpusID:220252273)Cited by:[§1](https://arxiv.org/html/2605.17137#S1.p3.1)\.
- \[12\]F\. Liu, X\. Tong, M\. Yuan, X\. Lin, F\. Luo, Z\. Wang, Z\. Lu, and Q\. Zhang\(2024\)Evolution of heuristics: towards efficient automatic algorithm design using large language model\.InInternational Conference on Machine Learning \(ICML\),External Links:[Link](https://arxiv.org/abs/2401.02051)Cited by:[§1](https://arxiv.org/html/2605.17137#S1.p2.1),[§4\.2](https://arxiv.org/html/2605.17137#S4.SS2.p1.2),[§4\.3](https://arxiv.org/html/2605.17137#S4.SS3.SSS0.Px5.p1.1),[§5](https://arxiv.org/html/2605.17137#S5.p5.1)\.
- \[13\]F\. Liu, R\. Zhang, Z\. Xie, R\. Sun, K\. Li, X\. Lin, Z\. Wang, Z\. Lu, and Q\. Zhang\(2024\)LLM4AD: a platform for algorithm design with large language model\.External Links:2412\.17287,[Link](https://arxiv.org/abs/2412.17287)Cited by:[§4\.2](https://arxiv.org/html/2605.17137#S4.SS2.p1.2)\.
- \[14\]M\. Macfarlane and C\. Bonnet\(2025\)Searching latent program spaces\.InThe Thirty\-ninth Annual Conference on Neural Information Processing Systems,External Links:[Link](https://openreview.net/forum?id=CsXKGIqZtr)Cited by:[§1](https://arxiv.org/html/2605.17137#S1.p3.1)\.
- \[15\]S\. Mandal, T\. A\. Anderson, J\. Turek, J\. Gottschlich, and A\. Muzahid\(2023\)Synthesizing programs with continuous optimization\.External Links:2211\.00828,[Link](https://arxiv.org/abs/2211.00828)Cited by:[§1](https://arxiv.org/html/2605.17137#S1.p3.1)\.
- \[16\]A\. Novikov, N\. Vũ, M\. Eisenberger, E\. Dupont, P\. Huang, A\. Z\. Wagner, S\. Shirobokov, B\. Kozlovskii, F\. J\. R\. Ruiz, A\. Mehrabian, M\. P\. Kumar, A\. See, S\. Chaudhuri, G\. Holland, A\. Davies, S\. Nowozin, P\. Kohli, and M\. Balog\(2025\)AlphaEvolve: a coding agent for scientific and algorithmic discovery\.External Links:2506\.13131,[Link](https://arxiv.org/abs/2506.13131)Cited by:[§1](https://arxiv.org/html/2605.17137#S1.p2.1)\.
- \[17\]G\. Papamakarios, E\. Nalisnick, D\. J\. Rezende, S\. Mohamed, and B\. Lakshminarayanan\(2021\-01\)Normalizing flows for probabilistic modeling and inference\.J\. Mach\. Learn\. Res\.22\(1\)\.External Links:ISSN 1532\-4435Cited by:[§3\.3\.1](https://arxiv.org/html/2605.17137#S3.SS3.SSS1.p1.5)\.
- \[18\]D\. Pisinger\(2005\)Where are the hard knapsack problems?\.Computers & Operations Research32\(9\),pp\. 2271–2284\.External Links:ISSN 0305\-0548,[Document](https://dx.doi.org/https%3A//doi.org/10.1016/j.cor.2004.03.002),[Link](https://www.sciencedirect.com/science/article/pii/S030505480400036X)Cited by:[§4\.3](https://arxiv.org/html/2605.17137#S4.SS3.SSS0.Px4.p1.1)\.
- \[19\]B\. Romera\-Paredes, M\. Barekatain, A\. Novikov, M\. Balog, M\. P\. Kumar, E\. Dupont, F\. J\. R\. Ruiz, J\. Ellenberg, P\. Wang, O\. Fawzi, P\. Kohli, and A\. Fawzi\(2023\)Mathematical discoveries from program search with large language models\.Nature\.External Links:[Document](https://dx.doi.org/10.1038/s41586-023-06924-6)Cited by:[§1](https://arxiv.org/html/2605.17137#S1.p2.1),[§4\.2](https://arxiv.org/html/2605.17137#S4.SS2.p1.2),[§4\.3](https://arxiv.org/html/2605.17137#S4.SS3.SSS0.Px5.p1.1),[§5](https://arxiv.org/html/2605.17137#S5.p5.1)\.
- \[20\]D\. Trivedi, J\. Zhang, S\. Sun, and J\. J\. Lim\(2021\)Learning to synthesize programs as interpretable and generalizable policies\.InAdvances in Neural Information Processing Systems,A\. Beygelzimer, Y\. Dauphin, P\. Liang, and J\. W\. Vaughan \(Eds\.\),External Links:[Link](https://openreview.net/forum?id=wP9twkexC3V)Cited by:[§1](https://arxiv.org/html/2605.17137#S1.p3.1)\.
- \[21\]T\. Vidal\(2022\)Hybrid genetic search for the cvrp: open\-source implementation and swap\* neighborhood\.Computers & Operations Research140,pp\. 105643\.External Links:ISSN 0305\-0548,[Document](https://dx.doi.org/https%3A//doi.org/10.1016/j.cor.2021.105643),[Link](https://www.sciencedirect.com/science/article/pii/S030505482100349X)Cited by:[§1](https://arxiv.org/html/2605.17137#S1.p1.1)\.
- \[22\]H\. Ye, J\. Wang, Z\. Cao, F\. Berto, C\. Hua, H\. Kim, J\. Park, and G\. Song\(2024\)ReEvo: large language models as hyper\-heuristics with reflective evolution\.InAdvances in Neural Information Processing Systems,Note:[https://github\.com/ai4co/reevo](https://github.com/ai4co/reevo)Cited by:[§1](https://arxiv.org/html/2605.17137#S1.p2.1),[§4\.2](https://arxiv.org/html/2605.17137#S4.SS2.p1.2),[§4\.3](https://arxiv.org/html/2605.17137#S4.SS3.SSS0.Px5.p1.1),[§5](https://arxiv.org/html/2605.17137#S5.p5.1)\.
- \[23\]Z\. Zheng, Z\. Xie, Z\. Wang, and B\. Hooi\(2025\)Monte carlo tree search for comprehensive exploration in LLM\-based automatic heuristic design\.InForty\-second International Conference on Machine Learning,External Links:[Link](https://openreview.net/forum?id=Do1OdZzYHr)Cited by:[§1](https://arxiv.org/html/2605.17137#S1.p2.1),[§4\.2](https://arxiv.org/html/2605.17137#S4.SS2.p1.2),[§4\.3](https://arxiv.org/html/2605.17137#S4.SS3.SSS0.Px5.p1.1)\.

## Appendix 0\.AImplementation Details

We provide a detailed breakdown of the model architectures, training objectives, and hyperparameters used in the LHS framework\. The pipeline consists of four primary components: the code encoder, the normalizing flow, the latent\-to\-prompt mapper, and the performance surrogate\.

### 0\.A\.1Latent Space Construction \(Encoder\)

To map discrete program syntax into a continuous manifold, we employ theQwen3\-Embedding\-0\.6Bmodel\. To ensure computational efficiency during gradient\-based search, we leverageMatryoshka Representation Learning\. Instead of utilizing the full high\-dimensional embedding 1024, we truncate the embedding vector to a compact dimensionality ofd=128d=128\. This 128\-dimensional vector serves as the inputzzfor all subsequent modules \(Flow, Mapper, and Surrogate\), significantly reducing the parameter count of downstream networks without sacrificing semantic fidelity\.

### 0\.A\.2Normalizing Flow Architecture

The normalizing flowFφ:ℝd→ℝdF\_\{\\varphi\}:\\mathbb\{R\}^\{d\}\\to\\mathbb\{R\}^\{d\}is implemented as a RealNVP \(Real\-valued Non\-Volume Preserving\) transformation\. Its primary role is to map the irregular latent distribution of the encoder to a standard Gaussian prior𝒩​\(0,I\)\\mathcal\{N\}\(0,I\)\.

- •Structure:The flow consists of44coupling layers\. We use an alternating binary mask strategy \(checkerboard pattern\) where parity alternates at each layer to ensure all dimensions are transformed\.
- •Coupling Networks:Each coupling layer utilizes two internal multi\-layer perceptrons \(MLPs\) to predict the scalessand translationttparameters\. These MLPs have a single hidden layer of size128128with LeakyReLU activation \(negative slope0\.20\.2\)\.
- •Stability:To prevent numerical instability during inversion, the output of the scale network is passed through atanh\\tanhactivation and multiplied by a learnable scale factor \(initialized to0\.80\.8\)\. Additionally,ActNorm\(Activation Normalization\) layers are inserted before each coupling layer to normalize feature statistics\.
- •Training:The flow is trained via maximum likelihood estimation \(minimizing NLL\) using the AdamW optimizer with a learning rate of1​e−31e^\{\-3\}and a weight decay of1​e−51e^\{\-5\}\.

### 0\.A\.3Latent\-to\-Prompt Mapper

To translate optimized latent vectors back into executable code, we introduce aLow\-Rank Mapperϕω\\phi\_\{\\omega\}\. Unlike standard MLP mappers which can be parameter\-heavy, our architecture is designed for parameter efficiency, preventing overfitting on the relatively small heuristic datasets\.

- •Objective:Map a latent vectorz∈ℝ128z\\in\\mathbb\{R\}^\{128\}to a sequence ofK=16K=16continuous soft\-prompt tokens compatible with the LLM decoder\.
- •Feature Expansion:The inputzzis first projected linearly to a tensor of shape\(K,128\)\(K,128\)\. Learnable positional embeddings are added to differentiate the tokens\.
- •Self\-Attention:To capture dependencies between prompt tokens, we apply a lightweight Multi\-Head Self\-Attention block \(2 heads, embedding dimension 128, dropout 0\.1\)\.
- •Shared Up\-Projection:To map the processed 128\-dimensional features to the high\-dimensional embedding space of the decoder \(e\.g\., 2560 dimensions for Qwen\-4B\), we use a shared MLP \(Pointwise Feed\-Forward Network\)\. This network projects each token individually, sharing weights across the sequence to minimize parameter count\.
- •Total Parameters:This low\-rank design results in approximately 1\.8M trainable parameters, compared to\>40\>40M for a standard full\-rank MLP mapper\.

### 0\.A\.4Performance Surrogate

We train a task\-specific differentiable surrogateRθ​\(u\)R\_\{\\theta\}\(u\)to predict heuristic quality from prior\-space vectorsuu\. The surrogate is a small MLP \(2 hidden layers, width 256\) with ReLU and dropout 0\.1, outputting a scalar score\. Supervision is provided via pairwise ranking: from evaluated programs we create preference pairs\(ui,uj\)\(u\_\{i\},u\_\{j\}\)wheneversi\>sjs\_\{i\}\>s\_\{j\}, and minimize a logistic ranking loss with temperatureτ\\tau:

ℒrank=−log⁡σ​\(τ​\(Rθ​\(ui\)−Rθ​\(uj\)\)\)\.\\mathcal\{L\}\_\{\\text\{rank\}\}=\-\\log\\sigma\(\\tau\(R\_\{\\theta\}\(u\_\{i\}\)\-R\_\{\\theta\}\(u\_\{j\}\)\)\)\.

### 0\.A\.5LLM for Code Generation

The decoding phase utilizes theQwen3\-4B\-InstructLarge Language Model\. The model is kept frozen during the training of the Mapper\. We employFlash Attention 2and gradient checkpointing to optimize memory usage during the soft\-prompt tuning phase\. During search we generate candidate code using a temperature of0\.70\.7and nucleus sampling \(t​o​p​\_​p=0\.9top\\\_p=0\.9\)\.

### 0\.A\.6Prompt augmentation

For mapper training, the prompt builder samples one of three prompt families per example: a task\-specific prompt, a broader problem\-class prompt, or a general coding prompt\. This improves robustness to prompt phrasing at generation time and reduces overfitting to a single instruction template\.

### 0\.A\.7Search Strategy and Hyperparameters

The search process operates iteratively in the structured prior space𝒰\\mathcal\{U\}\. In each iteration, we select high\-performing seed programs, optimize their representations via gradient ascent, and decode the results\.

- •Initialization \(Softmax Sampling\):To balance exploitation of known high\-quality heuristics with exploration of the manifold, we do not simply select the top\-kkprograms\. Instead, we compute selection probabilities based on the actual benchmark scores in our database\. For a pool of candidate programs, the probability of selecting programpip\_\{i\}as a seed is given by the softmax distribution: P​\(pi\)=exp⁡\(s​\(pi\)\)∑jexp⁡\(s​\(pj\)\)P\(p\_\{i\}\)=\\frac\{\\exp\(s\(p\_\{i\}\)\)\}\{\\sum\_\{j\}\\exp\(s\(p\_\{j\}\)\)\}\(15\)We sample 5 seeds without replacement from the top\-10 programs in the database for each search iteration\.
- •Gradient Ascent:We perform gradient ascent on the surrogatefθ​\(u\)f\_\{\\theta\}\(u\)forT=5T=5steps with a learning rate ofη=0\.001\\eta=0\.001\.
- •Decoding:The optimized prior vectorsu∗u^\{\*\}are mapped back to latent codesz∗=F−1​\(u∗\)z^\{\*\}=F^\{\-1\}\(u^\{\*\}\)and passed to the mapper\.

#### 0\.A\.7\.1Search\-time decoding prompts \(per task\)

During latent search, candidate latents are decoded by conditioning the frozen LLM with a task\-specific decoding prompt:

- •TSP \(constructive\): Based on the above context, write a Python function that implements a TSP heuristic\. The function should construct a tour by selecting the next city step\-by\-step\.
- •CVRP \(constructive\): Based on the above context, write a Python function that implements a CVRP heuristic\. The function should construct routes by selecting the next customer step\-by\-step\.
- •Knapsack \(constructive\): Based on the above context, write a Python function that implements a Knapsack heuristic\. The function should select items to include in the knapsack step\-by\-step\.
- •Online Bin Packing: Based on the above context, write a Python function that implements an online bin packing heuristic\. The function should assign each item to a bin as it arrives\.

## Appendix 0\.BData Generation and Augmentation

To train the normalizing flow and latent\-to\-prompt mapper, we require a diverse corpus of heuristic programs that covers a wide range of algorithmic behaviors and syntactic structures\. We constructed this dataset using a two\-stage process: initial seed generation followed by systematic data augmentation\.

### 0\.B\.1Initial Seed Generation

For each of the ten combinatorial optimization tasks, we generated an initial set of seed heuristics usingGPT\-5\.2\. The model was prompted with the task description and a template function signature to produce valid, executable Python code implementing constructive heuristics\. This process yielded approximately 100 unique seed programs per task\.

### 0\.B\.2Augmentation Strategy

To expand the seed corpus into a robust training dataset, we implemented a three\-tier augmentation pipeline designed to induce specific types of variation in the latent space\. We utilized the following strategies:

1. 1\.Syntactic & Structural Rewriting:Transformations that alter the code surface form without necessarily changing the algorithm’s logic\. This includes variable renaming \(e\.g\.,distance\_matrix→\\todist\_mat\), control flow modifications \(e\.g\., swappingforloops with list comprehensions\), and replacing mathematical operations with equivalent NumPy functions \(e\.g\.,x\*\*2→\\tonp\.square\(x\)\)\.
2. 2\.Parameter & Hyperparameter Tuning:Modifications to the numerical constants within the heuristic\. This involves perturbing weights in scoring functions \(e\.g\., changing0\.8⋅dist\+0\.2⋅cost0\.8\\cdot\\text\{dist\}\+0\.2\\cdot\\text\{cost\}to0\.6⋅dist\+0\.4⋅cost0\.6\\cdot\\text\{dist\}\+0\.4\\cdot\\text\{cost\}\) and adjusting threshold values for filtering operations \(e\.g\., top\-kkselection\)\.
3. 3\.Semantic\-Preserving Behavioral Diversity:Changes that introduce slight algorithmic variations while maintaining the core logic\. Techniques include adding deterministic noise for tie\-breaking, swapping aggregation methods \(e\.g\.,np\.meanvsnp\.median\), and replacing greedy selection \(argmin\) with softmin or random choice among top candidates\.

To ensure quality, all generated heuristics were validated for correct Python syntax and the presence of the required function signature\. We also filtered for structural uniqueness to remove near\-duplicate programs\.

### 0\.B\.3Final Dataset Composition

The final augmented dataset consists of 4,624 valid heuristic programs distributed across ten domains\. Table[6](https://arxiv.org/html/2605.17137#Pt0.A2.T6)summarizes the distribution of heuristics per task\.

Table 6:Distribution of Heuristics in the Augmented Dataset

## Appendix 0\.CBest Discovered Heuristics

We present the best heuristics discovered by each method for the Traveling Salesperson Problem \(TSP\)\.

### 0\.C\.1Traveling Salesperson Problem \(TSP\)

#### 0\.C\.1\.1Evolution of Heuristics \(EoH\)

importnumpyasnp

defselect\_next\_node\(current\_node:int,destination\_node:int,unvisited\_nodes:np\.ndarray,distance\_matrix:np\.ndarray\)\-\>int:

"""

Designanovelalgorithmtoselectthenextnodeineachstep\.

"""

iflen\(unvisited\_nodes\)==1:

returnunvisited\_nodes\[0\]

distances=distance\_matrix\[current\_node,unvisited\_nodes\]

centrality=np\.sum\(distance\_matrix\[:,unvisited\_nodes\],axis=0\)

centrality=centrality/np\.max\(centrality\)ifnp\.max\(centrality\)\>0elsenp\.ones\_like\(centrality\)

proximity\_score=1/\(distances\+1e\-6\)

direction\_vector=np\.array\(\[destination\_node\-current\_node\]\)

directional\_consistency=np\.abs\(unvisited\_nodes\-current\_node\)/\(np\.max\(unvisited\_nodes\)\-np\.min\(unvisited\_nodes\)\+1e\-6\)

directional\_consistency=directional\_consistency/np\.max\(directional\_consistency\)ifnp\.max\(directional\_consistency\)\>0elsenp\.ones\_like\(directional\_consistency\)

bias\_weight=0\.3

trajectory\_score=directional\_consistency\*bias\_weight

combined\_score=\(proximity\_score\*centrality\)\+\(trajectory\_score\*\(1\-centrality\)\)

entropy\_noise=np\.random\.dirichlet\(\[0\.5\]\*len\(unvisited\_nodes\)\)

combined\_score\+=entropy\_noise\*0\.05

next\_node\_idx=np\.argmax\(combined\_score\)

returnunvisited\_nodes\[next\_node\_idx\]

#### 0\.C\.1\.2FunSearch

importnumpyasnp

defselect\_next\_node\(current\_node:int,destination\_node:int,unvisited\_nodes:np\.ndarray,distance\_matrix:np\.ndarray\)\-\>int:

"""

Designanovelalgorithmtoselectthenextnodeineachstep\.

Args:

current\_node:IDofthecurrentnode\.

destination\_node:IDofthedestinationnode\.

unvisited\_nodes:ArrayofIDsofunvisitednodes\.

distance\_matrix:Distancematrixofnodes\.

Return:

IDofthenextnodetovisit\.

"""

iflen\(unvisited\_nodes\)==0:

return\-1

distances\_to\_unvisited=distance\_matrix\[current\_node,unvisited\_nodes\]

next\_node\_idx=np\.argmin\(distances\_to\_unvisited\)

next\_node=unvisited\_nodes\[next\_node\_idx\]

returnnext\_node

#### 0\.C\.1\.3MCTS\-AHD

importnumpyasnp

defselect\_next\_node\(current\_node:int,destination\_node:int,unvisited\_nodes:np\.ndarray,distance\_matrix:np\.ndarray\)\-\>int:

"""

Designanovelalgorithmtoselectthenextnodeineachstep\.

"""

iflen\(unvisited\_nodes\)==1:

returnunvisited\_nodes\[0\]

best\_next=unvisited\_nodes\[0\]

best\_score=float\(’inf’\)

current\_dists=\[distance\_matrix\[current\_node\]\[node\]fornodeinunvisited\_nodes\]

avg\_dist=np\.mean\(current\_dists\)

unvisited\_density=len\(unvisited\_nodes\)

base\_threshold=np\.percentile\(current\_dists,70\)

dynamic\_threshold=max\(base\_threshold,avg\_dist\*1\.1\)

neighbor\_counts=\[\]

fornodeinunvisited\_nodes:

shared\_neighbors=0

forotherinunvisited\_nodes:

ifother\!=node:

ifdistance\_matrix\[node\]\[other\]<dynamic\_threshold:

shared\_neighbors\+=1

neighbor\_counts\.append\(shared\_neighbors\)

edge\_frequency=np\.zeros\(len\(unvisited\_nodes\)\)

fori,nodeinenumerate\(unvisited\_nodes\):

freq=sum\(1forjinrange\(len\(unvisited\_nodes\)\)

ifj\!=ianddistance\_matrix\[node\]\[unvisited\_nodes\[j\]\]<dynamic\_threshold\)

edge\_frequency\[i\]=freq/len\(unvisited\_nodes\)iflen\(unvisited\_nodes\)\>0else0

dest\_factor=np\.array\(\[1/\(1\+distance\_matrix\[node\]\[destination\_node\]\*\*1\.5\)

fornodeinunvisited\_nodes\]\)

fori,nodeinenumerate\(unvisited\_nodes\):

dist\_cost=current\_dists\[i\]

conn\_gain=neighbor\_counts\[i\]/\(unvisited\_density\*1\.2\)

dest\_gain=dest\_factor\[i\]

temp\_penalty=edge\_frequency\[i\]\*0\.6

score=dist\_cost\+0\.35\*conn\_gain\+0\.3\*dest\_gain\-0\.25\*temp\_penalty

ifdist\_cost<avg\_dist\*0\.95:

ifscore<best\_score:

best\_score=score

best\_next=node

returnbest\_next

#### 0\.C\.1\.4ReEvo

importnumpyasnp

defselect\_next\_node\(current\_node:int,destination\_node:int,unvisited\_nodes:np\.ndarray,distance\_matrix:np\.ndarray\)\-\>int:

"""

Designanovelalgorithmtoselectthenextnodeineachstep\.

"""

iflen\(unvisited\_nodes\)==1:

returnunvisited\_nodes\[0\]

distances=distance\_matrix\[current\_node,unvisited\_nodes\]

node\_coords=np\.array\(\[

\(np\.sin\(i\*0\.3\+np\.random\.rand\(\)\*0\.1\)\*\(1\+0\.5\*np\.cos\(i\*0\.1\+np\.random\.rand\(\)\*0\.2\)\),

np\.cos\(i\*0\.3\+np\.random\.rand\(\)\*0\.1\)\*\(1\+0\.5\*np\.sin\(i\*0\.1\+np\.random\.rand\(\)\*0\.2\)\)\)

foriinrange\(len\(distance\_matrix\)\)

\]\)

current\_x,current\_y=node\_coords\[current\_node\]

unvisited\_coords=node\_coords\[unvisited\_nodes\]

dx=unvisited\_coords\[:,0\]\-current\_x

dy=unvisited\_coords\[:,1\]\-current\_y

direction\_magnitudes=np\.hypot\(dx,dy\)

direction\_magnitudes=np\.where\(direction\_magnitudes==0,1,direction\_magnitudes\)

direction\_angles=np\.arctan2\(dy,dx\)

alpha=0\.2

global\_trend\_angle=np\.mean\(direction\_angles\)

global\_trend\_angle=alpha\*np\.mean\(direction\_angles\)\+\(1\-alpha\)\*global\_trend\_angle

angle\_diffs=direction\_angles\-global\_trend\_angle

angle\_diffs=np\.abs\(angle\_diffs\)

angle\_diffs=np\.minimum\(angle\_diffs,2\*np\.pi\-angle\_diffs\)

directional\_score=1/\(1\+angle\_diffs\)

weighted\_distances=distances/\(1\+directional\_score\)

unvisited\_distance\_matrix=distance\_matrix\[unvisited\_nodes\[:,None\],unvisited\_nodes\]

radius\_threshold=0\.5

distances\_to\_neighbors=unvisited\_distance\_matrix

neighbor\_count=np\.sum\(distances\_to\_neighbors<=radius\_threshold,axis=1\)

clustering\_bonus=\(neighbor\_count\+1\)/\(len\(unvisited\_nodes\)\+1\)

weighted\_distances\*=clustering\_bonus

ifdestination\_nodeinunvisited\_nodes:

dest\_dist=distance\_matrix\[current\_node,destination\_node\]

dest\_penalty=np\.zeros\(len\(unvisited\_nodes\)\)

foridx,nodeinenumerate\(unvisited\_nodes\):

ifnode\!=destination\_node:

dist\_to\_dest=distance\_matrix\[current\_node,node\]

dist\_diff=abs\(dist\_to\_dest\-dest\_dist\)

penalty=np\.exp\(\-0\.5\*\(dist\_diff/1\.5\)\*\*2\)

dest\_penalty\[idx\]=penalty

else:

dest\_penalty\[idx\]=1\.0

weighted\_distances\+=dest\_penalty

next\_node\_idx=np\.argmin\(weighted\_distances\)

next\_node=unvisited\_nodes\[next\_node\_idx\]

returnnext\_node

#### 0\.C\.1\.5LHS \(Ours\)

importnumpyasnp

defselect\_next\_node\(current\_node:int,destination\_node:int,unvisited\_nodes:np\.ndarray,distance\_matrix:np\.ndarray\)\-\>int:

ifunvisited\_nodes\.size==0:

returnint\(destination\_node\)

cand=unvisited\_nodes\.astype\(int\)

r=distance\_matrix\[current\_node,cand\]

d=distance\_matrix\[np\.ix\_\(cand,cand\)\]\.copy\(\)

n=int\(np\.sum\(unvisited\_nodes==cand\)\)

d=np\.exp\(\-d/\(1e\-6\+np\.max\(d\)\+1e\-6\)\)

score=\(r\+\(d\.sum\(axis=1\)\*\(1\.0/\(n\+1e\-6\)\)\)\)\.clip\(0,None\)

noise=np\.arange\(len\(score\)\)\*1e\-9

score\+=noise

returnint\(cand\[int\(np\.argmin\(score\)\)\]\)

### 0\.C\.2Capacitated Vehicle Routing Problem \(CVRP\)

#### 0\.C\.2\.1Evolution of Heuristics \(EoH\)

importnumpyasnp

defselect\_next\_node\(current\_node:int,depot:int,unvisited\_nodes:np\.ndarray,rest\_capacity:np\.ndarray,demands:np\.ndarray,distance\_matrix:np\.ndarray\)\-\>int:

"""Designanovelalgorithmtoselectthenextnodeineachstep\."""

iflen\(unvisited\_nodes\)==0:

return\-1

distances=distance\_matrix\[current\_node\]\[unvisited\_nodes\]

demands\_arr=demands\[unvisited\_nodes\]

feasible\_mask=demands\_arr<=rest\_capacity

fromscipy\.statsimportentropy

fromnumpyimportlog

demand\_bins=np\.histogram\(demands\_arr,bins=5,range=\(0,np\.max\(demands\_arr\)\+1\)\)\[1\]

demand\_groups=np\.digitize\(demands\_arr,demand\_bins\)

prior\_prob=np\.bincount\(demand\_groups,minlength=5\)/len\(unvisited\_nodes\)

prior\_entropy=\-np\.sum\(prior\_prob\*np\.log\(prior\_prob\+1e\-6\)\)ifnp\.any\(prior\_prob\>0\)else0

node\_entropy\_change=\[\]

fori,node\_idinenumerate\(unvisited\_nodes\):

ifnotfeasible\_mask\[i\]:

continue

new\_demands=demands\_arr\.copy\(\)

new\_demands\[i\]=0

new\_demand\_groups=np\.digitize\(new\_demands,demand\_bins\)

new\_prob=np\.bincount\(new\_demand\_groups,minlength=5\)/\(len\(unvisited\_nodes\)\-1\)

new\_entropy=\-np\.sum\(new\_prob\*np\.log\(new\_prob\+1e\-6\)\)ifnp\.any\(new\_prob\>0\)else0

entropy\_reduction=prior\_entropy\-new\_entropy

node\_entropy\_change\.append\(entropy\_reduction\)

iflen\(node\_entropy\_change\)==0:

feasible\_mask=feasible\_mask&\(distances<=np\.max\(distances\)\*0\.9\)

ifnotnp\.any\(feasible\_mask\):

returnunvisited\_nodes\[np\.argmin\(distances\)\]

proximity\_score=1/\(distances\[feasible\_mask\]\+1e\-6\)

utilization\_score=np\.clip\(demands\_arr\[feasible\_mask\]/\(rest\_capacity\+1e\-6\),0\.2,0\.6\)

final\_score=proximity\_score\*0\.5\+utilization\_score\*0\.3

returnunvisited\_nodes\[np\.argmax\(final\_score\)\]

proximity\_score=1/\(distances\[feasible\_mask\]\+1e\-6\)

utilization\_score=np\.clip\(demands\_arr\[feasible\_mask\]/\(rest\_capacity\+1e\-6\),0\.2,0\.6\)

entropy\_score=np\.array\(node\_entropy\_change\)\[feasible\_mask\]

final\_score=proximity\_score\*0\.4\+utilization\_score\*0\.3\+entropy\_score\*0\.3

best\_idx=np\.argmax\(final\_score\)

returnunvisited\_nodes\[feasible\_mask\]\[best\_idx\]

#### 0\.C\.2\.2FunSearch

importnumpyasnp

defselect\_next\_node\(current\_node:int,depot:int,unvisited\_nodes:np\.ndarray,rest\_capacity:np\.ndarray,demands:np\.ndarray,distance\_matrix:np\.ndarray\)\-\>int:

"""Designanovelalgorithmtoselectthenextnodeineachstep\."""

best\_score=\-1

next\_node=\-1

unvisited\_distances=distance\_matrix\[current\_node\]\[unvisited\_nodes\]

unvisited\_demands=demands\[unvisited\_nodes\]

FLEX\_THRESHOLD=0\.1

DEPOT\_DISTANCE=distance\_matrix\[current\_node\]\[depot\]

DEMAND\_THRESHOLD=0\.5\*np\.max\(unvisited\_demands\)ifunvisited\_demands\.size\>0else0

fori,nodeinenumerate\(unvisited\_nodes\):

demand=unvisited\_demands\[i\]

distance=unvisited\_distances\[i\]

ifdemand\>rest\_capacity:

continue

ifdistance==0:

ifnode==depot:

score=0

else:

score=float\(’inf’\)

base\_score=score

else:

base\_score=demand/distance

future\_capacity=rest\_capacity\-demand

iffuture\_capacity<FLEX\_THRESHOLD\*demand:

base\_score\*=2\.5

ifdistance<DEPOT\_DISTANCE:

base\_score\*=1\.15

distance\_threshold=0\.7\*DEPOT\_DISTANCE

similar\_nodes=\[\]

forj,other\_nodeinenumerate\(unvisited\_nodes\):

ifother\_node\!=node:

d=distance\_matrix\[current\_node\]\[other\_node\]

ifd<distance\_thresholdandunvisited\_demands\[j\]\>0\.5\*demand:

similar\_nodes\.append\(other\_node\)

iflen\(similar\_nodes\)\>1:

base\_score\*=0\.85

ifdemand<DEMAND\_THRESHOLDanddistance<DEPOT\_DISTANCE\*0\.6:

base\_score\*=1\.2

dep\_dist=distance\_matrix\[node\]\[depot\]

ifdep\_dist\>0\.8\*DEPOT\_DISTANCE:

base\_score\*=1\.2

ifnext\_node==\-1:

closest\_dist=np\.inf

closest\_node=\-1

fornodeinunvisited\_nodes:

d=distance\_matrix\[current\_node\]\[node\]

ifd<closest\_dist:

closest\_dist=d

closest\_node=node

returnclosest\_nodeifclosest\_node\!=\-1elseunvisited\_nodes\[0\]

returnnext\_node

#### 0\.C\.2\.3MCTS\-AHD

importnumpyasnp

defselect\_next\_node\(current\_node:int,depot:int,unvisited\_nodes:np\.ndarray,rest\_capacity:np\.ndarray,demands:np\.ndarray,distance\_matrix:np\.ndarray\)\-\>int:

"""Designanovelalgorithmtoselectthenextnodeineachstep\."""

iflen\(unvisited\_nodes\)==0:

return\-1

n=len\(unvisited\_nodes\)

time\_to\_serve=np\.zeros\(n\)

capacity\_gap\_penalty=np\.zeros\(n\)

fori,nodeinenumerate\(unvisited\_nodes\):

dist=distance\_matrix\[current\_node\]\[node\]

demand\_util=demands\[node\]/\(rest\_capacity\+1e\-6\)

urgency=demands\[node\]\*1\.5

demand\_list=\[demands\[n\]forninunvisited\_nodes\]

median\_demand=np\.median\(demand\_list\)

capacity\_gap=max\(0,demands\[node\]\-median\_demand\)

capacity\_gap\_penalty\[i\]=capacity\_gap\*2\.0

time\_to\_serve\[i\]=dist\*\(1\+urgency\+capacity\_gap\_penalty\[i\]\)

threshold=np\.percentile\(time\_to\_serve,75\)

eligible\_nodes=unvisited\_nodes\[time\_to\_serve<threshold\]

iflen\(eligible\_nodes\)\>0:

best\_idx=np\.argmin\(\[distance\_matrix\[current\_node\]\[node\]fornodeineligible\_nodes\]\)

returneligible\_nodes\[best\_idx\]

best\_idx=np\.argmin\(\[distance\_matrix\[current\_node\]\[node\]fornodeinunvisited\_nodes\]\)

returnunvisited\_nodes\[best\_idx\]

#### 0\.C\.2\.4ReEvo

importnumpyasnp

defselect\_next\_node\(current\_node:int,depot:int,unvisited\_nodes:np\.ndarray,rest\_capacity:np\.ndarray,demands:np\.ndarray,distance\_matrix:np\.ndarray\)\-\>int:

"""Designanovelalgorithmtoselectthenextnodeineachstep\."""

best\_score=\-1

next\_node=\-1

fornodeinunvisited\_nodes:

demand=demands\[node\]

distance=distance\_matrix\[current\_node\]\[node\]

ifdistance==0:

continue

ifdemand\>rest\_capacity:

continue

capacity\_utilization=\(rest\_capacity\-demand\)/rest\_capacityifrest\_capacity\>0else0

efficiency\_score=demand/\(distance\+1e\-8\)

capacity\_bias=1\.5\*capacity\_utilization\+0\.3\*\(rest\_capacity/\(demand\+1e\-6\)\)

score=efficiency\_score\*\(1\.0\+capacity\_bias\)

ifscore\>best\_score:

best\_score=score

next\_node=node

returnnext\_nodeifnext\_node\!=\-1elseunvisited\_nodes\[0\]iflen\(unvisited\_nodes\)\>0elsedepot

#### 0\.C\.2\.5LHS \(Ours\)

importnumpyasnp

defselect\_next\_node\(current\_node:int,depot:int,unvisited\_nodes:np\.ndarray,rest\_capacity:float,demands:np\.ndarray,distance\_matrix:np\.ndarray\)\-\>int:

iflen\(unvisited\_nodes\)==0:

returndepot

cap=float\(rest\_capacity\)

feasible=np\.array\(\[nforninunvisited\_nodesifdemands\[n\]<=cap\],dtype=int\)

iflen\(feasible\)==0:

returndepot

d\_cur=distance\_matrix\[current\_node,feasible\]

d\_dep=distance\_matrix\[depot,feasible\]

d\_cur\_adj=d\_cur\-d\_dep

d\_cur\_adj=np\.clip\(d\_cur\_adj,\-np\.inf,1e12\)

d\_cur\_adj=np\.clip\(d\_cur\_adj,\-np\.inf,1e12\)

returnint\(feasible\[np\.argmin\(d\_cur\_adj\)\]\)

### 0\.C\.3Knapsack Problem \(KSP\)

#### 0\.C\.3\.1Evolution of Heuristics \(EoH\)

defselect\_next\_item\(remaining\_capacity:int,remaining\_items:List\[Tuple\[int,int,int\]\]\)\-\>Tuple\[int,int,int\]\|None:

"""

Selecttheitemwiththehighestvalue\-to\-weightratiothatfitsintheremainingcapacity\.

Args:

remaining\_capacity:Theremainingcapacityoftheknapsack\.

remaining\_items:Listoftuplescontaining\(weight,value,index\)ofremainingitems\.

Returns:

Theselecteditemasatuple\(weight,value,index\),orNoneifnoitemfits\.

"""

ifnotremaining\_items:

returnNone

best\_item=None

max\_ratio=\-1

forweight,value,indexinremaining\_items:

ifweight<=remaining\_capacity:

ratio=value/weight

ifratio\>max\_ratio:

max\_ratio=ratio

best\_item=\(weight,value,index\)

returnbest\_itemifbest\_itemelseNone

#### 0\.C\.3\.2FunSearch

defselect\_next\_item\(remaining\_capacity:int,remaining\_items:List\[Tuple\[int,int,int\]\]\)\-\>Tuple\[int,int,int\]\|None:

"""

Selecttheitemwiththehighestvalue\-to\-weightratiothatfitsintheremainingcapacity\.

Args:

remaining\_capacity:Theremainingcapacityoftheknapsack\.

remaining\_items:Listoftuplescontaining\(weight,value,index\)ofremainingitems\.

Returns:

Theselecteditemasatuple\(weight,value,index\),orNoneifnoitemfits\.

"""

best\_item=None

best\_ratio=\-1

foriteminremaining\_items:

weight,value,index=item

ifweight<=remaining\_capacity:

ratio=value/weight

ifratio\>best\_ratio:

best\_ratio=ratio

best\_item=item

returnbest\_item

#### 0\.C\.3\.3MCTS\-AHD

defselect\_next\_item\(remaining\_capacity:int,remaining\_items:List\[Tuple\[int,int,int\]\]\)\-\>Tuple\[int,int,int\]\|None:

"""

Selecttheitemwiththehighestvalue\-to\-weightratiothatfitsintheremainingcapacity\.

Args:

remaining\_capacity:Theremainingcapacityoftheknapsack\.

remaining\_items:Listoftuplescontaining\(weight,value,index\)ofremainingitems\.

Returns:

Theselecteditemasatuple\(weight,value,index\),orNoneifnoitemfits\.

"""

ifnotremaining\_items:

returnNone

best\_item=None

best\_ratio=\-1

adjustment\_factor=1\.0

ifremaining\_capacity<=10:

adjustment\_factor=1\.5

elifremaining\_capacity<=20:

adjustment\_factor=1\.2

forweight,value,indexinremaining\_items:

ifweight<=remaining\_capacity:

ratio=\(value/weight\)\*adjustment\_factor\+\(value/100\)

ifratio\>best\_ratio:

best\_ratio=ratio

best\_item=\(weight,value,index\)

returnbest\_item

#### 0\.C\.3\.4ReEvo

defselect\_next\_item\(remaining\_capacity:int,remaining\_items:List\[Tuple\[int,int,int\]\]\)\-\>Tuple\[int,int,int\]\|None:

"""

Selecttheitemwiththehighestvalue\-to\-weightratiothatfitsintheremainingcapacity\.

Args:

remaining\_capacity:Theremainingcapacityoftheknapsack\.

remaining\_items:Listoftuplescontaining\(weight,value,index\)ofremainingitems\.

Returns:

Theselecteditemasatuple\(weight,value,index\),orNoneifnoitemfits\.

"""

ifnotremaining\_items:

returnNone

item\_freq=\{\}

foriteminremaining\_items:

index=item\[2\]

item\_freq\[index\]=item\_freq\.get\(index,0\)\+1

underrepresented=\[itemforiteminremaining\_itemsifitem\_freq\[item\[2\]\]<2\]

total\_items=len\(remaining\_items\)

capacity\_threshold=remaining\_capacity

iftotal\_items\>5andcapacity\_threshold\>30:

best\_item=None

best\_ratio=\-1

best\_value=\-1

foriteminremaining\_items:

weight,value,index=item

ifweight<=remaining\_capacity:

ratio=value/weight

if\(ratio\>best\_ratioor

\(ratio==best\_ratioandvalue\>best\_value\)\):

best\_ratio=ratio

best\_value=value

best\_item=item

returnbest\_item

eliftotal\_items\>3andcapacity\_threshold\>10:

candidates=\[\]

foriteminremaining\_items:

weight,value,index=item

ifweight<=remaining\_capacity:

weight\_penalty=0\.0

ifweight\>remaining\_capacity\*0\.6:

weight\_penalty=1\.5

diversity\_factor=1\.0ifitem\_freq\.get\(index,0\)<2else0\.7

ratio=value/weight

score=\(ratio\*value\*diversity\_factor\)\-weight\_penalty

candidates\.append\(\(score,weight,value,index\)\)

ifnotcandidates:

returnNone

candidates\.sort\(key=lambdax:x\[0\],reverse=True\)

selected=candidates\[0\]

return\(selected\[1\],selected\[2\],selected\[3\]\)

else:

candidates=\[\]

foriteminremaining\_items:

weight,value,index=item

ifweight<=remaining\_capacity:

weight\_penalty=0\.0

ifweight\>remaining\_capacity\*0\.4:

weight\_penalty=2\.0

diversity\_factor=1\.5ifitem\_freq\.get\(index,0\)<2else0\.5

ratio=value/weight

score=\(ratio\*value\*diversity\_factor\)\-weight\_penalty

candidates\.append\(\(score,weight,value,index\)\)

ifnotcandidates:

returnNone

candidates\.sort\(key=lambdax:x\[0\],reverse=True\)

selected=candidates\[0\]

return\(selected\[1\],selected\[2\],selected\[3\]\)

#### 0\.C\.3\.5LHS \(Ours\)

importnumpyasnp

fromtypingimportList,Tuple,Optional

defselect\_next\_item\(remaining\_capacity:int,remaining\_items:List\[Tuple\[int,int,int\]\]\)\-\>Optional\[Tuple\[int,int,int\]\]:

fit=\[\(w,v,idx\)for\(w,v,idx\)inremaining\_itemsifw<=remaining\_capacity\]

ifnotfit:

returnNone

forw,v,idxinfit:

ifw==0:

return\(w,v,idx\)

densities=\[float\(v\)/max\(1e\-12,float\(w\)\)for\(w,v,idx\)infit\]

median=np\.median\(densities\)

mean=np\.mean\(densities\)

best\_score=\-np\.inf

best\_item=None

forw,v,idxinfit:

ifw==0:

score=v

else:

score=\(v/max\(1e\-12,float\(w\)\)\)\*\(1\.0\+0\.5\*\(w/remaining\_capacity\)\)

ifscore\>best\_score:

best\_score=score

best\_item=\(w,v,idx\)

returnbest\_item

### 0\.C\.4Online Bin Packing \(OBP\)

#### 0\.C\.4\.1Evolution of Heuristics \(EoH\)

importnumpyasnp

defpriority\(item:float,bins:np\.ndarray\)\-\>np\.ndarray:

"""Returnsprioritywithwhichwewanttoadditemtoeachbin\.

Args:

item:Sizeofitemtobeaddedtothebin\.

bins:Arrayofcapacitiesforeachbin\.

Return:

Arrayofsamesizeasbinswithpriorityscoreofeachbin\.

"""

return1/\(bins\+1e\-8\)/\(bins/\(bins\+item\)\)

#### 0\.C\.4\.2FunSearch

importnumpyasnp

defpriority\(item:float,bins:np\.ndarray\)\-\>np\.ndarray:

"""Returnsprioritywithwhichwewanttoadditemtoeachbin\.

Args:

item:Sizeofitemtobeaddedtothebin\.

bins:Arrayofcapacitiesforeachbin\.

Return:

Arrayofsamesizeasbinswithpriorityscoreofeachbin\.

"""

returnitem\-bins

#### 0\.C\.4\.3MCTS\-AHD

importnumpyasnp

defpriority\(item:float,bins:np\.ndarray\)\-\>np\.ndarray:

"""Returnsprioritywithwhichwewanttoadditemtoeachbin\.

Args:

item:Sizeofitemtobeaddedtothebin\.

bins:Arrayofcapacitiesforeachbin\.

Return:

Arrayofsamesizeasbinswithpriorityscoreofeachbin\.

"""

ifitem<=0:

returnnp\.ones\_like\(bins\)\*1e\-8

scaled\_bins=bins/item

max\_ratio=1\.4

min\_ratio=0\.7

base\_width=0\.4

sensitivity\_factor=1\.0\+0\.5\*np\.log1p\(1\.0/\(1\.0\+np\.exp\(\-2\*\(item\-1\.0\)\)\)\)

width=base\_width\*sensitivity\_factor

lorentzian\_kernel=1\.0/\(1\.0\+\(\(scaled\_bins\-1\.0\)/width\)\*\*2\)

threshold\_factor=0\.8

lorentzian\_kernel=lorentzian\_kernel\*threshold\_factor

phase\_shift=np\.pi\*\(scaled\_bins\-0\.5\)

cosine\_modulation=np\.cos\(phase\_shift\+np\.pi/6\)\*0\.15

modulated=lorentzian\_kernel\+cosine\_modulation

median\_ratio=np\.median\(scaled\_bins\)

entropy\_term=np\.clip\(1\.0/\(1\.0\+0\.4\*\(np\.abs\(scaled\_bins\-median\_ratio\)/0\.4\)\*\*2\),0\.4,1\.0\)

final=modulated\*entropy\_term

returnnp\.clip\(final\+1e\-8,1e\-8,1\.0\)

#### 0\.C\.4\.4ReEvo

importnumpyasnp

defpriority\(item:float,bins:np\.ndarray\)\-\>np\.ndarray:

"""Returnsprioritywithwhichwewanttoadditemtoeachbin\.

Args:

item:Sizeofitemtobeaddedtothebin\.

bins:Arrayofcapacitiesforeachbin\.

Return:

Arrayofsamesizeasbinswithpriorityscoreofeachbin\.

"""

residual\_capacity=bins\-item

feasible\_residues=np\.clip\(residual\_capacity,0,None\)

ifnp\.any\(feasible\_residues\>0\):

inverse\_residues=1\.0/\(feasible\_residues\+1e\-8\)

normalized\_priority=inverse\_residues/inverse\_residues\.sum\(\)

priority=normalized\_priority

else:

priority=np\.zeros\_like\(bins\)

returnpriority

#### 0\.C\.4\.5LHS \(Ours\)

importnumpyasnp

defpriority\(item:float,bins:np\.ndarray\)\-\>np\.ndarray:

a=bins\.copy\(\)\.astype\(float\)

a\[a<item\]=\-np\.inf

return\-np\.clip\(a,0,None\)