CEC 2006 g21
The problem
The CEC 2006 special session on constrained optimization (Liang et al., 2006) collected 24 test problems, g01 to g24, with their best known solutions and rules for comparing algorithms. g21 is the twenty-first, the report's equations 43 and 44 (page 13). The report takes it from Epperly's collection of global optimization test problems with solutions (the report's reference 6). It has 7 variables:
x1 in [0, 1000] x2, x3 in [0, 40] x4 in [100, 300] x5 in [6.3, 6.7]
x6 in [5.9, 6.4] x7 in [4.5, 6.25]
Minimize f(x) = x1, subject to one inequality, g(x) ≤ 0, and five equalities, each h(x) = 0:
g1 = −x1 + 35 x2^0.6 + 35 x3^0.6
h1 = −300 x3 + 7500 x5 − 7500 x6 − 25 x4 x5 + 25 x4 x6 + x3 x4
h2 = 100 x2 + 155.365 x4 + 2500 x7 − x2 x4 − 25 x4 x7 − 15536.5
h3 = −x5 + ln(−x4 + 900)
h4 = −x6 + ln(x4 + 300)
h5 = −x7 + ln(−2 x4 + 700)
The report counts an equality as met when |h(x)| ≤ 0.0001, and so does genoxide's G21. The best
known value is f = 193.724510070035, at the report's x:
x = (193.724510070035, 5.6e−27, 17.3191887294085, 100.047897801387, 6.68445185362378,
5.99168428444265, 6.21451648886070)
where g1 is active and every equality is at the edge of the tolerance, |h| = 0.0001. The value is the best known, not proven optimal.
One free variable
The equalities leave less freedom than they seem to. h3, h4 and h5 give x5, x6 and x7 from x4. h1 factors as (x4 − 300) (x3 − 25 (x5 − x6)) = 0, so x3 = 25 (x5 − x6) for x4 below 300. h2 factors as (100 − x4) (x2 − 155.365 + 25 x7) = 0, so x2 = 155.365 − 25 x7 for x4 above 100, and any x2 at x4 = 100. And f = x1 is least with g1 active, x1 = 35 x2^0.6 + 35 x3^0.6. With the equalities met exactly, f is a function of x4 alone.
That function has a minimum at each end. At x4 = 100, x2 = 0 and x3 = 25 ln 2 meet every constraint exactly, with f = 35 (25 ln 2)^0.6 = 193.788. It then rises, to 330.6 near x4 = 293, and falls again to 325.1 at x4 = 299.53, where x2 reaches its upper bound 40: past it, h2 would need x2 above 40. The tolerance of 0.0001 lets the best known solution go 0.064 below 193.788.
What makes it hard
The feasible region has no volume: random points never meet five equalities. The report's table 3 gives a feasible share of 0.0000 %. A search starts outside, guided only by the violation, and must reach a thin layer around a curve in 7 dimensions, the curve that x4 traces.
Then the curve has two ends, and a search that reaches the layer near x4 = 300 finds the other minimum, f = 324.70, 131 above f*, with x2 at 40 and x3 near 0. Moving along the curve to the good end means passing the maximum near x4 = 293: from there, every step toward it is worse, and Deb's rules never take a worse feasible step. With Deb's rules alone, where a search first reaches the layer decides where it ends.
And the best known solution has x2 = 5.6·10⁻²⁷. f grows like 35 x2^0.6, whose slope is infinite at 0: an error below 1e-8 needs x2 below about 10⁻¹⁶.
Representation
A Real genome of 7 genes, x1 to x7, within the bounds above. genoxide's problems::cec2006::G21
is the fitness: the value f(x) and the total constraint violation, the sum of max(0, g1(x)) and of
max(0, |h(x)| − 0.0001) over the equalities, 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. Before the first feasible solution, the search is a minimization of the violation, which leads to the region.
Algorithm
SHADE (Tanabe and Fukunaga, 2013, IEEE CEC 2013: 71-78) is genoxide's default differential evolution: current-to-pbest/1 mutation with an archive, and a memory of the scale factor F and the crossover rate CR that worked. It uses genoxide's defaults, but a population of 30 instead of 100, with restarts when the population has converged or stagnated, and a trial outside the bounds brought back halfway between its parent and the bound, which lets x2 approach 0 geometrically.
SHADE compares solutions with Deb's rules at an ε level, the ε constrained method of Takahama and
Sakai (2006, "Constrained optimization by the ε constrained differential evolution with
gradient-based mutation and feasible elites", IEEE CEC 2006), whose εDE won the CEC 2006
competition: a violation up to ε counts as none. Two solutions within ε compare by value, and the
others as Deb's rules have it. ε starts at 1,000 and follows Takahama and Sakai's schedule,
ε(t) = 1000 (1 − t / 150,000)⁵ after t evaluations, down to 0 at 150,000. The example applies it
with genoxide's features: the fitness function returns the violation beyond ε, and the engine's
control lowers ε every 10 generations and scores the population again (reevaluate), since the
old violations no longer compare with the new ones. The re-evaluations count in the budget. From
150,000 evaluations on, the rules are Deb's, and the problem the report's. While ε is large, the
value leads the population along the whole curve, toward the good end, before the layer is thin.
The run has the report's budget of 500,000 evaluations, and stops once ε is 0 and its best solution is feasible with an absolute error f(x) − f* of at most 1e-8. The report counts a run as successful with an error of at most 1e-4; the example asks for more.
With Deb's rules alone, no algorithm of genoxide solves g21 every time: each run reaches one end of the curve or the other. With 100 seeds, SHADE with a population of 50 met the target on 65, after a median of 45,350 evaluations (from 37,500 to 196,150); the other 35 ended at the minimum of 324.70. With the default population of 100, it met the target on 66, after a median of 93,500; L-SHADE, whose population shrinks over the budget, on 61. The restarts don't change this: with the default population, runs met the target on the same 16 of 25 seeds with and without them. CMA-ES with IPOP restarts (Hansen and Ostermeier, 2001; Auger and Hansen, 2005) came within 3·10⁻⁷ of f* on 15 of 25 seeds, but within 1e-8 on none: it converges pressed against the bound x2 = 0 and the edges of five tolerances at once.
At an ε level, SHADE met the target on 96 to 98 % of the runs with populations of 30 to 50 and ε starting at 100 or 1,000, and on 92 % with a population of 20, in 1,000 runs each. It needs its restarts: without them, no run met the target.
Output
The first line names the run. The second gives what stopped it, after how many evaluations, 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, by the ε level, was first feasible without it, and when its error first met the report's criterion of success. The fourth compares f(x) with f, to 6 significant digits. The fifth gives the solution, with the genes near 0 in scientific notation, and the last the constraints: for g1, "active" on the boundary (|g| ≤ 1e-6), else the value of g, negative when it's satisfied; for h1 to h5, "active" when the equality is met within the tolerance, else by how much |h| exceeds it. In Python, the fitness function evaluates the problem in Rust, a generation at a time, so both versions print the same.
The page's plot shows each variable on its range, and each constraint's state: violated, active or satisfied (an equality met within the tolerance shows as active). Its curve shows the error f − f* of the best solution, and of the population's median, on a log scale, measured without the ε level. The best's curve begins once the best is feasible, and the median's once half the population is: after the end of the ε schedule.
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. The runs that fail end at the other minimum, f = 324.70, 131 above f, or on their way along the curve.
The recorded run's best, by the ε level, uses the tolerance ε allows until the end of the schedule. It is first feasible without ε after 166,259 evaluations, within the report's criterion after 187,289, and meets its target after 193,229. The solution is the report's x* to 6 digits, with x2 = 7.8·10⁻¹⁸, g1 active and every equality met.
With seeds 1 to 1,000, 973 runs meet the target, after 175,919 to 486,899 evaluations (183,809 for half of them), and one more ends within 1e-4 of f*. Of the other 26, 17 end at the other minimum, 8 between 204.24 and 298.26, still moving along the curve when the budget runs out, and one at 322.77, near the other minimum but still infeasible.
Known optimum: 193.724510070035 (best known, with the equalities met within 0.0001)
Source: examples/cec2006_g21
Interactive run: tachsin.gr/projects/genoxide/examples/cec2006-g21
cargo run --release --example cec2006_g21
//! CEC 2006 g21: the linear function x1 of 7 variables under 1 nonlinear inequality and 5
//! nonlinear equality constraints, from the CEC 2006 special session on constrained optimization
//! (Liang et al., 2006). The best known value is 193.724510070035, with the equalities met within
//! the report's tolerance of 0.0001.
//!
//! genoxide's `G21` gives the value of a solution and its constraint violation. SHADE, genoxide's
//! default differential evolution, with a population of 30, compares them with Deb's feasibility
//! rules at an ε level (Takahama and Sakai, 2006): a violation up to ε counts as none. ε starts at
//! 1,000 and falls to 0 over the first 150,000 evaluations, and the population is scored again
//! each time it falls. The run has the report's budget of 500,000 evaluations, and stops once ε is
//! 0 and 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_g21
//! ```
mod trace;
use genoxide::prelude::*;
use genoxide::problems::Problem;
use genoxide::problems::cec2006::{EQUALITY_TOLERANCE, G21};
use std::sync::Arc;
use std::sync::atomic::{AtomicU64, Ordering};
// 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;
// the ε level at the start: a violation up to it counts as none
const EPSILON: f64 = 1000.0;
// the evaluations after which ε is 0
const CONTROL: u64 = 150_000;
// ε falls, and the population is scored again, every this many generations
const EVERY: u64 = 10;
// a constraint within this of its boundary is active
const ACTIVE: f64 = 1e-6;
// the ε level after `evaluations`: EPSILON (1 - evaluations / CONTROL)^5, then 0
fn epsilon(evaluations: u64) -> f64 {
if evaluations >= CONTROL {
return 0.0;
}
let rest = 1.0 - evaluations as f64 / CONTROL as f64;
EPSILON * rest * rest * rest * rest * rest
}
// the ε level in use, shared by the fitness function, the control and the stop condition
#[derive(Clone)]
struct Level(Arc<AtomicU64>);
impl Level {
fn get(&self) -> f64 {
f64::from_bits(self.0.load(Ordering::Relaxed))
}
fn set(&self, epsilon: f64) {
self.0.store(epsilon.to_bits(), Ordering::Relaxed);
}
}
fn main() -> Result<()> {
// the report's equality tolerance, EQUALITY_TOLERANCE
let problem = G21::default();
let optimum = problem.optimum().expect("known");
let f_star = optimum.value();
// a third of the default population: see the README
let shade = De::builder(problem.representation())
.population_size(30)
.minimize()
.seed(1)
.build()?;
let level = Level(Arc::new(AtomicU64::new(EPSILON.to_bits())));
// the value, and the violation beyond ε
let fitness = |x: &Reals| {
let (value, violation) = problem.evaluate(x);
(value, (violation - level.get()).max(0.0))
};
// the run's target, once ε is 0: a feasible best within ERROR of f*
let target = {
let level = level.clone();
Stop::custom(move |progress| {
level.get() == 0.0
&& progress.best().is_some_and(|best| {
best.is_feasible() && best.score().is_some_and(|value| value <= f_star + ERROR)
})
})
};
// 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);
let outcome = Engine::new(shade, fitness)
.stop_when(target.or(Stop::evaluations(BUDGET)))
.on_generation(|snapshot| {
let progress = snapshot.progress();
// the best by the ε level, measured without it
let (value, violation) = problem.evaluate(snapshot.best().genome());
if violation == 0.0 {
feasible.get_or_insert(progress.evaluations());
if value - f_star <= SUCCESS {
success.get_or_insert(progress.evaluations());
}
}
trace.record(snapshot);
})
.control(|shade, progress| {
// every EVERY generations, and once it reaches 0, ε follows its schedule
let next = epsilon(progress.evaluations());
let due = progress.generation() % EVERY == 0 || next == 0.0;
if progress.generation() > 0 && due && next != level.get() {
level.set(next);
shade.reevaluate()?;
}
Ok(())
})
.run()?;
let best = outcome.best_fitness();
let value = best.score().expect("valid");
let x = outcome.best_genome();
println!("SHADE, a population of 30, with Deb's rules at an epsilon level on g21, seed 1");
let (stop, error) = if outcome.stop_reason() == StopReason::Custom {
("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: 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} {}", gene(xi)))
.collect();
println!("{}", genes.join(", "));
// g1, then h1 to h5, each equality as its excess over the tolerance, 0 when it's met
let constraints = problem.constraints(x);
let inequalities = (1..)
.zip(constraints.inequalities())
.map(|(i, &g)| format!("g{i} {}", state(g)));
let equalities = (1..)
.zip(constraints.equalities())
.map(|(i, &h)| format!("h{i} {}", state((h.abs() - EQUALITY_TOLERANCE).max(0.0))));
let constraints: Vec<String> = inequalities.chain(equalities).collect();
println!("{}", constraints.join(", "));
trace.write();
Ok(())
}
// a constraint g(x) <= 0: "active" on its boundary, else its value
fn state(g: f64) -> String {
if g.abs() <= ACTIVE {
"active".to_string()
} else {
significant(g, 4)
}
}
// a gene: in scientific notation when it's that close to 0, e.g. 3.0e-12, else to 6 significant
// digits
fn gene(value: f64) -> String {
if value != 0.0 && value.abs() < 1e-4 {
format!("{value:.1e}")
} else {
significant(value, 6)
}
}
// 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_g21/main.py
"""CEC 2006 g21: the linear function x1 of 7 variables under 1 nonlinear inequality and 5 nonlinear
equality constraints, from the CEC 2006 special session on constrained optimization (Liang et al.,
2006). The best known value is 193.724510070035, with the equalities met within the report's
tolerance of 0.0001.
genoxide's ``G21`` gives the value of a solution and its constraint violation. SHADE, genoxide's
default differential evolution, with a population of 30, compares them with Deb's feasibility
rules at an ε level (Takahama and Sakai, 2006): a violation up to ε counts as none. ε starts at
1,000 and falls to 0 over the first 150,000 evaluations, and the population is scored again each
time it falls. The run has the report's budget of 500,000 evaluations, and stops once ε is 0 and
the error f(x) − f* is at most 1e-8. The example prints the best solution and its constraints.
The fitness function evaluates the problem in Rust, a generation at a time.
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_g21/main.py
"""
import math
import genoxide as gx
import numpy as np
from genoxide.problems.cec2006 import EQUALITY_TOLERANCE
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
# the ε level at the start: a violation up to it counts as none
EPSILON = 1000.0
# the evaluations after which ε is 0
CONTROL = 150_000
# ε falls, and the population is scored again, every this many generations
EVERY = 10
def epsilon(evaluations):
"""The ε level after ``evaluations``: EPSILON (1 - evaluations / CONTROL)^5, then 0."""
if evaluations >= CONTROL:
return 0.0
rest = 1.0 - evaluations / CONTROL
return EPSILON * rest * rest * rest * rest * rest
# 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 state(g):
"""A constraint g(x) <= 0: "active" on its boundary, else its value."""
return "active" if abs(g) <= ACTIVE else significant(g, 4)
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 gene(value):
"""A gene: in scientific notation when it's that close to 0, e.g. 3.0e-12, else to 6
significant digits."""
return scientific(value, 1) if value != 0 and abs(value) < 1e-4 else significant(value, 6)
def count(evaluations):
"""The evaluations, or "never"."""
return "never" if evaluations is None else str(evaluations)
# the report's equality tolerance, EQUALITY_TOLERANCE
problem = gx.problems.cec2006.G21()
optimum = problem.optimum
f_star = optimum.value
# a third of the default population: see the README
shade = gx.De(problem.genome, objective=problem.objective, population_size=30, seed=1)
# the ε level in use
level = {"epsilon": EPSILON}
def fitness(genomes):
"""The values of a generation, and their violations beyond ε."""
values, violations = problem.evaluate(genomes)
return values, np.maximum(violations - level["epsilon"], 0.0)
# 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
first = {"feasible": None, "success": None}
def on_generation(progress):
# the best by the ε level, measured without it
values, violations = problem.evaluate(progress.best_genome[np.newaxis])
feasible = violations[0] == 0.0
if feasible:
error = values[0] - f_star
if first["feasible"] is None:
first["feasible"] = progress.evaluations
if first["success"] is None and error <= SUCCESS:
first["success"] = progress.evaluations
trace.record(progress)
# the run's target, once ε is 0: a feasible best within ERROR of f*
return not (level["epsilon"] == 0.0 and feasible and values[0] <= f_star + ERROR)
def control(shade, progress):
# every EVERY generations, and once it reaches 0, ε follows its schedule
following = epsilon(progress.evaluations)
due = progress.generation % EVERY == 0 or following == 0.0
if progress.generation > 0 and due and following != level["epsilon"]:
level["epsilon"] = following
shade.reevaluate()
result = shade.run(
fitness,
batch=True,
evaluations=BUDGET,
on_generation=on_generation,
control=control,
)
value = result.best_fitness
print("SHADE, a population of 30, with Deb's rules at an epsilon level on g21, seed 1")
if result.stop_reason == "aborted":
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: 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} {gene(xi)}" for i, xi in enumerate(x, 1)))
# g1, then h1 to h5, each equality as its excess over the tolerance, 0 when it's met
g1, *h = problem.constraints(result.best_genome).tolist()
constraints = [f"g1 {state(g1)}"]
constraints += [
f"h{i} {state(max(abs(hj) - EQUALITY_TOLERANCE, 0.0))}" for i, hj in enumerate(h, 1)
]
print(", ".join(constraints))
trace.write()
What it prints, from a seeded run:
SHADE, a population of 30, with Deb's rules at an epsilon level on g21, seed 1
stopped by the target after 193229 evaluations: f(x) - f* < 1e-8, feasible
first feasible after 166259 evaluations, f(x) - f* <= 1e-4 after 187289
f(x) 193.725, f* 193.725 (best known)
x1 193.725, x2 7.8e-18, x3 17.3192, x4 100.048, x5 6.68445, x6 5.99168, x7 6.21452
g1 active, h1 active, h2 active, h3 active, h4 active, h5 active