CEC 2006 g18
The problem
The CEC 2006 special session on constrained optimization (Liang et al., 2006) collected 24 problems, g01 to g24, from the literature, with their best known solutions and rules for comparing algorithms. g18 is the eighteenth, the report's equations 37 and 38 (page 11). The report takes it from Himmelblau (1972, Applied Nonlinear Programming, McGraw-Hill). The source maximizes; the report, and genoxide, minimize the negated value.
There are 9 variables, x1 to x8 in [−10, 10] and x9 in [0, 20]. The problem is
minimize f(x) = −0.5 (x1 x4 − x2 x3 + x3 x9 − x5 x9 + x5 x8 − x6 x7)
subject to g1 = x3² + x4² − 1 ≤ 0 g8 = (x3 − x7)² + (x4 − x8)² − 1 ≤ 0
g2 = x9² − 1 ≤ 0 g9 = x7² + (x8 − x9)² − 1 ≤ 0
g3 = x5² + x6² − 1 ≤ 0 g10 = x2 x3 − x1 x4 ≤ 0
g4 = x1² + (x2 − x9)² − 1 ≤ 0 g11 = −x3 x9 ≤ 0
g5 = (x1 − x5)² + (x2 − x6)² − 1 ≤ 0 g12 = x5 x9 ≤ 0
g6 = (x1 − x7)² + (x2 − x8)² − 1 ≤ 0 g13 = x6 x7 − x5 x8 ≤ 0
g7 = (x3 − x5)² + (x4 − x6)² − 1 ≤ 0
The formulas have a geometric reading, derived here from the definition. Take six points in the plane: A = (x1, x2), B = (x3, x4), C = (0, x9), D = (x5, x6), E = (x7, x8) and the origin O. Then −f is the shoelace formula for the area of the hexagon A-B-C-D-E-O. g1 to g9 keep the distance between two of the points at most 1: O and B, O and C, O and D, A and C, A and D, A and E, B and D, B and E, C and E. These are the hexagon's nine diagonals, the pairs of vertices that aren't neighbors. g10 to g13 keep the triangles O-A-B, O-B-C, O-C-D and O-D-E turning the same way.
The best known value is f = −0.866025403784439, −√3/2 to its digits, at the report's x. It isn't proven optimal. g1, g3, g4, g6, g7 and g9 are active there.
What makes it hard
No constraint bounds the hexagon's six sides, and the shoelace formula doesn't require a hexagon whose sides don't cross. The best known solutions use this. Their active constraints are the sides of two equilateral triangles with sides of 1, A-C-E and B-D-O. The hexagon A-B-C-D-E-O jumps from one triangle to the other, its sides cross five times, and the formula adds up the two triangles' areas: 2 · √3/4 = √3/2. No figure whose points are all within 1 of each other has that much area: even a disk of diameter 1 has π/4 ≈ 0.785.
A proper hexagon competes with it. At f = −0.674981, 0.191 above f, lies a convex hexagon whose sides don't cross and whose 15 distances are all at most 1. It's a local optimum, far from the best known solutions: x9 is 1 there, and 0.6 at the report's x. Of 25 runs of CMA-ES without restarts, 9 converged to it and stayed there.
The best known solutions are also not unique. The two triangles can shift against each other: at the report's x*, A, B and C are 0.6 from D, E and O, and at the run's solution 0.83. The feasible region is small too: of 10 million random points in the box, none was feasible.
Representation
A Real genome of 9 genes, x1 to x9, within the report's bounds. genoxide's
problems::cec2006::G18 is the fitness: the value f(x) and the total constraint violation, the
sum of max(0, g(x)) over the 13 inequalities, 0 for a feasible solution.
genoxide compares fitnesses with Deb's feasibility rules (Deb, 2000, Computer Methods in Applied Mechanics and Engineering 186: 311-338): a feasible solution beats an infeasible one, two feasible ones compare by value, and two infeasible ones by violation. The rules need no penalty weights.
Algorithm
CMA-ES (Hansen and Ostermeier, 2001, Evolutionary Computation 9(2): 159-195) samples a population from a normal distribution, and adapts its mean, step size and covariance matrix. It uses genoxide's defaults: a population of 4 + ⌊3 ln 9⌋ = 10, a step size of 0.3 of each gene's range and a random start. With IPOP restarts (Auger and Hansen, 2005, IEEE CEC 2005: 1769-1776), a run that has converged starts again from a random point, with twice the population. A sample outside the bounds is drawn again, up to 100 times, and then clipped to them. Deb's rules rank the samples.
The run has the report's budget of 500,000 evaluations, and stops once its best solution is feasible with an error f(x) − f* of at most 1e-8, an absolute error. The report counts a run as successful with an error of at most 1e-4; the example asks for more.
Why IPOP restarts: they let a run that converged to the local optimum start again elsewhere, with a larger population each time. With seeds 1 to 25, IPOP-CMA-ES met the target on all 25 runs, after a median of 18,480 evaluations (12,050 to 192,960). Without restarts, CMA-ES met it on 15; of the other 10, 9 stayed at the local optimum and one at −0.5, 0.366 above f*. L-SHADE (Tanabe and Fukunaga, 2014, IEEE CEC 2014: 1658-1665), genoxide's differential evolution with a population that shrinks over the budget, also met the target on all 25, after a median of 86,540 evaluations and at most 99,405: slower, but more even. SHADE (Tanabe and Fukunaga, 2013, IEEE CEC 2013: 71-78) met it on 24 of 25, after a median of 74,700; the other run stayed at the local optimum. The example shows CMA-ES, whose run with seed 1 misses the trap and needs no restart.
Output
The first line names the run. The second gives what stopped it, after how many evaluations and
restarts, the error f(x) − f and whether the best solution is feasible: "< 1e-8" means the run met
its target. The third gives when the best solution was first feasible, and when its error first
met the report's criterion of success. The fourth compares f(x) with f, to 6 significant digits,
and the fifth gives the solution, to 4 digits. The last names the active constraints, those with
|g| ≤ 1e-6, and counts the others. In Python, run evaluates the problem in Rust, so both versions
print the same.
The page shows each variable on its range, and each of the 13 constraints' state. Its curve shows the error f − f* of the best feasible solution, and of the population's median, on a log scale. The best's curve begins at the first feasible solution, and the median's once half the population is feasible.
The project page plays this run back.
Good results
A good run is feasible and ends within 1e-4 of f, the report's success. A value below f is possible, since f* is only the best known, but the runs here end just above it.
The recorded run finds its first feasible solution after 460 evaluations. Its first CMA-ES run, of 10 samples, finds the two triangles without coming near the local optimum, where x9 is 1: from then on, its best solution's x9 stays below 0.86. It is within 1e-4 of f after 5,310 evaluations, stays about 5.2e-5 above f from about 7,000 to 11,000 evaluations, and meets its target after 20,250, without a restart.
The solution is x = (−0.08706, −0.1688, 0.8193, −0.5734, −0.08696, −0.9962, 0.8192, 0.2539, 0.8274), another best known solution than the report's, with the same six active constraints: g1, g3, g4, g6, g7 and g9. A and D, B and E, C and O are 0.83 apart, further than at the report's x*: the two triangles are shifted against each other along C-O.
Known optimum: −0.866025403784439, −√3/2 (best known)
Source: examples/cec2006_g18
Interactive run: tachsin.gr/projects/genoxide/examples/cec2006-g18
cargo run --release --example cec2006_g18
//! CEC 2006 g18: a quadratic in 9 variables under 13 nonlinear inequality constraints, the area
//! of a hexagon whose diagonals are at most 1, maximized and negated, from the CEC 2006 special
//! session on constrained optimization (Liang et al., 2006). The best known value is
//! −0.866025403784439, −√3/2.
//!
//! genoxide's `G18` gives the value of a solution and its constraint violation, which Deb's
//! feasibility rules compare: a feasible solution beats an infeasible one. CMA-ES with IPOP
//! restarts searches the 9 variables within the report's budget of 500,000 evaluations, and stops
//! once the error f(x) − f* is at most 1e-8. The example prints the best solution and its
//! constraints.
//!
//! With `GENOXIDE_TRACE=<file>`, it also writes a trace of its run for the plot on the example's
//! page, with `trace.rs`.
//!
//! ```text
//! cargo run --release --example cec2006_g18
//! ```
mod trace;
use genoxide::prelude::*;
use genoxide::problems::Problem;
use genoxide::problems::cec2006::G18;
// the CEC 2006 report's budget of evaluations per run
const BUDGET: u64 = 500_000;
// the run stops once its best is feasible with an error f(x) - f* at most this
const ERROR: f64 = 1e-8;
// the report counts a run as successful once its error is at most this
const SUCCESS: f64 = 1e-4;
// a constraint within this of its boundary is active
const ACTIVE: f64 = 1e-6;
fn main() -> Result<()> {
let problem = G18;
let optimum = problem.optimum().expect("known");
let f_star = optimum.value();
let cmaes = Cmaes::builder(problem.representation())
.restarts(cmaes::Restarts::Ipop)
.minimize()
.seed(1)
.build()?;
// with GENOXIDE_TRACE=<file>, a trace of the run for the plot on the example's page
let mut trace = trace::Trace::from_env();
// the evaluations when the best is first feasible, and when its error first meets the
// report's criterion of success
let (mut feasible, mut success) = (None, None);
// the restarts, each of which doubles the population
let (mut restarts, mut size) = (0, 0);
let outcome = Engine::new(cmaes, problem)
.stop_when(Stop::target(f_star + ERROR).or(Stop::evaluations(BUDGET)))
.on_generation(|snapshot| {
let progress = snapshot.progress();
let best = progress.best().filter(|best| best.is_feasible());
let error = best.and_then(Fitness::score).map(|value| value - f_star);
if error.is_some() {
feasible.get_or_insert(progress.evaluations());
}
if error.is_some_and(|error| error <= SUCCESS) {
success.get_or_insert(progress.evaluations());
}
let population = snapshot.population().len();
if size != 0 && population > size {
restarts += 1;
}
size = population;
trace.record(snapshot);
})
.run()?;
let best = outcome.best_fitness();
let value = best.score().expect("valid");
let x = outcome.best_genome();
println!("IPOP-CMA-ES with Deb's feasibility rules on g18, seed 1");
let (stop, error) = if outcome.stop_reason() == StopReason::Target {
("stopped by the target", format!("< {ERROR:.0e}"))
} else {
("stopped", format!("{:.1e}", value - f_star))
};
let evaluations = outcome.evaluations();
let feasibility = if best.is_feasible() {
"feasible"
} else {
"infeasible"
};
println!(
"{stop} after {evaluations} evaluations and {restarts} restarts: f(x) - f* {error}, \
{feasibility}"
);
println!(
"first feasible after {} evaluations, f(x) - f* <= 1e-4 after {}",
count(feasible),
count(success)
);
println!(
"f(x) {}, f* {} ({})",
significant(value, 6),
significant(f_star, 6),
if optimum.is_proven() {
"proven"
} else {
"best known"
}
);
let genes: Vec<String> = (1..)
.zip(&x[..])
.map(|(i, xi)| format!("x{i} {}", significant(*xi, 4)))
.collect();
println!("{}", genes.join(", "));
// g1 to g13: the active ones, and how many of the others are satisfied or violated
let constraints = problem.constraints(x);
let g = constraints.inequalities();
let active: Vec<String> = (1..)
.zip(g)
.filter(|(_, g)| g.abs() <= ACTIVE)
.map(|(i, _)| format!("g{i}"))
.collect();
let satisfied = g.iter().filter(|&&g| g < -ACTIVE).count();
let violated = g.iter().filter(|&&g| g > ACTIVE).count();
println!(
"active (|g| <= 1e-6): {}; {satisfied} satisfied, {violated} violated",
active.join(", ")
);
trace.write();
Ok(())
}
// the evaluations, or "never"
fn count(evaluations: Option<u64>) -> String {
evaluations.map_or("never".to_string(), |evaluations| evaluations.to_string())
}
// `digits` significant digits, e.g. 29.9953 or -30665.5 for 6
fn significant(value: f64, digits: i32) -> String {
let magnitude = value.abs().log10().floor() as i32;
let decimals = (digits - 1 - magnitude).max(0) as usize;
format!("{value:.decimals$}")
}
python examples/cec2006_g18/main.py
"""CEC 2006 g18: a quadratic in 9 variables under 13 nonlinear inequality constraints, the area of
a hexagon whose diagonals are at most 1, maximized and negated, from the CEC 2006 special session
on constrained optimization (Liang et al., 2006). The best known value is −0.866025403784439,
−√3/2.
genoxide's ``G18`` gives the value of a solution and its constraint violation, which Deb's
feasibility rules compare: a feasible solution beats an infeasible one. CMA-ES with IPOP restarts
searches the 9 variables within the report's budget of 500,000 evaluations, and stops once the
error f(x) − f* is at most 1e-8. The example prints the best solution and its constraints.
``run`` evaluates the problem in Rust.
With ``GENOXIDE_TRACE=<file>``, it also writes a trace of its run for the plot on the example's
page, with trace.py.
python examples/cec2006_g18/main.py
"""
import math
import genoxide as gx
import numpy as np
from trace import Trace
# the CEC 2006 report's budget of evaluations per run
BUDGET = 500_000
# the run stops once its best is feasible with an error f(x) - f* at most this
ERROR = 1e-8
# the report counts a run as successful once its error is at most this
SUCCESS = 1e-4
# a constraint within this of its boundary is active
ACTIVE = 1e-6
def significant(value, digits):
"""``digits`` significant digits, e.g. 29.9953 or -30665.5 for 6."""
magnitude = math.floor(math.log10(abs(value)))
return f"{value:.{max(digits - 1 - magnitude, 0)}f}"
def scientific(value, decimals):
"""Scientific notation as Rust writes it, e.g. 1.2e-5 for 1 decimal."""
mantissa, exponent = f"{value:.{decimals}e}".split("e")
return f"{mantissa}e{int(exponent)}"
def count(evaluations):
"""The evaluations, or "never"."""
return "never" if evaluations is None else str(evaluations)
problem = gx.problems.cec2006.G18()
optimum = problem.optimum
f_star = optimum.value
cmaes = gx.Cmaes(problem.genome, restarts="ipop", objective=problem.objective, seed=1)
# with GENOXIDE_TRACE=<file>, a trace of the run for the plot on the example's page
trace = Trace(problem)
# the evaluations when the best is first feasible, and when its error first meets the report's
# criterion of success; the restarts, each of which doubles the population
first = {"feasible": None, "success": None, "restarts": 0, "size": 0}
def on_generation(progress):
_, violations = problem.evaluate(progress.best_genome[np.newaxis])
if violations[0] == 0.0:
error = progress.best_fitness - f_star
if first["feasible"] is None:
first["feasible"] = progress.evaluations
if first["success"] is None and error <= SUCCESS:
first["success"] = progress.evaluations
size = len(progress.scores)
if first["size"] and size > first["size"]:
first["restarts"] += 1
first["size"] = size
trace.record(progress)
result = cmaes.run(
problem, target=f_star + ERROR, evaluations=BUDGET, on_generation=on_generation
)
value = result.best_fitness
print("IPOP-CMA-ES with Deb's feasibility rules on g18, seed 1")
if result.stop_reason == "target":
stop, error = "stopped by the target", f"< {scientific(ERROR, 0)}"
else:
stop, error = "stopped", scientific(value - f_star, 1)
feasibility = "feasible" if result.violation == 0.0 else "infeasible"
print(
f"{stop} after {result.evaluations} evaluations and {first['restarts']} restarts: "
f"f(x) - f* {error}, {feasibility}"
)
print(
f"first feasible after {count(first['feasible'])} evaluations, f(x) - f* <= 1e-4 after "
f"{count(first['success'])}"
)
proven = "proven" if optimum.proven else "best known"
print(f"f(x) {significant(value, 6)}, f* {significant(f_star, 6)} ({proven})")
x = result.best_genome.tolist()
print(", ".join(f"x{i} {significant(xi, 4)}" for i, xi in enumerate(x, 1)))
# g1 to g13: the active ones, and how many of the others are satisfied or violated
g = problem.constraints(result.best_genome).tolist()
active = [f"g{i}" for i, gi in enumerate(g, 1) if abs(gi) <= ACTIVE]
satisfied = sum(gi < -ACTIVE for gi in g)
violated = sum(gi > ACTIVE for gi in g)
print(f"active (|g| <= 1e-6): {', '.join(active)}; {satisfied} satisfied, {violated} violated")
trace.write()
What it prints, from a seeded run:
IPOP-CMA-ES with Deb's feasibility rules on g18, seed 1
stopped by the target after 20250 evaluations and 0 restarts: f(x) - f* < 1e-8, feasible
first feasible after 460 evaluations, f(x) - f* <= 1e-4 after 5310
f(x) -0.866025, f* -0.866025 (best known)
x1 -0.08706, x2 -0.1688, x3 0.8193, x4 -0.5734, x5 -0.08696, x6 -0.9962, x7 0.8192, x8 0.2539, x9 0.8274
active (|g| <= 1e-6): g1, g3, g4, g6, g7, g9; 7 satisfied, 0 violated