swarmkit 0.1.0

Composable particle swarm optimization with nested searches
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122

//! `cargo run --example elliptical` — gbest PSO finding points on an
//! implicit ellipse. Uses [`Floats<2>`] as the particle unit so the
//! example doesn't have to hand-write `Add`/`Sub`/`Mul<f64>` boilerplate.

#![allow(
    clippy::print_stdout,
    clippy::unwrap_used,
    reason = "stand-alone demo: stdout is the user-facing surface; .unwrap() on a finite iterator is harmless"
)]

use rand::rngs::ThreadRng;
use rand::{Rng, RngExt as _};
use swarmkit::{
    Boundary, Contextful, Evolution, FieldwiseClamp as _, FitCalc, Floats, GBestMover,
    IntoGBestSearcher as _, PSOCoeffs, ParticleInit, ParticleMover as _, Searcher as _,
};

type Vec2 = Floats<2>;

#[derive(Copy, Clone, Default, Debug)]
struct Bounds {
    min: Vec2,
    max: Vec2,
}

#[derive(Copy, Clone, Default, Debug)]
struct RectBoundary {
    bounds: Bounds,
}

impl Contextful for RectBoundary {
    type TContext = ();
}

impl Boundary for RectBoundary {
    type T = Vec2;
    fn handle(&self, pos: Vec2) -> Vec2 {
        pos.clamp(self.bounds.min, self.bounds.max)
    }
}

struct RandomInit {
    particle_count: usize,
    bounds: Bounds,
}

impl ParticleInit for RandomInit {
    type T = Vec2;

    fn init_pos<R: Rng>(&self, rng: &mut R) -> Vec<Vec2> {
        (0..self.particle_count)
            .map(|_| sample(&self.bounds, rng))
            .collect()
    }

    fn init_vel<R: Rng>(&self, rng: &mut R) -> Vec<Vec2> {
        (0..self.particle_count)
            .map(|_| sample(&self.bounds, rng) * 0.5)
            .collect()
    }
}

fn sample<R: Rng>(b: &Bounds, rng: &mut R) -> Vec2 {
    Floats([
        rng.random_range(b.min[0]..b.max[0]),
        rng.random_range(b.min[1]..b.max[1]),
    ])
}

struct EllipticalFitCalc;

impl Contextful for EllipticalFitCalc {
    type TContext = ();
}

impl FitCalc for EllipticalFitCalc {
    type T = Vec2;

    fn calculate_fit(&self, v: Vec2) -> f64 {
        // exp(-|f|) maps "on ellipse" to fitness 1, smoothly falling off
        // outside; avoids the 1/|f| singularity.
        std::f64::consts::E.powf(-0.001 * elliptical_eq(v).abs())
    }
}

fn elliptical_eq(v: Vec2) -> f64 {
    let (x, y) = (v[0], v[1]);
    (x + y) * (x + y) - x * y - 100.0
}

fn main() {
    let bounds = Bounds {
        min: Floats([-100.0, -100.0]),
        max: Floats([100.0, 100.0]),
    };
    let fit_calc = EllipticalFitCalc;
    let init = RandomInit {
        particle_count: 30,
        bounds,
    };
    let coeffs = PSOCoeffs::new(0.2, 2.0, 1.0);

    let mut group = init.init(&mut ThreadRng::default());
    let mut searcher = GBestMover::<Vec2>::new(coeffs)
        .bounded_by(RectBoundary { bounds })
        .into_gbest_searcher(fit_calc);

    let mut evolution: Evolution<Vec2> = Evolution::default();
    let best = searcher
        .iter(80, &mut group, Some(&mut evolution))
        .last()
        .unwrap();

    let residual = elliptical_eq(best.best_pos).abs();
    println!(
        "best position : ({:.4}, {:.4})",
        best.best_pos[0], best.best_pos[1]
    );
    println!("best fitness  : {:.6}", best.best_fit);
    println!("|f(x, y)|     : {residual:.4}  (0 = on ellipse)");
    println!("iterations    : {}", evolution.frames().len());
}