From 0f5f2903a05e93a8a3e2fd3560f3a33eb1583afc Mon Sep 17 00:00:00 2001 From: Rutherther Date: Sun, 2 Nov 2025 18:44:11 +0100 Subject: [PATCH] feat: add a lot of algorithm showcases --- .envrc | 2 +- .gitignore | 1 + codes/Cargo.toml | 2 +- codes/eoa_lib/src/evolution.rs | 6 +- codes/eoa_lib/src/perturbation/mod.rs | 25 +- codes/tsp_hw01/src/initializers.rs | 9 +- codes/tsp_hw01/src/main.rs | 628 ++++++++++++++++++++++-- codes/tsp_plotter/Cargo.toml | 10 + codes/tsp_plotter/src/main.rs | 669 ++++++++++++++++++++++++++ 9 files changed, 1302 insertions(+), 50 deletions(-) create mode 100644 codes/tsp_plotter/Cargo.toml create mode 100644 codes/tsp_plotter/src/main.rs diff --git a/.envrc b/.envrc index e5dd41d79c2977d29dba6dd66aae50e7063d8b69..d1e88736a8b61366ff553b0d7e31f7dfe03f3013 100644 --- a/.envrc +++ b/.envrc @@ -1,3 +1,3 @@ use guix SHELL=$(guix shell -m manifest.scm bash-minimal -- bash -c 'echo $GUIX_ENVIRONMENT') -export LD_LIBRARY_PATH=$SHELL/lib +#export LD_LIBRARY_PATH=$SHELL/lib diff --git a/.gitignore b/.gitignore index 411f97b9f66025c5502bde5f756090effa7c1ae4..a366c60ccf4aae548d7a3be7ea17a73f32910232 100644 --- a/.gitignore +++ b/.gitignore @@ -1,2 +1,3 @@ *.png +*.svg target diff --git a/codes/Cargo.toml b/codes/Cargo.toml index 28906b229fec4eca6cffd594d58df2a225d694ed..6e6cedb45c2e6cfe418bc75ea5e5df5707026b52 100644 --- a/codes/Cargo.toml +++ b/codes/Cargo.toml @@ -1,3 +1,3 @@ [workspace] resolver = "3" -members = ["eoa_lib", "tsp_hw01"] +members = ["eoa_lib", "tsp_hw01", "tsp_plotter"] diff --git a/codes/eoa_lib/src/evolution.rs b/codes/eoa_lib/src/evolution.rs index aefd767d96425dc34c7e7836b213dd95f632faa7..0cad35bdc4a8d23736e0ec29773de42b23187576 100644 --- a/codes/eoa_lib/src/evolution.rs +++ b/codes/eoa_lib/src/evolution.rs @@ -26,7 +26,8 @@ pub struct EvolutionResult { pub population: EvaluatedPopulation, pub stats: EvolutionStats, pub best_candidate: EvaluatedChromosome, - pub iterations: usize + pub iterations: usize, + pub evaluations: usize } pub fn evolution_algorithm @@ -163,7 +164,8 @@ pub fn evolution_algorithm population: current_population, best_candidate, stats, - iterations + iterations, + evaluations: current_evaluation }) } diff --git a/codes/eoa_lib/src/perturbation/mod.rs b/codes/eoa_lib/src/perturbation/mod.rs index d8d9dd9429293f6aeca137f64db96d01133019b0..6e6d0e9d14d7e89133a94cc058b1c689779d1f31 100644 --- a/codes/eoa_lib/src/perturbation/mod.rs +++ b/codes/eoa_lib/src/perturbation/mod.rs @@ -1,7 +1,7 @@ use std::{any::Any, borrow::{Borrow, BorrowMut}, marker::PhantomData}; use nalgebra::{allocator::Allocator, DefaultAllocator, Dim, SVector}; -use rand::{distr::Distribution, Rng, RngCore}; +use rand::{distr::Distribution, Rng, RngCore, prelude::IteratorRandom}; use rand_distr::{uniform, Normal, NormalError, Uniform}; use crate::binary_string::BinaryString; @@ -432,6 +432,29 @@ impl PerturbationOperator for CombinedPerturbation { } } +pub struct OneOfPerturbation { + perturbations: Vec>> +} + +impl OneOfPerturbation { + pub fn new(perturbations: Vec>>) -> Self { + Self { + perturbations + } + } +} + +impl PerturbationOperator for OneOfPerturbation { + type Chromosome = T; + + fn perturb(&self, chromosome: &mut Self::Chromosome, rng: &mut dyn RngCore) { + let chosen = (0..self.perturbations.len()).choose(rng); + if let Some(chosen) = chosen { + self.perturbations[chosen].perturb(chromosome, rng); + } + } +} + #[cfg(test)] pub mod tests { use crate::binary_string::BinaryString; diff --git a/codes/tsp_hw01/src/initializers.rs b/codes/tsp_hw01/src/initializers.rs index 2cc6ddc44786971f039e01bbcdd4c3d3566c86bc..b180133574cced13dadd6c70ba39f00f58c34a2b 100644 --- a/codes/tsp_hw01/src/initializers.rs +++ b/codes/tsp_hw01/src/initializers.rs @@ -125,7 +125,7 @@ where { fn initialize_single(&self, size: D, rng: &mut dyn RngCore) -> NodePermutation { let starting_node = rng.random_range(0..size.value()); - self.initialize_from(size, starting_node, rng.random_range(0.0..=0.7), rng) + self.initialize_from(size, starting_node, 0.0, rng) } fn initialize( @@ -170,7 +170,12 @@ where // 3. Initialize randomly with random p_choose_second (call initialize_single) for _ in 0..count { - population.push(self.initialize_single(size, rng)) + population.push(self.initialize_from( + size, + rng.random_range(0..size.value()), + rng.random_range(0.0..0.5), + rng + )); } population diff --git a/codes/tsp_hw01/src/main.rs b/codes/tsp_hw01/src/main.rs index dc68cc9113e3cd9dc795b470a2fe74d28da110b7..1c47320698fba49cdde35185ebfa44a3bae4ec6a 100644 --- a/codes/tsp_hw01/src/main.rs +++ b/codes/tsp_hw01/src/main.rs @@ -13,7 +13,7 @@ use perturbations::{MovePerturbation, Random2OptPerturbation, ReverseSubsequence use binary_string_representation::TSPBinaryStringWrapper; use nalgebra::{Dim, Dyn}; use eoa_lib::{ - binary_string::BinaryString, comparison::MinimizingOperator, crossover::BinaryNPointCrossover, evolution::{evolution_algorithm, EvolutionStats}, initializer::{Initializer, RandomInitializer}, local_search::{local_search_first_improving, LocalSearchStats}, pairing::AdjacentPairing, perturbation::{apply_to_perturbations, BinaryStringBitPerturbation, BinaryStringFlipNPerturbation, BinaryStringSingleBitPerturbation, CombinedPerturbation, MutationPerturbation, PerturbationOperator}, replacement::{BestReplacement, TournamentReplacement}, selection::{BestSelection, RouletteWheelSelection, TournamentSelection}, terminating::MaximumCyclesTerminatingCondition + binary_string::BinaryString, comparison::MinimizingOperator, crossover::BinaryNPointCrossover, evolution::{evolution_algorithm, EvolutionStats}, evolutionary_strategy::IdentityStrategy, initializer::{Initializer, RandomInitializer}, local_search::{local_search_first_improving, LocalSearchStats}, pairing::AdjacentPairing, perturbation::{apply_to_perturbations, BinaryStringBitPerturbation, BinaryStringFlipNPerturbation, BinaryStringSingleBitPerturbation, CombinedPerturbation, IdentityPerturbation, MutationPerturbation, OneOfPerturbation, PerturbationOperator}, random_search::random_search, replacement::{BestReplacement, TournamentReplacement}, selection::{BestSelection, RouletteWheelSelection, TournamentSelection}, terminating::MaximumCyclesTerminatingCondition }; use rand::rng; use std::env; @@ -22,6 +22,14 @@ use std::io::{BufRead, BufReader, Write}; use flate2::read::GzDecoder; use chrono::{DateTime, Local}; +// Algorithm iteration/cycle constants +const EA_MAX_ITERATIONS: usize = 5000; +const LS_MAX_CYCLES: usize = 250 * 5000 + 500; + +// EA population constants +const EA_POPULATION_SIZE: usize = 500; +const EA_PARENTS_COUNT: usize = 250; + fn load_tsp_instance(filename: &str) -> Result, Box> { let file = File::open(filename)?; let reader: Box = if filename.ends_with(".gz") { @@ -79,8 +87,6 @@ fn extract_evolution_data( final_solution: &NodePermutation, final_evaluation: f64, final_iteration: usize, - initial_population_size: usize, - offspring_count: usize, ) -> PlotData { let mut iterations = Vec::new(); let mut evaluations = Vec::new(); @@ -90,11 +96,6 @@ fn extract_evolution_data( evaluations.push(candidate.evaluated_chromosome.evaluation); } - // Add final result - let final_fitness_evaluations = initial_population_size + final_iteration * offspring_count; - iterations.push(final_fitness_evaluations); - evaluations.push(final_evaluation); - PlotData { best_solution: final_solution.clone(), iterations, @@ -110,8 +111,6 @@ fn extract_binary_evolution_data( final_solution: &NodePermutation, final_evaluation: f64, final_iteration: usize, - initial_population_size: usize, - offspring_count: usize, ) -> PlotData { let mut iterations = Vec::new(); let mut evaluations = Vec::new(); @@ -121,11 +120,6 @@ fn extract_binary_evolution_data( evaluations.push(candidate.evaluated_chromosome.evaluation); } - // Add final result - let final_fitness_evaluations = initial_population_size + final_iteration * offspring_count; - iterations.push(final_fitness_evaluations); - evaluations.push(final_evaluation); - PlotData { best_solution: final_solution.clone(), iterations, @@ -231,24 +225,151 @@ fn run_evolution_algorithm(instance: &TSPInstance) -> Result) -> Result> { + let mut rng = rng(); + let initializer = MinimumSpanningTreeInitializer::new(instance); + let dimension = instance.dimension(); + + // Create combined perturbation with two mutations wrapped in MutationPerturbation + let move_mutation = MutationPerturbation::new(Box::new(MovePerturbation::new()), 0.1); + let swap_mutation = MutationPerturbation::new(Box::new(SwapPerturbation::new()), 0.1); + let reverse_mutation = MutationPerturbation::new(Box::new(ReverseSubsequencePerturbation::new()), 0.1); + let mut combined_perturbation = CombinedPerturbation::new(vec![ + Box::new(move_mutation), + Box::new(swap_mutation), + Box::new(reverse_mutation), + ]); + + // Set up other components + let mut crossover = EdgeRecombinationCrossover::new(); + let mut selection = RouletteWheelSelection::new(); + let mut replacement = BestReplacement::new(); + let mut pairing = AdjacentPairing::new(); + let better_than_operator = MinimizingOperator::new(); + + // Create initial population + let population_size = EA_POPULATION_SIZE; + let initial_population = initializer.initialize(dimension, population_size, &mut rng); + + let initial_population = eoa_lib::replacement::Population::from_vec(initial_population); + + // Run evolution algorithm + let parents_count = EA_PARENTS_COUNT; + let result = evolution_algorithm( + initial_population.clone(), + parents_count, + instance, + &mut selection, + &mut pairing, + &mut crossover, + &mut combined_perturbation, + &mut replacement, + &better_than_operator, + EA_MAX_ITERATIONS, // max iterations + &mut rng, + |iteration, stats, _, _, _, _, perturbation, _| { + let iters_till_end = EA_MAX_ITERATIONS - iteration + 1; + let iters_since_better = + iteration - stats.best_candidates.last().map(|c| c.iteration).unwrap_or(0); + MutationPerturbation::apply_to_mutations( + perturbation, + &mut |p| { + p.probability = (0.5 * (1.0 + (iters_since_better as f64 / iters_till_end as f64))).min(1.0); + } + ); + } + )?; + + // Extract plotting data + let plot_data = extract_evolution_data( + &result.stats, + &result.best_candidate.chromosome, + result.best_candidate.evaluation, + result.iterations, + ); + + Ok(plot_data) +} + +fn run_evolution_algorithm_nn(instance: &TSPInstance) -> Result> { + let mut rng = rng(); + let initializer = NearestNeighborInitializer::new(instance); + let dimension = instance.dimension(); + + // Create combined perturbation with two mutations wrapped in MutationPerturbation + let move_mutation = MutationPerturbation::new(Box::new(MovePerturbation::new()), 0.1); + let swap_mutation = MutationPerturbation::new(Box::new(SwapPerturbation::new()), 0.1); + let reverse_mutation = MutationPerturbation::new(Box::new(ReverseSubsequencePerturbation::new()), 0.1); + let mut combined_perturbation = CombinedPerturbation::new(vec![ + Box::new(move_mutation), + Box::new(swap_mutation), + Box::new(reverse_mutation), + ]); + + // Set up other components + let mut crossover = EdgeRecombinationCrossover::new(); + let mut selection = RouletteWheelSelection::new(); + let mut replacement = BestReplacement::new(); + let mut pairing = AdjacentPairing::new(); + let better_than_operator = MinimizingOperator::new(); + + // Create initial population + let population_size = EA_POPULATION_SIZE; + let initial_population = initializer.initialize(dimension, population_size, &mut rng); let initial_population = eoa_lib::replacement::Population::from_vec(initial_population); // Run evolution algorithm - let parents_count = 250; + let parents_count = EA_PARENTS_COUNT; let result = evolution_algorithm( initial_population.clone(), parents_count, @@ -259,10 +380,10 @@ fn run_evolution_algorithm(instance: &TSPInstance) -> Result) -> Result) -> Result> { + let mut rng = rng(); + let initializer = TSPRandomInitializer::new(); + let dimension = instance.dimension(); + + // Create combined perturbation with lower probabilities and no adaptive changes + let move_mutation = MutationPerturbation::new(Box::new(MovePerturbation::new()), 0.001); + let swap_mutation = MutationPerturbation::new(Box::new(SwapPerturbation::new()), 0.001); + let reverse_mutation = MutationPerturbation::new(Box::new(ReverseSubsequencePerturbation::new()), 0.001); + let mut combined_perturbation = CombinedPerturbation::new(vec![ + Box::new(move_mutation), + Box::new(swap_mutation), + Box::new(reverse_mutation), + ]); + + // Set up other components with Cycle Crossover + let mut crossover = CycleCrossover::new(); + let mut selection = RouletteWheelSelection::new(); + let mut replacement = BestReplacement::new(); + let mut pairing = AdjacentPairing::new(); + let better_than_operator = MinimizingOperator::new(); + + // Create initial population + let population_size = EA_POPULATION_SIZE; + let initial_population = initializer.initialize(dimension, population_size, &mut rng); + + let initial_population = eoa_lib::replacement::Population::from_vec(initial_population); + + // Run evolution algorithm + let parents_count = EA_PARENTS_COUNT; + let result = evolution_algorithm( + initial_population.clone(), + parents_count, + instance, + &mut selection, + &mut pairing, + &mut crossover, + &mut combined_perturbation, + &mut replacement, + &better_than_operator, + EA_MAX_ITERATIONS, // max iterations + &mut rng, + |_iteration, _stats, _, _, _, _, _perturbation, _| { + // No adaptive perturbation changes - keep probabilities constant + } + )?; + + // Extract plotting data + let plot_data = extract_evolution_data( + &result.stats, + &result.best_candidate.chromosome, + result.best_candidate.evaluation, + result.iterations, + ); + + Ok(plot_data) +} + +fn run_evolution_algorithm_pmx(instance: &TSPInstance) -> Result> { + let mut rng = rng(); + let initializer = TSPRandomInitializer::new(); + let dimension = instance.dimension(); + + // Create combined perturbation with lower probabilities and no adaptive changes + let move_mutation = MutationPerturbation::new(Box::new(MovePerturbation::new()), 0.001); + let swap_mutation = MutationPerturbation::new(Box::new(SwapPerturbation::new()), 0.001); + let reverse_mutation = MutationPerturbation::new(Box::new(ReverseSubsequencePerturbation::new()), 0.001); + let mut combined_perturbation = CombinedPerturbation::new(vec![ + Box::new(move_mutation), + Box::new(swap_mutation), + Box::new(reverse_mutation), + ]); + + // Set up other components with Partially Mapped Crossover + let mut crossover = PartiallyMappedCrossover::new(); + let mut selection = RouletteWheelSelection::new(); + let mut replacement = BestReplacement::new(); + let mut pairing = AdjacentPairing::new(); + let better_than_operator = MinimizingOperator::new(); + + // Create initial population + let population_size = EA_POPULATION_SIZE; + let initial_population = initializer.initialize(dimension, population_size, &mut rng); + + let initial_population = eoa_lib::replacement::Population::from_vec(initial_population); + + // Run evolution algorithm + let parents_count = EA_PARENTS_COUNT; + let result = evolution_algorithm( + initial_population.clone(), + parents_count, + instance, + &mut selection, + &mut pairing, + &mut crossover, + &mut combined_perturbation, + &mut replacement, + &better_than_operator, + EA_MAX_ITERATIONS, // max iterations + &mut rng, + |_iteration, _stats, _, _, _, _, _perturbation, _| { + // No adaptive perturbation changes - keep probabilities constant + } + )?; + + // Extract plotting data + let plot_data = extract_evolution_data( + &result.stats, + &result.best_candidate.chromosome, + result.best_candidate.evaluation, + result.iterations, + ); + + Ok(plot_data) +} + +fn run_evolution_algorithm_erx(instance: &TSPInstance) -> Result> { + let mut rng = rng(); + let initializer = TSPRandomInitializer::new(); + let dimension = instance.dimension(); + + // Create combined perturbation with lower probabilities and no adaptive changes + let move_mutation = MutationPerturbation::new(Box::new(MovePerturbation::new()), 0.001); + let swap_mutation = MutationPerturbation::new(Box::new(SwapPerturbation::new()), 0.001); + let reverse_mutation = MutationPerturbation::new(Box::new(ReverseSubsequencePerturbation::new()), 0.001); + let mut combined_perturbation = CombinedPerturbation::new(vec![ + Box::new(move_mutation), + Box::new(swap_mutation), + Box::new(reverse_mutation), + ]); + + // Set up other components with Edge Recombination Crossover + let mut crossover = EdgeRecombinationCrossover::new(); + let mut selection = RouletteWheelSelection::new(); + let mut replacement = BestReplacement::new(); + let mut pairing = AdjacentPairing::new(); + let better_than_operator = MinimizingOperator::new(); + + // Create initial population + let population_size = EA_POPULATION_SIZE; + let initial_population = initializer.initialize(dimension, population_size, &mut rng); + + let initial_population = eoa_lib::replacement::Population::from_vec(initial_population); + + // Run evolution algorithm + let parents_count = EA_PARENTS_COUNT; + let result = evolution_algorithm( + initial_population.clone(), + parents_count, + instance, + &mut selection, + &mut pairing, + &mut crossover, + &mut combined_perturbation, + &mut replacement, + &better_than_operator, + EA_MAX_ITERATIONS, // max iterations + &mut rng, + |_iteration, _stats, _, _, _, _, _perturbation, _| { + // No adaptive perturbation changes - keep probabilities constant + } + )?; + + // Extract plotting data + let plot_data = extract_evolution_data( + &result.stats, + &result.best_candidate.chromosome, + result.best_candidate.evaluation, + result.iterations, ); Ok(plot_data) @@ -315,14 +606,14 @@ fn run_evolution_algorithm_binary(instance: &TSPInstance) -> Result) -> Result) -> Result) -> Result> { +fn run_local_search_reverse(instance: &TSPInstance) -> Result> { let mut rng = rng(); let initializer = TSPRandomInitializer::new(); let dimension = instance.dimension(); @@ -388,9 +675,42 @@ fn run_local_search(instance: &TSPInstance) -> Result) -> Result> { + let mut rng = rng(); + let initializer = MinimumSpanningTreeInitializer::new(instance); + let dimension = instance.dimension(); + + // Create a MST-based initial solution + let initial_solution = initializer.initialize_single(dimension, &mut rng); + // Run local search let mut perturbation = MovePerturbation::new(); - let mut terminating_condition = MaximumCyclesTerminatingCondition::new(250 * 5000 + 500); + let mut terminating_condition = MaximumCyclesTerminatingCondition::new(LS_MAX_CYCLES); let better_than_operator = MinimizingOperator::new(); let result = local_search_first_improving( @@ -413,13 +733,191 @@ fn run_local_search(instance: &TSPInstance) -> Result) -> Result> { + let mut rng = rng(); + let initializer = NearestNeighborInitializer::new(instance); + let dimension = instance.dimension(); + + // Create a nearest neighbor-based initial solution + let initial_solution = initializer.initialize_single(dimension, &mut rng); + + // Run local search + let mut perturbation = MovePerturbation::new(); + let mut terminating_condition = MaximumCyclesTerminatingCondition::new(LS_MAX_CYCLES); + let better_than_operator = MinimizingOperator::new(); + + let result = local_search_first_improving( + instance, + &mut terminating_condition, + &mut perturbation, + &better_than_operator, + &initial_solution, + &mut rng, + )?; + + // Extract plotting data + let plot_data = extract_local_search_data( + &result.stats, + &result.best_candidate.pos, + result.best_candidate.fit, + result.cycles, + ); + + Ok(plot_data) +} + +fn run_local_search_swap(instance: &TSPInstance) -> Result> { + let mut rng = rng(); + let initializer = TSPRandomInitializer::new(); + let dimension = instance.dimension(); + + // Create a random initial solution + let initial_solution = initializer.initialize_single(dimension, &mut rng); + + // Run local search + let mut perturbation = SwapPerturbation::new(); + let mut terminating_condition = MaximumCyclesTerminatingCondition::new(LS_MAX_CYCLES); + let better_than_operator = MinimizingOperator::new(); + + let result = local_search_first_improving( + instance, + &mut terminating_condition, + &mut perturbation, + &better_than_operator, + &initial_solution, + &mut rng, + )?; + + // Extract plotting data + let plot_data = extract_local_search_data( + &result.stats, + &result.best_candidate.pos, + result.best_candidate.fit, + result.cycles, + ); + + Ok(plot_data) +} + +fn run_local_search_move(instance: &TSPInstance) -> Result> { + let mut rng = rng(); + let initializer = TSPRandomInitializer::new(); + let dimension = instance.dimension(); + + // Create a random initial solution + let initial_solution = initializer.initialize_single(dimension, &mut rng); + + // Run local search + let mut perturbation = MovePerturbation::new(); + let mut terminating_condition = MaximumCyclesTerminatingCondition::new(LS_MAX_CYCLES); + let better_than_operator = MinimizingOperator::new(); + + let result = local_search_first_improving( + instance, + &mut terminating_condition, + &mut perturbation, + &better_than_operator, + &initial_solution, + &mut rng, + )?; + + // Extract plotting data + let plot_data = extract_local_search_data( + &result.stats, + &result.best_candidate.pos, + result.best_candidate.fit, + result.cycles, + ); + + Ok(plot_data) +} + +fn run_local_search_mix(instance: &TSPInstance) -> Result> { + let mut rng = rng(); + let initializer = TSPRandomInitializer::new(); + let dimension = instance.dimension(); + + // Create a random initial solution + let initial_solution = initializer.initialize_single(dimension, &mut rng); + + // Run local search with mixed perturbations + let mut perturbation = OneOfPerturbation::new(vec![ + Box::new(MovePerturbation::new()), + Box::new(SwapPerturbation::new()), + Box::new(ReverseSubsequencePerturbation::new()), + ]); + let mut terminating_condition = MaximumCyclesTerminatingCondition::new(LS_MAX_CYCLES); + let better_than_operator = MinimizingOperator::new(); + + let result = local_search_first_improving( + instance, + &mut terminating_condition, + &mut perturbation, + &better_than_operator, + &initial_solution, + &mut rng, + )?; + + // Extract plotting data + let plot_data = extract_local_search_data( + &result.stats, + &result.best_candidate.pos, + result.best_candidate.fit, + result.cycles, + ); + + Ok(plot_data) +} + +fn run_random_search(instance: &TSPInstance) -> Result> { + let mut rng = rng(); + let initializer = TSPRandomInitializer::new(); + let dimension = instance.dimension(); + + // Run random search + let mut terminating_condition = MaximumCyclesTerminatingCondition::new(LS_MAX_CYCLES); + let better_than_operator = MinimizingOperator::new(); + + let result = random_search::< + Dyn, + NodePermutation, + f64, + TSPInstance, + MaximumCyclesTerminatingCondition, + IdentityPerturbation>, + MinimizingOperator, + IdentityStrategy, + TSPRandomInitializer + >( + instance, + &mut terminating_condition, + &better_than_operator, + &initializer, + dimension, + &mut rng, + )?; + + // Extract plotting data + let mut plot_data = extract_local_search_data( + &result.stats, + &result.best_candidate.pos, + result.best_candidate.fit, + result.cycles, + ); + + // Update algorithm name for random search + plot_data.algorithm_name = "Random Search".to_string(); + + Ok(plot_data) +} + fn main() -> Result<(), Box> { let args: Vec = env::args().collect(); if args.len() != 3 { eprintln!("Usage: {} ", args[0]); eprintln!(" instance_name: e.g., kroA100, berlin52, eil51"); - eprintln!(" algorithm: ea (evolution algorithm) or ls (local search)"); + eprintln!(" algorithm: ea, ea_mst, ea_nn, ea_cx, ea_pmx, ea_erx, ls_reverse, ls_swap, ls_move, ls_mix, ls_mst, ls_nn, rs, or ea_binary"); std::process::exit(1); } @@ -449,16 +947,60 @@ fn main() -> Result<(), Box> { println!("Running Evolution Algorithm..."); run_evolution_algorithm(&instance)? }, + "ea_mst" => { + println!("Running Evolution Algorithm with MST initialization..."); + run_evolution_algorithm_mst(&instance)? + }, + "ea_nn" => { + println!("Running Evolution Algorithm with Nearest Neighbor initialization..."); + run_evolution_algorithm_nn(&instance)? + }, + "ea_cx" => { + println!("Running Evolution Algorithm with Cycle Crossover..."); + run_evolution_algorithm_cx(&instance)? + }, + "ea_pmx" => { + println!("Running Evolution Algorithm with Partially Mapped Crossover..."); + run_evolution_algorithm_pmx(&instance)? + }, + "ea_erx" => { + println!("Running Evolution Algorithm with Edge Recombination Crossover..."); + run_evolution_algorithm_erx(&instance)? + }, "ea_binary" => { println!("Running Evolution Algorithm (Binary)..."); run_evolution_algorithm_binary(&instance)? }, - "ls" => { - println!("Running Local Search..."); - run_local_search(&instance)? + "ls_reverse" => { + println!("Running Local Search with Reverse Subsequence perturbation..."); + run_local_search_reverse(&instance)? + }, + "ls_swap" => { + println!("Running Local Search with Swap perturbation..."); + run_local_search_swap(&instance)? + }, + "ls_move" => { + println!("Running Local Search with Move perturbation..."); + run_local_search_move(&instance)? + }, + "ls_mix" => { + println!("Running Local Search with mixed perturbations..."); + run_local_search_mix(&instance)? + }, + "ls_mst" => { + println!("Running Local Search with MST initialization..."); + run_local_search_mst(&instance)? + }, + "ls_nn" => { + println!("Running Local Search with Nearest Neighbor initialization..."); + run_local_search_nn(&instance)? + }, + "rs" => { + println!("Running Random Search..."); + run_random_search(&instance)? }, _ => { - eprintln!("Unknown algorithm: {}. Use 'ea', 'ea_binary', or 'ls'", algorithm); + eprintln!("Unknown algorithm: {}. Use 'ea', 'ea_mst', 'ea_nn', 'ea_cx', 'ea_pmx', 'ea_erx', 'ls_reverse', 'ls_swap', 'ls_move', 'ls_mix', 'ls_mst', 'ls_nn', 'rs', or 'ea_binary'", algorithm); std::process::exit(1); } }; diff --git a/codes/tsp_plotter/Cargo.toml b/codes/tsp_plotter/Cargo.toml new file mode 100644 index 0000000000000000000000000000000000000000..00f2a0887f4e192f9fabf211dafb64fce09b96bc --- /dev/null +++ b/codes/tsp_plotter/Cargo.toml @@ -0,0 +1,10 @@ +[package] +name = "tsp_plotter" +version = "0.1.0" +edition = "2024" + +[dependencies] +plotters = "0.3" +csv = "1.3" +serde = { version = "1.0", features = ["derive"] } +serde_json = "1.0" diff --git a/codes/tsp_plotter/src/main.rs b/codes/tsp_plotter/src/main.rs new file mode 100644 index 0000000000000000000000000000000000000000..1054ba4d742c3fcd3a7a0f37257daad8d4ed9c73 --- /dev/null +++ b/codes/tsp_plotter/src/main.rs @@ -0,0 +1,669 @@ +use csv::Reader; +use plotters::prelude::*; +use std::collections::HashMap; +use std::fs; +use std::path::Path; +use serde::{Deserialize, Serialize}; + +#[derive(Debug, Clone, Serialize, Deserialize)] +#[serde(rename_all = "snake_case")] +enum PlotType { + FitnessEvolution, + SuccessProbability, +} + +#[derive(Debug)] +struct DataPoint { + evaluations: u32, + fitness: f64, + percentage_deviation: f64, +} + +#[derive(Debug, Clone, Serialize, Deserialize)] +struct PlotConfig { + instances: Vec, + algorithms: Vec, + group_by_algorithm: bool, + base_path: String, + output_path: String, + optimal_solutions: HashMap, + targets: Vec, + plot_type: PlotType, + average_targets: bool, + algorithm_labels: Option>, +} + +impl Default for PlotConfig { + fn default() -> Self { + let mut optimal_solutions = HashMap::new(); + optimal_solutions.insert("eil51".to_string(), 426.0); + optimal_solutions.insert("kroA100".to_string(), 21282.0); + + Self { + instances: vec!["eil51".to_string()], + algorithms: vec!["ea".to_string(), "ls".to_string(), "rs".to_string()], + group_by_algorithm: true, + base_path: "../tsp_hw01/solutions".to_string(), + output_path: "comparison_eil51.svg".to_string(), + optimal_solutions, + targets: vec![1.0, 5.0, 10.0], + plot_type: PlotType::FitnessEvolution, + average_targets: false, + algorithm_labels: None, + } + } +} + +fn calculate_percentage_deviation(fitness: f64, optimal: f64) -> f64 { + ((fitness - optimal) / optimal) * 100.0 +} + +fn get_algorithm_label(algorithm: &str, config: &PlotConfig) -> String { + config + .algorithm_labels + .as_ref() + .and_then(|labels| labels.get(algorithm)) + .map(|label| label.clone()) + .unwrap_or_else(|| algorithm.to_uppercase()) +} + +fn create_step_function(data: Vec) -> Vec { + if data.is_empty() { + return data; + } + + let mut result = Vec::new(); + + for i in 0..data.len() { + let current_point = &data[i]; + + // Add the actual data point (the vertical part of the step) + result.push(DataPoint { + evaluations: current_point.evaluations, + fitness: current_point.fitness, + percentage_deviation: current_point.percentage_deviation, + }); + + // If this is not the last point, add a horizontal step to the next evaluation + if i + 1 < data.len() { + let next_evaluation = data[i + 1].evaluations; + // Add a point just before the next evaluation with current fitness value + // This creates the horizontal step line + result.push(DataPoint { + evaluations: next_evaluation, + fitness: current_point.fitness, + percentage_deviation: current_point.percentage_deviation, + }); + } + } + + result +} + +#[derive(Debug)] +struct ProbabilityPoint { + evaluations: u32, + probability: f64, +} + +fn calculate_success_probability( + algorithm_data: &HashMap>, + target_percentage: f64, +) -> Vec { + if algorithm_data.is_empty() { + return Vec::new(); + } + + // Collect all unique evaluation points across all runs + let mut all_evaluations = std::collections::BTreeSet::new(); + for (_, points) in algorithm_data { + for point in points { + all_evaluations.insert(point.evaluations); + } + } + + let mut probability_points = Vec::new(); + + for &evaluation in &all_evaluations { + let total_runs = algorithm_data.len(); + let mut successful_runs = 0; + + // For each run, check if it has achieved the target at this evaluation + for (_, points) in algorithm_data { + // Find the best performance achieved up to this evaluation + let best_percentage = points + .iter() + .filter(|p| p.evaluations <= evaluation) + .map(|p| p.percentage_deviation) + .fold(f64::INFINITY, f64::min); + + if best_percentage <= target_percentage { + successful_runs += 1; + } + } + + let probability = successful_runs as f64 / total_runs as f64; + probability_points.push(ProbabilityPoint { + evaluations: evaluation, + probability, + }); + } + + probability_points +} + +fn calculate_averaged_success_probability( + algorithm_data: &HashMap>, + targets: &[f64], +) -> Vec { + if algorithm_data.is_empty() || targets.is_empty() { + return Vec::new(); + } + + // Calculate probability for each target + let target_probabilities: Vec> = targets + .iter() + .map(|&target| calculate_success_probability(algorithm_data, target)) + .collect(); + + if target_probabilities.is_empty() { + return Vec::new(); + } + + // Collect all unique evaluation points across all targets + let mut all_evaluations = std::collections::BTreeSet::new(); + for prob_data in &target_probabilities { + for point in prob_data { + all_evaluations.insert(point.evaluations); + } + } + + let mut averaged_points = Vec::new(); + + for &evaluation in &all_evaluations { + let mut total_probability = 0.0; + let mut count = 0; + + // Average probabilities across all targets for this evaluation + for prob_data in &target_probabilities { + if let Some(point) = prob_data.iter().find(|p| p.evaluations == evaluation) { + total_probability += point.probability; + count += 1; + } + } + + if count > 0 { + let averaged_probability = total_probability / count as f64; + averaged_points.push(ProbabilityPoint { + evaluations: evaluation, + probability: averaged_probability, + }); + } + } + + averaged_points +} + +fn read_csv_file(file_path: &Path, optimal_solution: f64) -> Result, Box> { + let mut reader = Reader::from_path(file_path)?; + let mut data = Vec::new(); + + for result in reader.records() { + let record = result?; + if record.len() >= 2 { + let evaluations: u32 = record[0].parse()?; + let fitness: f64 = record[1].parse()?; + let percentage_deviation = calculate_percentage_deviation(fitness, optimal_solution); + data.push(DataPoint { evaluations, fitness, percentage_deviation }); + } + } + + // Add a point at evaluation 1 with the same fitness as the first point + if !data.is_empty() { + let first_fitness = data[0].fitness; + let first_percentage = data[0].percentage_deviation; + data.insert(0, DataPoint { + evaluations: 1, + fitness: first_fitness, + percentage_deviation: first_percentage + }); + } + + Ok(data) +} + +fn read_all_csv_files(directory: &Path, optimal_solution: f64) -> Result>, Box> { + let mut all_data = HashMap::new(); + + if directory.is_dir() { + for entry in fs::read_dir(directory)? { + let entry = entry?; + let path = entry.path(); + + if path.extension().and_then(|s| s.to_str()) == Some("csv") { + let file_name = path.file_stem() + .and_then(|s| s.to_str()) + .unwrap_or("unknown") + .to_string(); + + match read_csv_file(&path, optimal_solution) { + Ok(data) => { + all_data.insert(file_name, data); + } + Err(e) => { + eprintln!("Error reading {}: {}", path.display(), e); + } + } + } + } + } + + Ok(all_data) +} + +#[derive(Debug)] +struct PlotData { + data: HashMap>>, +} + +fn read_plot_data(config: &PlotConfig) -> Result> { + let mut data = HashMap::new(); + let base_path = Path::new(&config.base_path); + + for algorithm in &config.algorithms { + let mut algorithm_data = HashMap::new(); + + for instance in &config.instances { + let path = base_path.join(format!("{}/{}", algorithm, instance)); + + if path.exists() { + let optimal_solution = config.optimal_solutions.get(instance) + .copied() + .unwrap_or_else(|| { + eprintln!("Warning: No optimal solution found for instance '{}', using default value 1000.0", instance); + 1000.0 + }); + + let instance_data = read_all_csv_files(&path, optimal_solution)?; + + // Apply step function to create proper step visualization + let mut step_data = HashMap::new(); + for (run_key, points) in instance_data { + let step_points = create_step_function(points); + step_data.insert(run_key, step_points); + } + + algorithm_data.extend(step_data); + } + } + + if !algorithm_data.is_empty() { + data.insert(algorithm.clone(), algorithm_data); + } + } + + Ok(PlotData { data }) +} + +fn get_color_palette() -> Vec { + vec![ + BLUE, RED, GREEN, CYAN, MAGENTA, + RGBColor(255, 165, 0), // Orange + RGBColor(128, 0, 128), // Purple + RGBColor(255, 192, 203), // Pink + RGBColor(165, 42, 42), // Brown + RGBColor(128, 128, 128), // Gray + ] +} + +fn create_plot(plot_data: &PlotData, config: &PlotConfig) -> Result<(), Box> { + // Create plots directory if it doesn't exist + fs::create_dir_all("plots")?; + + let output_path = Path::new("plots").join(&config.output_path); + + let root = SVGBackend::new(&output_path, (1024, 768)).into_drawing_area(); + root.fill(&WHITE)?; + + let mut min_evaluations = u32::MAX; + let mut max_evaluations = 0u32; + + let mut min_percentage = f64::MAX; + let mut max_percentage = f64::MIN; + + // Find min/max across all data + for (_, algorithm_data) in &plot_data.data { + for (_, points) in algorithm_data { + for point in points { + min_evaluations = min_evaluations.min(point.evaluations); + max_evaluations = max_evaluations.max(point.evaluations); + min_percentage = min_percentage.min(point.percentage_deviation); + max_percentage = max_percentage.max(point.percentage_deviation); + } + } + } + + println!("Percentage deviation range: {:.2}% to {:.2}%", min_percentage, max_percentage); + + // Use actual data range with some padding for better visualization + let padding_factor = 0.1; + let range = max_percentage - min_percentage; + let y_min = (min_percentage - range * padding_factor).max(0.1); // Ensure y_min > 0 for log scale + let y_max = max_percentage + range * padding_factor; + + let mut chart = ChartBuilder::on(&root) + .margin(10) + .x_label_area_size(50) + .y_label_area_size(90) + .build_cartesian_2d( + (min_evaluations as f64..max_evaluations as f64).log_scale(), + (y_min..y_max).log_scale(), + )?; + + chart + .configure_mesh() + .x_desc("Evaluations") + .y_desc("Percentage deviation from optimal (%)") + .y_labels(10) + .x_label_formatter(&|x| { + let power = x.log10().round() as i32; + match power { + 0 => "10⁰".to_string(), + 1 => "10¹".to_string(), + 2 => "10²".to_string(), + 3 => "10³".to_string(), + 4 => "10⁴".to_string(), + 5 => "10⁵".to_string(), + 6 => "10⁶".to_string(), + 7 => "10⁷".to_string(), + 8 => "10⁸".to_string(), + 9 => "10⁹".to_string(), + _ => format!("10^{}", power), + } + }) + .y_label_formatter(&|y| { + if *y >= 100.0 { + format!("{:.0}%", *y) + } else if *y >= 10.0 { + format!("{:.1}%", *y) + } else { + format!("{:.2}%", *y) + } + }) + .x_max_light_lines(8) + .y_max_light_lines(12) + .axis_desc_style(("sans-serif", 18)) + .label_style(("sans-serif", 16)) + .draw()?; + + let colors = get_color_palette(); + let mut color_index = 0; + let mut legend_added = HashMap::new(); + + if config.group_by_algorithm { + // Group by algorithm: EA, LS, RS each get one color + for algorithm in &config.algorithms { + if let Some(algorithm_data) = plot_data.data.get(algorithm) { + let color = colors[color_index % colors.len()]; + color_index += 1; + + for (_, points) in algorithm_data { + let series = chart + .draw_series(LineSeries::new( + points.iter().map(|p| (p.evaluations as f64, p.percentage_deviation)), + &color, + ))?; + + if !legend_added.contains_key(algorithm) { + let label = get_algorithm_label(algorithm, config); + series.label(label).legend(move |(x, y)| PathElement::new(vec![(x, y), (x + 10, y)], &color)); + legend_added.insert(algorithm.clone(), true); + } + } + } + } + } else { + // Group by instance: Each algorithm-instance combination gets different color + for algorithm in &config.algorithms { + if let Some(algorithm_data) = plot_data.data.get(algorithm) { + for instance in &config.instances { + for (run_key, points) in algorithm_data { + if run_key.contains(instance) || config.instances.len() == 1 { + let color = colors[color_index % colors.len()]; + color_index += 1; + + let series = chart + .draw_series(LineSeries::new( + points.iter().map(|p| (p.evaluations as f64, p.percentage_deviation)), + &color, + ))?; + + let algo_label = get_algorithm_label(algorithm, config); + let label = format!("{} {}", algo_label, instance); + if !legend_added.contains_key(&label) { + series.label(&label).legend(move |(x, y)| PathElement::new(vec![(x, y), (x + 10, y)], &color)); + legend_added.insert(label, true); + } + } + } + } + } + } + } + + // Draw horizontal target lines + for &target in &config.targets { + if target >= y_min && target <= y_max { + let series = chart + .draw_series(LineSeries::new( + vec![(min_evaluations as f64, target), (max_evaluations as f64, target)], + BLACK.stroke_width(2), + ))?; + + series + .label(&format!("{}% target", target)) + .legend(move |(x, y)| PathElement::new(vec![(x, y), (x + 10, y)], BLACK)); + } + } + + chart.configure_series_labels() + .background_style(&WHITE) + .border_style(&BLACK) + .position(SeriesLabelPosition::LowerLeft) + .draw()?; + root.present()?; + + println!("Plot saved to: {}", output_path.display()); + Ok(()) +} + +fn create_success_probability_plot(plot_data: &PlotData, config: &PlotConfig) -> Result<(), Box> { + // Create plots directory if it doesn't exist + fs::create_dir_all("plots")?; + + let output_path = Path::new("plots").join(&config.output_path); + + let root = SVGBackend::new(&output_path, (1024, 768)).into_drawing_area(); + root.fill(&WHITE)?; + + let mut min_evaluations = u32::MAX; + let mut max_evaluations = 0u32; + + // Find min/max evaluations across all data + for (_, algorithm_data) in &plot_data.data { + for (_, points) in algorithm_data { + for point in points { + min_evaluations = min_evaluations.min(point.evaluations); + max_evaluations = max_evaluations.max(point.evaluations); + } + } + } + + let mut chart = ChartBuilder::on(&root) + .margin(10) + .x_label_area_size(50) + .y_label_area_size(90) + .build_cartesian_2d( + (min_evaluations as f64..max_evaluations as f64).log_scale(), + 0.0f64..1.0f64, + )?; + + chart + .configure_mesh() + .x_desc("Evaluations") + .y_desc("Probability of Success") + .y_labels(11) + .x_label_formatter(&|x| { + let power = x.log10().round() as i32; + match power { + 0 => "10⁰".to_string(), + 1 => "10¹".to_string(), + 2 => "10²".to_string(), + 3 => "10³".to_string(), + 4 => "10⁴".to_string(), + 5 => "10⁵".to_string(), + 6 => "10⁶".to_string(), + 7 => "10⁷".to_string(), + 8 => "10⁸".to_string(), + 9 => "10⁹".to_string(), + _ => format!("10^{}", power), + } + }) + .y_label_formatter(&|y| { + format!("{:.1}", *y) + }) + .x_max_light_lines(8) + .y_max_light_lines(0) + .axis_desc_style(("sans-serif", 18)) + .label_style(("sans-serif", 16)) + .draw()?; + + let colors = get_color_palette(); + let mut color_index = 0; + + // Plot probability curves + if config.average_targets { + // Plot one averaged curve per algorithm + for algorithm in &config.algorithms { + if let Some(algorithm_data) = plot_data.data.get(algorithm) { + let probability_data = calculate_averaged_success_probability(algorithm_data, &config.targets); + + if !probability_data.is_empty() { + let color = colors[color_index % colors.len()]; + color_index += 1; + + let series = chart + .draw_series(LineSeries::new( + probability_data.iter().map(|p| (p.evaluations as f64, p.probability)), + &color, + ))?; + + let algo_label = get_algorithm_label(algorithm, config); + let label = format!("{} (avg)", algo_label); + series + .label(&label) + .legend(move |(x, y)| PathElement::new(vec![(x, y), (x + 10, y)], &color)); + } + } + } + } else { + // Plot individual curves for each algorithm and target combination + for algorithm in &config.algorithms { + if let Some(algorithm_data) = plot_data.data.get(algorithm) { + for &target in &config.targets { + let probability_data = calculate_success_probability(algorithm_data, target); + + if !probability_data.is_empty() { + let color = colors[color_index % colors.len()]; + color_index += 1; + + let series = chart + .draw_series(LineSeries::new( + probability_data.iter().map(|p| (p.evaluations as f64, p.probability)), + &color, + ))?; + + let algo_label = get_algorithm_label(algorithm, config); + let label = format!("{} - {}%", algo_label, target); + series + .label(&label) + .legend(move |(x, y)| PathElement::new(vec![(x, y), (x + 10, y)], &color)); + } + } + } + } + } + + chart.configure_series_labels().background_style(&WHITE).border_style(&BLACK).draw()?; + root.present()?; + + println!("Probability plot saved to: {}", output_path.display()); + Ok(()) +} + +fn load_config(config_path: &str) -> Result> { + + if Path::new(config_path).exists() { + let config_str = fs::read_to_string(config_path)?; + let config: PlotConfig = serde_json::from_str(&config_str)?; + println!("Loaded configuration from {}", config_path); + Ok(config) + } else { + let config = PlotConfig::default(); + + // Create default config file + let config_str = serde_json::to_string_pretty(&config)?; + fs::write(config_path, config_str)?; + println!("Created default configuration file: {}", config_path); + + Ok(config) + } +} + +fn main() -> Result<(), Box> { + let args: Vec = std::env::args().collect(); + + if args.len() != 2 { + eprintln!("Usage: {} ", args[0]); + eprintln!("Example: {} plot_config.json", args[0]); + std::process::exit(1); + } + + let config_path = &args[1]; + let config = load_config(config_path)?; + + println!("Configuration:"); + println!(" Instances: {:?}", config.instances); + println!(" Algorithms: {:?}", config.algorithms); + println!(" Group by algorithm: {}", config.group_by_algorithm); + println!(" Base path: {}", config.base_path); + println!(" Output path: {}", config.output_path); + println!(" Plot type: {:?}", config.plot_type); + println!(" Average targets: {}", config.average_targets); + println!(); + + let base_directory = Path::new(&config.base_path); + println!("Reading CSV files from: {}", base_directory.display()); + + let plot_data = read_plot_data(&config)?; + + if plot_data.data.is_empty() { + eprintln!("No CSV files found for any configured algorithms and instances"); + return Ok(()); + } + + for (algorithm, algorithm_data) in &plot_data.data { + println!("Found {} CSV files for algorithm: {}", algorithm_data.len(), algorithm); + } + + match config.plot_type { + PlotType::FitnessEvolution => { + create_plot(&plot_data, &config)?; + } + PlotType::SuccessProbability => { + create_success_probability_plot(&plot_data, &config)?; + } + } + + Ok(()) +}