Powell
The problem
Powell's singular function is a sum of two squares and two fourth powers of four genes, to minimize; in more dimensions, the sum of the same over blocks of four:
f(x) = Σ over blocks (x₁, x₂, x₃, x₄) of
(x₁ + 10x₂)² + 5 (x₃ − x₄)² + (x₂ − 2x₃)⁴ + 10 (x₁ − x₄)⁴, each xᵢ in [−4, 5]
Its minimum is 0, at the origin. Here n = 24, six blocks.
Powell (1962) defined it in 4 dimensions, with the start (3, −1, 0, 1), where f = 215, and no bounds; the paper couldn't be read. genoxide takes the formula from Steihaug and Suleiman (2013, Journal of Global Optimization 56(3): 845-853), who restate it from Powell, and the extension to blocks of four and the bounds from Laguna and Martí (2005, function 36, with n = 24). Jamil and Yang (2013, function 91) print (x₂ − x₃)⁴ for (x₂ − 2x₃)⁴, and give the start as the minimizer.
The origin is the only point where every term is 0: x₁ = −10x₂, x₃ = x₄, x₂ = 2x₃ and x₁ = x₄ leave
only 0. genoxide's problems::Powell gives it as proven.
What makes it hard
The function is convex, but its Hessian at the minimum is singular: the fourth powers have no curvature at 0, so in each block, two directions (along which x₁ + 10x₂ and x₃ − x₄ stay 0) are flat to second order. There the function falls only as the fourth power of the distance, so an error of 1e-8 needs the genes within about 1e-2 of the minimum along those directions, and within 1e-4 along the others. Methods that rely on a quadratic model converge only linearly, as Steihaug and Suleiman show for Newton's method.
The squares also couple the genes of each block with weights up to 10, so the level sets are long and oblique. The blocks don't interact: a search that learns one block's shape can use it for all six.
Representation
A Real genome of 24 genes, each in [−4, 5]: the point x itself. The fitness is f(x), to minimize.
The function is genoxide's problems::Powell, which brings its bounds and its minimum.
Algorithm
Five algorithms, each with a budget of 10,000 evaluations per dimension, 240,000 in all, and a target of 1e-8, from seed 1:
- CMA-ES (Hansen and Ostermeier, 2001, Evolutionary Computation 9(2): 159-195), which samples a population of 13 from a normal distribution and adapts its mean, its step size and its covariance matrix, from a step size of 0.3 of each gene's range and a random start;
- sep-CMA-ES (Ros and Hansen, 2008, PPSN X: 296-305), the same with a diagonal covariance matrix: a scale per gene but no correlations;
- differential evolution with genoxide's defaults, SHADE (Tanabe and Fukunaga, CEC 2013), with a population of 100;
- particle swarm optimization (Kennedy and Eberhart, 1995), 40 particles with Clerc and Kennedy's constriction coefficients and a global topology;
- a real-coded genetic algorithm: a population of 100, tournaments of 3, simulated binary crossover (Deb and Agrawal, 1995) with η = 15 and polynomial mutation with η = 20 at a rate of 1/24 per gene.
Output
The first two lines give the dimension and the budget. Then a row per algorithm: the evaluations it
had used when its best error first reached 1, 1e-2, 1e-4, 1e-6 and 1e-8, and the best error it
found, to two significant digits. A dash is an error not reached. The function has no sin, cos
or exp, so the runs are the same on every platform, and in Python, run evaluates the function
in Rust, so both versions print the same.
The project page plays back another run: CMA-ES on the function in 4 dimensions, Powell's own, each gene on its range with the minimum marked.
Good results
The minimum is 0. CMA-ES reaches 1e-8 after 19,864 evaluations: 1,716 to get down to an error of 1, then about 2,300 per decade. Its full covariance matrix learns the oblique squares, and its step size keeps shrinking along the flat directions.
sep-CMA-ES is faster at first, reaching 1e-2 after 2,665 evaluations, but then slows: it reaches 1e-6 only after 104,260 evaluations and ends at 1.7e-7, since a diagonal matrix can't line up with the oblique squares. SHADE reaches 1e-8 after 79,400 evaluations. PSO reaches 1e-4 after 59,640 and ends at 4.7e-6, and the genetic algorithm reaches 1e-2 after 119,908 and ends at 3.3e-3: along the flat directions, a small error still means genes far from 0, and their steps don't shrink with the error.
Known optimum: 0 (at the origin)
Source: examples/powell
Interactive run: tachsin.gr/projects/genoxide/examples/powell
cargo run --release --example powell
//! Powell: minimize Powell's singular function in 24 dimensions, six blocks of four genes, convex,
//! with a singular Hessian at the minimum.
//!
//! Compares how fast CMA-ES, with a full and with a diagonal covariance matrix (sep-CMA-ES),
//! differential evolution, particle swarm optimization and a real-coded genetic algorithm close in
//! on the minimum, 0 at the origin: the evaluations each takes until its error is at most 1, 1e-2,
//! 1e-4, 1e-6 and 1e-8. The function is genoxide's `problems::Powell`.
//!
//! With `GENOXIDE_TRACE=<file>`, it also writes a trace of a run for the plot on the example's
//! page, with `trace.rs`.
//!
//! ```text
//! cargo run --release --example powell
//! ```
mod trace;
use genoxide::observer::Snapshot;
use genoxide::prelude::*;
use genoxide::problems::{Powell, Problem};
const DIMENSIONS: usize = 24;
const BUDGET: u64 = 10_000 * DIMENSIONS as u64;
// the errors at which the table gives each run's evaluations
const ERRORS: [f64; 5] = [1e0, 1e-2, 1e-4, 1e-6, 1e-8];
const COLUMNS: [&str; 5] = ["1", "1e-2", "1e-4", "1e-6", "1e-8"];
fn main() -> Result<()> {
let problem = Powell::new(DIMENSIONS);
let minimum = problem.optimum().expect("known").value();
let stop = || Stop::target(minimum + 1e-8).or(Stop::evaluations(BUDGET));
println!("Powell in {DIMENSIONS} dimensions, {BUDGET} evaluations at most");
println!("Evaluations until the error is at most");
print!("{:<10}", "algorithm");
COLUMNS.iter().for_each(|column| print!("{column:>9}"));
println!("{:>9}", "best");
for (name, covariance) in [
("CMA-ES", cmaes::Covariance::Full),
("sep-CMA-ES", cmaes::Covariance::Diagonal),
] {
let cmaes = Cmaes::builder(problem.representation())
.covariance(covariance)
.minimize()
.seed(1)
.build()?;
let mut reached = Reached::new(minimum);
let outcome = Engine::new(cmaes, problem)
.stop_when(stop())
.on_generation(|snapshot| reached.record(snapshot))
.run()?;
reached.print(name, &outcome);
}
let de = De::builder(problem.representation())
.minimize()
.seed(1)
.build()?;
let mut reached = Reached::new(minimum);
let outcome = Engine::new(de, problem)
.stop_when(stop())
.on_generation(|snapshot| reached.record(snapshot))
.run()?;
reached.print("DE", &outcome);
let pso = Pso::builder(problem.representation())
.population_size(40)
.minimize()
.seed(1)
.build()?;
let mut reached = Reached::new(minimum);
let outcome = Engine::new(pso, problem)
.stop_when(stop())
.on_generation(|snapshot| reached.record(snapshot))
.run()?;
reached.print("PSO", &outcome);
let ga = Ga::builder(problem.representation())
.population_size(100)
.select(Tournament::new(3)?)
.crossover(SimulatedBinaryCrossover::new(15.0)?)
.mutate(PolynomialMutation::per_gene(1.0 / DIMENSIONS as f64, 20.0)?)
.minimize()
.seed(1)
.build()?;
let mut reached = Reached::new(minimum);
let outcome = Engine::new(ga, problem)
.stop_when(stop())
.on_generation(|snapshot| reached.record(snapshot))
.run()?;
reached.print("GA", &outcome);
// with GENOXIDE_TRACE=<file>, a trace for the plot on the example's page, of a separate
// run in 4 dimensions, Powell's own
trace::record_small()?;
Ok(())
}
// the evaluations after the first generation whose best error was at most each of ERRORS, for a
// function whose minimum is `minimum`
struct Reached {
minimum: f64,
evaluations: [Option<u64>; 5],
}
impl Reached {
fn new(minimum: f64) -> Self {
let evaluations = [None; 5];
Self {
minimum,
evaluations,
}
}
fn record(&mut self, snapshot: &Snapshot<'_, Reals>) {
let progress = snapshot.progress();
let Some(best) = progress.best().and_then(Fitness::score) else {
return;
};
let error = best - self.minimum;
for (reached, bound) in self.evaluations.iter_mut().zip(ERRORS) {
if reached.is_none() && error <= bound {
*reached = Some(progress.evaluations());
}
}
}
// a row of the table: the evaluations, "-" for an error not reached, and the best error
fn print(&self, name: &str, outcome: &Outcome<Reals>) {
print!("{name:<10}");
for reached in self.evaluations {
let reached = reached.map_or("-".to_string(), |evaluations| evaluations.to_string());
print!("{reached:>9}");
}
// rounding can put a solution a few ulps below the minimum
let best = outcome.best_fitness().score().expect("valid");
let error = (best - self.minimum).max(0.0);
println!("{:>9}", format!("{error:.1e}"));
}
}
python examples/powell/main.py
"""Powell: minimize Powell's singular function in 24 dimensions, six blocks of four genes, convex,
with a singular Hessian at the minimum.
Compares how fast CMA-ES, with a full and with a diagonal covariance matrix (sep-CMA-ES),
differential evolution, particle swarm optimization and a real-coded genetic algorithm close in on
the minimum, 0 at the origin: the evaluations each takes until its error is at most 1, 1e-2, 1e-4,
1e-6 and 1e-8. The function is genoxide's `problems::Powell`.
With ``GENOXIDE_TRACE=<file>``, it also writes a trace of a run for the plot on the example's page,
with trace.py.
python examples/powell/main.py
"""
import genoxide as gx
from trace import record_small
DIMENSIONS = 24
BUDGET = 10_000 * DIMENSIONS
# the errors at which the table gives each run's evaluations
ERRORS = [1e0, 1e-2, 1e-4, 1e-6, 1e-8]
COLUMNS = ["1", "1e-2", "1e-4", "1e-6", "1e-8"]
class Reached:
"""The evaluations after the first generation whose best error was at most each of ERRORS,
for a function whose minimum is ``minimum``."""
def __init__(self, minimum):
self.minimum = minimum
self.evaluations = [None] * len(ERRORS)
def record(self, progress):
if progress.best_fitness is None:
return
error = progress.best_fitness - self.minimum
for i, bound in enumerate(ERRORS):
if self.evaluations[i] is None and error <= bound:
self.evaluations[i] = progress.evaluations
def print(self, name, result):
"""A row of the table: the evaluations, "-" for an error not reached, and the best
error."""
cells = ["-" if reached is None else str(reached) for reached in self.evaluations]
# rounding can put a solution a few ulps below the minimum
error = max(result.best_fitness - self.minimum, 0.0)
mantissa, exponent = f"{error:.1e}".split("e")
cells.append(f"{mantissa}e{int(exponent)}")
print(f"{name:<10}" + "".join(f"{cell:>9}" for cell in cells))
problem = gx.problems.Powell(DIMENSIONS)
minimum = problem.optimum.value
print(f"Powell in {DIMENSIONS} dimensions, {BUDGET} evaluations at most")
print("Evaluations until the error is at most")
print(f"{'algorithm':<10}" + "".join(f"{column:>9}" for column in COLUMNS) + f"{'best':>9}")
for name, algorithm in (
("CMA-ES", gx.Cmaes(problem.genome, objective="minimize", seed=1)),
("sep-CMA-ES", gx.Cmaes(problem.genome, covariance="diagonal", objective="minimize", seed=1)),
("DE", gx.De(problem.genome, objective="minimize", seed=1)),
("PSO", gx.Pso(problem.genome, population_size=40, objective="minimize", seed=1)),
(
"GA",
gx.Ga(
problem.genome,
population_size=100,
select=gx.Tournament(3),
crossover=gx.SimulatedBinaryCrossover(15.0),
mutation=gx.PolynomialMutation(20.0, rate=1 / DIMENSIONS),
objective="minimize",
seed=1,
),
),
):
reached = Reached(minimum)
result = algorithm.run(
problem, target=minimum + 1e-8, evaluations=BUDGET, on_generation=reached.record
)
reached.print(name, result)
# with GENOXIDE_TRACE=<file>, a trace for the plot on the example's page, of a separate run in
# 4 dimensions, Powell's own
record_small()
What it prints, from a seeded run:
Powell in 24 dimensions, 240000 evaluations at most
Evaluations until the error is at most
algorithm 1 1e-2 1e-4 1e-6 1e-8 best
CMA-ES 1716 4004 7943 13546 19864 9.9e-9
sep-CMA-ES 1404 2665 11427 104260 - 1.7e-7
DE 9300 14800 44100 60900 79400 9.7e-9
PSO 6520 13520 59640 - - 4.7e-6
GA 10901 119908 - - - 3.3e-3