Skip to content

Sum of different powers

The problem

The sum of different powers raises each gene's absolute value to a power one more than its index, to minimize:

f(x) = Σ |xᵢ|^(i+1),   i from 1 to n, each xᵢ in [−1, 1]

Its minimum is 0, at the origin. Here n = 30, so the powers run from 2 to 31. Its origin is unknown: genoxide takes the definition and the bounds from Molga and Smutnicki (2005, section 2.8), and they're still to be checked against an original (issue #168).

What makes it hard

It's unimodal and separable, and each term is smallest at 0. But the terms differ widely in how much they matter. Near the minimum, the first gene's term is a parabola, while the thirtieth's, |x|³¹, is flat: at x₃₀ = 0.5 it's 5·10⁻¹⁰, and an error of 1e-8 allows x₃₀ up to 0.55. A search reaches small values long before the later genes are near 0, and the flatter terms give it little to follow.

Shifted and rotated, every direction mixes steep and flat terms, and the shapes that the steps must learn are no longer along the axes.

Representation

A Real genome of 30 genes, each in [−1, 1]: the point x itself. The fitness is f(x), to minimize. The function is genoxide's problems::SumOfDifferentPowers, which brings its bounds and its minimum.

Algorithm

Five algorithms, each with a budget of 10,000 evaluations per dimension, 300,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 14 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/30 per gene.

The second table runs the same algorithms on the function shifted and rotated, with genoxide's problems::Shifted and problems::Rotated and seed 1: the minimum moves to a random point in the middle 80% of the box, and an orthogonal matrix, drawn from normal numbers made orthonormal by Gram-Schmidt as BBOB draws its rotations, turns the function about it. That's how the CEC and BBOB suites use the function, with their own data; genoxide generates its instances instead.

Output

The first line gives the dimension and the budget. Then two tables, the function as it is and shifted and rotated: a row per algorithm, the evaluations it had used when its best error first reached each value of the heading, and the best error it found, to two significant digits. A dash is an error not reached. The function is evaluated with genoxide's portable math, so the runs are the same on every platform, and in Python, run evaluates it in Rust, so both versions print the same.

The project page plays back another run: CMA-ES on the function in 2 dimensions, |x₁|² + |x₂|³, so that the population can be drawn on its contour. It meets the target after 204 evaluations.

Good results

The minimum is 0. On the function as it is, sep-CMA-ES reaches 1e-8 first, after 2,464 evaluations, then PSO (4,880), SHADE (6,500), CMA-ES (13,006) and the genetic algorithm (15,922): every one gets there. The function is separable, and the methods that work a gene at a time, or learn one scale per gene, are the fastest.

Shifted and rotated, CMA-ES takes about as long as before, 12,334 evaluations: its full covariance matrix learns the rotation. SHADE needs five times as many, 32,200, and PSO 243,400. sep-CMA-ES reaches 1e-6 after 10,220 evaluations but ends at 1.4e-8, and the genetic algorithm at 5.9e-8.

Reference: Molga, M. and Smutnicki, C. (2005). Test functions for optimization needs.

Known optimum: 0 (at the origin)

Source: examples/sum_of_different_powers

Interactive run: tachsin.gr/projects/genoxide/examples/sum-of-different-powers

cargo run --release --example sum_of_different_powers
//! Sum of different powers: minimize the sum of the genes' absolute values to powers from 2 to 31, in 30
//! dimensions: the later the gene, the flatter the function near 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::SumOfDifferentPowers`. Then the same on the function shifted and rotated, with genoxide's `problems::Shifted`
//! and `problems::Rotated`, as the CEC and BBOB suites transform it.
//!
//! 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 sum_of_different_powers
//! ```

mod trace;

use genoxide::observer::Snapshot;
use genoxide::prelude::*;
use genoxide::problems::{Problem, Rotated, Shifted, SumOfDifferentPowers};

const DIMENSIONS: usize = 30;
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<()> {
    println!("Sum of different powers in {DIMENSIONS} dimensions, {BUDGET} evaluations at most");
    compare(
        "Sum of different powers",
        &SumOfDifferentPowers::new(DIMENSIONS),
    )?;
    // the same function, shifted and rotated, as the CEC and BBOB suites transform theirs
    let rotated = Rotated::new(Shifted::new(SumOfDifferentPowers::new(DIMENSIONS), 1), 1);
    compare("Shifted and rotated (seed 1)", &rotated)?;
    // with GENOXIDE_TRACE=<file>, a trace for the plot on the example's page, of a separate
    // run in 2 dimensions: the plot is the function's contour
    trace::record_small()?;
    Ok(())
}

// the table of the five algorithms on `problem`, after a line that names it
fn compare<P>(name: &str, problem: &P) -> Result<()>
where
    P: Problem<Representation = Real> + FitnessFunction<Reals, Output = f64> + Clone,
{
    let minimum = problem.optimum().expect("known").value();
    let stop = || Stop::target(minimum + 1e-8).or(Stop::evaluations(BUDGET));
    println!("{name}: 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.clone())
            .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.clone())
        .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.clone())
        .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.clone())
        .stop_when(stop())
        .on_generation(|snapshot| reached.record(snapshot))
        .run()?;
    reached.print("GA", &outcome);
    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/sum_of_different_powers/main.py
"""Sum of different powers: minimize the sum of the genes' absolute values to powers from 2 to 31,
in 30 dimensions: the later the gene, the flatter the function near 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::SumOfDifferentPowers`. Then the same on the
function shifted and rotated, with genoxide's `problems::Shifted` and `problems::Rotated`, as the
CEC and BBOB suites transform it.

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/sum_of_different_powers/main.py
"""

import genoxide as gx

from trace import record_small

DIMENSIONS = 30
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"]


def error_text(error):
    """An error to two significant digits, as Rust writes it: 9.9e-9."""
    mantissa, exponent = f"{error:.1e}".split("e")
    return f"{mantissa}e{int(exponent)}"


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
        cells.append(error_text(max(result.best_fitness - self.minimum, 0.0)))
        print(f"{name:<10}" + "".join(f"{cell:>9}" for cell in cells))


def compare(name, problem):
    """The table of the five algorithms on ``problem``, after a line that names it."""
    minimum = problem.optimum.value
    print(f"{name}: evaluations until the error is at most")
    print(f"{'algorithm':<10}" + "".join(f"{column:>9}" for column in COLUMNS) + f"{'best':>9}")
    genome = problem.genome
    for label, algorithm in (
        ("CMA-ES", gx.Cmaes(genome, objective="minimize", seed=1)),
        ("sep-CMA-ES", gx.Cmaes(genome, covariance="diagonal", objective="minimize", seed=1)),
        ("DE", gx.De(genome, objective="minimize", seed=1)),
        ("PSO", gx.Pso(genome, population_size=40, objective="minimize", seed=1)),
        (
            "GA",
            gx.Ga(
                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(label, result)


print(f"Sum of different powers in {DIMENSIONS} dimensions, {BUDGET} evaluations at most")
compare("Sum of different powers", gx.problems.SumOfDifferentPowers(DIMENSIONS))
# the same function, shifted and rotated, as the CEC and BBOB suites transform theirs
rotated = gx.problems.Rotated(
    gx.problems.Shifted(gx.problems.SumOfDifferentPowers(DIMENSIONS), seed=1), seed=1
)
compare("Shifted and rotated (seed 1)", rotated)

# with GENOXIDE_TRACE=<file>, a trace for the plot on the example's page, of a separate run in
# 2 dimensions: the plot is the function's contour
record_small()

What it prints, from a seeded run:

Sum of different powers in 30 dimensions, 300000 evaluations at most
Sum of different powers: evaluations until the error is at most
algorithm         1     1e-2     1e-4     1e-6     1e-8     best
CMA-ES           28      294     1624     5306    13006   9.1e-9
sep-CMA-ES       28      252      742     1540     2464   9.2e-9
DE              100     1100     3100     4800     6500  7.4e-10
PSO              40      760     1680     3120     4880   7.2e-9
GA              100     1517     3781     6915    15922   8.7e-9
Shifted and rotated (seed 1): evaluations until the error is at most
algorithm         1     1e-2     1e-4     1e-6     1e-8     best
CMA-ES          168      476     1666     4718    12334   9.0e-9
sep-CMA-ES      182      322      784    10220        -   1.4e-8
DE              600     2100     5100    16200    32200   7.5e-9
PSO             280      920     2560    14280   243400   9.9e-9
GA              293     1609     3316    27087        -   5.9e-8