Algorithms

BCA(;N, n, K, η_max, resize_population)

Bilevel Centers Algorithm — a physics-inspired metaheuristic for single-objective bilevel optimisation. It uses a nested scheme where both the upper and lower levels employ a center-of-mass variation operator.

How it works

At each iteration, for every individual in the upper-level population:

1. A random subset U of size K is selected.
2. The center of mass c of U is computed, weighted by the combined fitness Q = F + f (i.e., both level objectives).
3. A new candidate p = x_i + η · (c - u_worst) is generated, where u_worst is the worst element in U and η is a random step size.
4. The lower-level problem is solved for p (using the same center-of-mass strategy), producing one or more optimal lower-level responses.
5. The new pair (p, y_opt) replaces a worse member of the population.

The population size may be dynamically reduced during the run (resize_population), shifting from exploration to exploitation.

Parameters

• N — upper-level population size (auto‑computed from K × D_ul if 0).
• n — lower-level population size (auto‑computed if 0).
• K — number of solutions used to compute the center of mass (default 7). Larger K → faster convergence (exploitation); smaller K → more exploration.
• η_max — maximum step size for the variation operator (default 2.0).
• resize_population — if true, the population shrinks linearly over the run (default true).

Reference

Mejía-de-Dios, J. A., & Mezura-Montes, E. (2018, November). A physics-inspired algorithm for bilevel optimization. In 2018 IEEE International Autumn Meeting on Power, Electronics and Computing (ROPEC) (pp. 1–6). IEEE.

QBCA(;N, K, η_ul_max, η_ll_max, α, β, autodiff)

Quasi-Newton Bilevel Centers Algorithm — combines the center-of-mass variation (upper level) with a BFGS quasi-Newton solver (lower level) and Tikhonov regularisation for improved lower-level accuracy.

How it works

The upper-level search follows the same center-of-mass strategy as BCA. The key difference is in the lower-level solver:

1. For a given upper-level x, the lower level is first explored with a lightweight ECA centre-of-mass search to get an approximate y.
2. That y is then refined by a BFGS quasi-Newton method, minimising a Tikhonov- regularised objective f(x, y) + α·‖y‖² + β·‖y - y₀‖², which stabilises the lower-level solution and prevents overfitting to a single x.
3. Both the upper and lower levels compute their own centers of mass, and the lower-level step also uses a separate step size η_ll_max.

This hybrid approach reduces the number of lower-level function evaluations compared to BCA, especially on problems where the lower level is smooth and unimodal.

Parameters

• N — upper-level population size (auto‑computed if 0).
• K — number of solutions to generate centers (default 3).
• η_ul_max — upper-level step size (default 2.0).
• η_ll_max — lower-level step size (default 1/η_ul_max).
• α, β — Tikhonov regularisation weights (default 0.05 each).
• autodiff — differentiation mode for BFGS: :finite (finite differences) or :forward (forward-mode AD) (default :finite).

Reference

Mejía-de-Dios, J. A., & Mezura-Montes, E. (2019, June). A metaheuristic for bilevel optimization using Tikhonov regularization and the quasi-Newton method. In 2019 IEEE Congress on Evolutionary Computation (CEC) (pp. 3134–3141). IEEE.

QBCA2(;N, K, η_max, δ1, δ2, ε1, ε2, λ, t_reevaluate)

Improved Quasi-Newton BCA — an enhanced version of QBCA that explicitly detects and avoids pseudo-feasible solutions.

How it works

A pseudo-feasible solution (x, y) is one where:

• y is optimal for the lower level given x, and
• F(x, y) is close to the upper-level optimum, but
• the lower-level optimum is not unique (multiple different y give the same f value).

Such solutions can mislead the upper-level search because the leader observes a good F value but the follower's response is not stable. QBCA2 detects pseudo-feasible pairs using thresholds (δ1, δ2, ε1, ε2) and avoids replacing solutions with unstable ones.

The upper-level search uses the center-of-mass operator, while the lower level is solved with BFGS (optionally preceded by Nelder-Mead). A surrogate model (SECA) can optionally be enabled to further reduce lower-level evaluations.

Parameters

• N — upper-level population size (auto‑computed if 0).
• K — number of solutions to generate centers (default 3).
• η_max — maximum step size (default 2.0).
• δ1, δ2 — position difference thresholds for pseudo-feasibility detection (default 1e-2).
• ε1, ε2 — objective difference thresholds for pseudo-feasibility detection (default 1e-2).
• t_reevaluate — frequency (in iterations) for re-evaluating the lower level of the elite solution (default 10).
• λ — regularisation for the optional surrogate model (default 1e-12).
• autodiff — differentiation mode for BFGS (:finite or :forward).
• use_surrogate_model — if true, a kernel-interpolation surrogate (SECA) assists the lower-level search (default false).

Reference

Mejía-de-Dios, J. A., Mezura-Montes, E., & Toledo-Hernández, P. (2022). Pseudo-feasible solutions in evolutionary bilevel optimization: Test problems and performance assessment. Applied Mathematics and Computation, 412, 126577.

SABO(;N, K, η_max, δ1, δ2, ε1, ε2, λ, t_reevaluate)

Surrogate-Assisted Bilevel Optimization — reduces the cost of lower-level optimisation by employing a kernel-interpolation surrogate model (BiApprox) to approximate the lower-level optimal response.

How it works

SABO extends the QBCA2 framework with a surrogate strategy:

1. An elite surrogate solution is maintained, tracking the most promising pair (x, y).
2. At each iteration, the surrogate model predicts the lower-level response for candidate upper-level points, reducing expensive lower-level evaluations.
3. Every t_reevaluate iterations, the entire population is re-evaluated with the true lower-level solver (BFGS) to correct surrogate drift.
4. The same pseudo-feasibility detection from QBCA2 is used to avoid unstable solutions.

The surrogate model is a kernel interpolation (radial basis function with Gaussian kernel) trained on previously evaluated (x, y) pairs. This makes SABO particularly effective when the lower level is expensive to evaluate.

Parameters

• N — upper-level population size (auto‑computed if 0).
• K — number of solutions to generate centers (default 3).
• η_max — maximum step size (default 2.0).
• δ1, δ2, ε1, ε2 — thresholds for pseudo-feasibility detection (default 1e-2).
• λ — regularisation for the kernel interpolation surrogate (default 1e-5).
• t_reevaluate — re-evaluation period (default 5 iterations).
• autodiff — differentiation mode (:finite or :forward).

Reference

Mejía-de-Dios, J. A., & Mezura-Montes, E. (2020, June). A surrogate-assisted metaheuristic for bilevel optimization. In Proceedings of the 2020 Genetic and Evolutionary Computation Conference (GECCO) (pp. 629–635).

BLEMO(; ul = NSGA2(), ll = NSGA2())

Bilevel Evolutionary Multi-Objective Optimization — a nested evolutionary algorithm for multi-objective bilevel problems.

BLEMO maintains multiple subpopulations, each corresponding to a distinct upper-level decision x. For each x, the lower-level NSGA-II finds a set of optimal responses. New upper-level candidates are generated via SBX crossover and polynomial mutation on the current elite set, and the lower level re-optimises from the inherited subpopulations. Non-dominated sorting and crowding distance are used for environmental selection at both levels.

The Pareto archive of non-dominated bilevel solutions is stored in method.parameters.archive.

Parameters (keyword constructor)

• ul — upper-level NSGA2 configuration.
• ll — lower-level NSGA2 configuration.

SMS_MOBO(;
    ul = Metaheuristics.SMS_EMOA(;N = 100),  # upper level optimizer
    ll = Metaheuristics.NSGA2(;N = 50),      # lower level optimizer
    ul_offsprings = 10,
    options_ul = Options(iterations = 100, f_calls_limit = Inf),
    options_ll = Options(iterations = 50, f_calls_limit = Inf),
)

S-metric-selection Multi-Objective Bilevel Optimization — described in Mejía-de-Dios & Mezura-Montes (2022).

SMS_MOBO extends the S-metric selection idea (SMS-EMOA) to the bilevel setting. Both the upper and lower levels solve multi-objective problems:

• The upper level uses SMS_EMOA (hypervolume-based indicator) to drive the search toward a well-distributed Pareto front.
• The lower level uses a multi-objective evolutionary algorithm (NSGA2, SPEA2, or SMS_EMOA) to find the optimal response set for each upper-level decision.
• A family concept groups an upper-level x with its associated lower-level optimal y solutions, and a super-rank criterion combines dominance information from both levels for environmental selection.

The non-dominated upper-level solutions are stored in method.parameters.archive and are also accessible via res.population after optimisation.

Parameters

• ul — upper-level optimizer (SMS_EMOA is the only valid choice).
• ll — lower-level optimizer (NSGA2, SPEA2, or SMS_EMOA).
• ul_offsprings — number of new upper-level candidate solutions generated per generation (default 10).

Example

using BilevelHeuristics

function F(x, y)   # two upper-level objectives
    [y[1] - x[1], y[2]], [-1.0 - sum(y)], [0.0]
end

function f(x, y)   # two lower-level objectives
    y, [-x[1]^2 + sum(y .^ 2)], [0.0]
end

bounds_ul = [0.0 1.0]'
bounds_ll = [-1 -1; 1 1.0]

res = optimize(F, f, bounds_ul, bounds_ll, SMS_MOBO())
# Access the Pareto archive:
archive = SMS_MOBO().parameters.archive

Nested(; ul, ll)

A flexible nested framework for bilevel optimisation that accepts any Metaheuristics.jl algorithm as the upper- and lower-level solver.

Unlike the specialised algorithms (BCA, QBCA, etc.) which use tightly-coupled strategies, Nested decouples the two levels entirely:

1. The upper-level optimizer proposes candidate x values.
2. For each x, the lower-level optimizer runs to completion, finding optimal y.
3. The joint solution (x, y) is evaluated and the upper-level population is updated.

This makes Nested suitable for prototyping — you can experiment with different combinations of algorithms (e.g. GA + DE, PSO + BFGS) without implementing any bilevel-specific logic.

Parameters

• ul — a configured Metaheuristics.AbstractAlgorithm for the upper level.
• ll — a configured Metaheuristics.AbstractAlgorithm for the lower level.

Example

using BilevelHeuristics, Metaheuristics

ul = GA(N = 50)          # Genetic Algorithm at upper level
ll = DE(N = 50)          # Differential Evolution at lower level
res = optimize(F, f, bounds_ul, bounds_ll, Nested(; ul, ll))

DBMA(; Nu, Nl, F, CR)

Differential Evolution Bilevel Multi-objective Algorithm — an ε-constrained DE framework for bilevel problems where the upper level is multi-objective and the lower level is single-objective (or multi-objective).

DBMA extends the ε-constrained Differential Evolution (εDE) to both levels. The ε-level control gradually relaxes constraint violations, allowing the algorithm to explore infeasible regions early and converge to feasible Pareto-optimal solutions later.

Internally, DBMA inherits from the Nested framework and specialises the truncation and decision-making steps to handle multi-objective upper-level populations using ε-dominance.

Parameters

• Nu — upper-level population size (default 50).
• Nl — lower-level population size (default 50).
• F — DE mutation scaling factor (default 0.7).
• CR — DE crossover rate (default 0.5).

Example

using BilevelHeuristics

F(x, y) = [y[1] - x[1], y[2]], [-1.0 - sum(y)], [0.0]   # two objectives
f(x, y) = sum((x - y).^2) + y[1]^2, [0.0], [0.0]        # single objective

bounds_ul = [0.0 1.0]'
bounds_ll = [-ones(5) ones(5)]

res = optimize(F, f, bounds_ul, bounds_ll, DBMA())