From d16e0379a34646f1fc6ebe8f29c96bf7dc2d6734 Mon Sep 17 00:00:00 2001 From: Stephen Waits Date: Tue, 5 May 2026 08:00:47 -0600 Subject: [PATCH] feat(algorithms): add MOPSO (Multi-Objective Particle Swarm) Coello, Pulido & Lechuga 2004 MOPSO: PSO adapted for multi-objective optimization via an external Pareto archive used as the source of swarm leaders. Each generation: - Evaluate every particle's current position - Insert non-dominated members into the archive (using ParetoArchive) - For each particle, pick a leader from the archive (uniform random among archive members) - Update velocity using inertia + cognitive (toward pbest) + social (toward leader) - Update positions, clamp to bounds - Refresh personal bests using Pareto comparison: pbest is replaced only when the new position dominates it; on non-dominated, keep with 50/50 random tiebreak Vec decisions only. Truncates the archive to `archive_size` via the existing simple-tail truncation. Tests: produces a non-empty front on Schaffer N.1, deterministic reruns, panic on single-objective. --- src/algorithms/mod.rs | 2 + src/algorithms/mopso.rs | 229 ++++++++++++++++++++++++++++++++++++++++ src/prelude.rs | 2 +- 3 files changed, 232 insertions(+), 1 deletion(-) create mode 100644 src/algorithms/mopso.rs diff --git a/src/algorithms/mod.rs b/src/algorithms/mod.rs index e79d4fe..517159e 100644 --- a/src/algorithms/mod.rs +++ b/src/algorithms/mod.rs @@ -5,6 +5,7 @@ pub mod differential_evolution; pub mod genetic_algorithm; pub mod hill_climber; pub mod moead; +pub mod mopso; pub mod nsga2; pub mod nsga3; pub mod paes; @@ -20,6 +21,7 @@ pub use differential_evolution::*; pub use genetic_algorithm::*; pub use hill_climber::*; pub use moead::*; +pub use mopso::*; pub use nsga2::*; pub use nsga3::*; pub use paes::*; diff --git a/src/algorithms/mopso.rs b/src/algorithms/mopso.rs new file mode 100644 index 0000000..8076216 --- /dev/null +++ b/src/algorithms/mopso.rs @@ -0,0 +1,229 @@ +//! `Mopso` — Coello, Pulido & Lechuga 2004 Multi-Objective Particle Swarm. + +use rand::Rng as _; +use rand::seq::IndexedRandom; + +use crate::algorithms::parallel_eval::evaluate_batch; +use crate::core::candidate::Candidate; +use crate::core::population::Population; +use crate::core::problem::Problem; +use crate::core::result::OptimizationResult; +use crate::core::rng::rng_from_seed; +use crate::operators::real::RealBounds; +use crate::pareto::archive::ParetoArchive; +use crate::pareto::dominance::{Dominance, pareto_compare}; +use crate::pareto::front::{best_candidate, pareto_front}; +use crate::traits::Optimizer; + +/// Configuration for [`Mopso`]. +#[derive(Debug, Clone)] +pub struct MopsoConfig { + /// Number of particles in the swarm. + pub swarm_size: usize, + /// Number of generations. + pub generations: usize, + /// External Pareto archive size cap (simple-tail truncation). + pub archive_size: usize, + /// Inertia weight `w`. + pub inertia: f64, + /// Cognitive coefficient `c_1` (toward personal best). + pub cognitive: f64, + /// Social coefficient `c_2` (toward archive leader). + pub social: f64, + /// Seed for the deterministic RNG. + pub seed: u64, +} + +impl Default for MopsoConfig { + fn default() -> Self { + Self { + swarm_size: 40, + generations: 200, + archive_size: 100, + inertia: 0.7, + cognitive: 1.5, + social: 1.5, + seed: 42, + } + } +} + +/// Multi-objective particle swarm with an external Pareto archive. +/// +/// `Vec` decisions only. Each particle maintains a personal best (the +/// last position that was Pareto-non-dominated by any later position). The +/// social leader is sampled uniformly from the external archive each step. +#[derive(Debug, Clone)] +pub struct Mopso { + /// Algorithm configuration. + pub config: MopsoConfig, + /// Per-variable bounds — used both to seed the swarm and to clamp positions. + pub bounds: RealBounds, +} + +impl Mopso { + /// Construct a `Mopso`. + pub fn new(config: MopsoConfig, bounds: RealBounds) -> Self { + Self { config, bounds } + } +} + +impl

Optimizer

for Mopso +where + P: Problem> + Sync, +{ + fn run(&mut self, problem: &P) -> OptimizationResult { + assert!(self.config.swarm_size >= 1, "Mopso swarm_size must be >= 1"); + assert!(self.config.archive_size >= 1, "Mopso archive_size must be >= 1"); + let objectives = problem.objectives(); + assert!( + objectives.is_multi_objective(), + "Mopso requires multi-objective problems (use ParticleSwarm for single-objective)", + ); + let dim = self.bounds.bounds.len(); + let n = self.config.swarm_size; + let mut rng = rng_from_seed(self.config.seed); + + let mut positions: Vec> = { + use crate::traits::Initializer as _; + self.bounds.initialize(n, &mut rng) + }; + let mut velocities: Vec> = (0..n) + .map(|_| { + self.bounds + .bounds + .iter() + .map(|&(lo, hi)| 0.1 * (hi - lo) * (rng.random::() * 2.0 - 1.0)) + .collect() + }) + .collect(); + let v_max: Vec = self.bounds.bounds.iter().map(|&(lo, hi)| hi - lo).collect(); + + let initial_pop = evaluate_batch(problem, positions.clone()); + let mut evaluations = initial_pop.len(); + + // Personal bests start at initial positions. + let mut pbest_decisions: Vec> = positions.clone(); + let mut pbest_evals: Vec = + initial_pop.iter().map(|c| c.evaluation.clone()).collect(); + + // External archive seeded with the non-dominated subset. + let mut archive = ParetoArchive::new(objectives.clone()); + for c in initial_pop { + archive.insert(c); + } + archive.truncate(self.config.archive_size); + + for _ in 0..self.config.generations { + // --- Phase 1: serial position/velocity updates (uses RNG) --- + for i in 0..n { + let leader = archive + .members() + .choose(&mut rng) + .map(|c| c.decision.clone()) + .unwrap_or_else(|| positions[i].clone()); + #[allow(clippy::needless_range_loop)] // body indexes velocities/positions/bounds. + for j in 0..dim { + let r1: f64 = rng.random(); + let r2: f64 = rng.random(); + let cognitive_term = + self.config.cognitive * r1 * (pbest_decisions[i][j] - positions[i][j]); + let social_term = self.config.social * r2 * (leader[j] - positions[i][j]); + let mut v = self.config.inertia * velocities[i][j] + + cognitive_term + + social_term; + if v > v_max[j] { + v = v_max[j]; + } else if v < -v_max[j] { + v = -v_max[j]; + } + velocities[i][j] = v; + let (lo, hi) = self.bounds.bounds[j]; + positions[i][j] = (positions[i][j] + v).clamp(lo, hi); + } + } + + // --- Phase 2: parallel-friendly batch evaluation --- + let evaluated = evaluate_batch(problem, positions.clone()); + evaluations += evaluated.len(); + + // --- Phase 3: serial pbest + archive updates --- + for (i, cand) in evaluated.iter().enumerate() { + let dominance = + pareto_compare(&cand.evaluation, &pbest_evals[i], &objectives); + let replace = match dominance { + Dominance::Dominates => true, + Dominance::DominatedBy => false, + Dominance::Equal | Dominance::NonDominated => rng.random_bool(0.5), + }; + if replace { + pbest_decisions[i] = cand.decision.clone(); + pbest_evals[i] = cand.evaluation.clone(); + } + } + for c in evaluated { + archive.insert(c); + } + archive.truncate(self.config.archive_size); + } + + let members = archive.into_vec(); + let front = pareto_front(&members, &objectives); + let best = best_candidate(&members, &objectives); + OptimizationResult::new( + Population::new(members), + front, + best, + evaluations, + self.config.generations, + ) + } +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::tests_support::{SchafferN1, Sphere1D}; + + fn make_optimizer(seed: u64) -> Mopso { + Mopso::new( + MopsoConfig { + swarm_size: 30, + generations: 30, + archive_size: 30, + inertia: 0.7, + cognitive: 1.5, + social: 1.5, + seed, + }, + RealBounds::new(vec![(-5.0, 5.0)]), + ) + } + + #[test] + fn produces_pareto_front() { + let mut opt = make_optimizer(1); + let r = opt.run(&SchafferN1); + assert!(!r.pareto_front.is_empty()); + } + + #[test] + fn deterministic_with_same_seed() { + let mut a = make_optimizer(99); + let mut b = make_optimizer(99); + let ra = a.run(&SchafferN1); + let rb = b.run(&SchafferN1); + let oa: Vec> = + ra.pareto_front.iter().map(|c| c.evaluation.objectives.clone()).collect(); + let ob: Vec> = + rb.pareto_front.iter().map(|c| c.evaluation.objectives.clone()).collect(); + assert_eq!(oa, ob); + } + + #[test] + #[should_panic(expected = "multi-objective")] + fn single_objective_panics() { + let mut opt = make_optimizer(0); + let _ = opt.run(&Sphere1D); + } +} diff --git a/src/prelude.rs b/src/prelude.rs index 1b61849..5061d7a 100644 --- a/src/prelude.rs +++ b/src/prelude.rs @@ -24,7 +24,7 @@ pub use crate::operators::{ pub use crate::algorithms::{ CmaEs, CmaEsConfig, DifferentialEvolution, DifferentialEvolutionConfig, GeneticAlgorithm, GeneticAlgorithmConfig, HillClimber, HillClimberConfig, Moead, - MoeadConfig, Nsga2, Nsga2Config, Nsga3, Nsga3Config, Paes, PaesConfig, ParticleSwarm, + MoeadConfig, Mopso, MopsoConfig, Nsga2, Nsga2Config, Nsga3, Nsga3Config, Paes, PaesConfig, ParticleSwarm, ParticleSwarmConfig, RandomSearch, RandomSearchConfig, SimulatedAnnealing, SimulatedAnnealingConfig, Spea2, Spea2Config, TabuSearch, TabuSearchConfig, };