Skip to main content

graphmind_optimization/algorithms/
rao.rs

1use crate::common::{Individual, OptimizationResult, Problem, SolverConfig};
2use ndarray::Array1;
3use rand::prelude::*;
4use rayon::prelude::*;
5
6#[derive(Debug, Clone, Copy)]
7pub enum RaoVariant {
8    Rao1,
9    Rao2,
10    Rao3,
11}
12
13pub struct RaoSolver {
14    pub config: SolverConfig,
15    pub variant: RaoVariant,
16}
17
18impl RaoSolver {
19    pub fn new(config: SolverConfig, variant: RaoVariant) -> Self {
20        Self { config, variant }
21    }
22
23    pub fn solve<P: Problem>(&self, problem: &P) -> OptimizationResult {
24        let mut rng = thread_rng();
25        let dim = problem.dim();
26        let (lower, upper) = problem.bounds();
27
28        // Initialize population
29        let mut population: Vec<Individual> = (0..self.config.population_size)
30            .map(|_| {
31                let mut vars = Array1::zeros(dim);
32                for i in 0..dim {
33                    vars[i] = rng.gen_range(lower[i]..upper[i]);
34                }
35                let fitness = problem.fitness(&vars);
36                Individual::new(vars, fitness)
37            })
38            .collect();
39
40        let mut history = Vec::with_capacity(self.config.max_iterations);
41
42        for iter in 0..self.config.max_iterations {
43            if iter % 10 == 0 {
44                println!(
45                    "Rao Solver: Iteration {}/{}",
46                    iter, self.config.max_iterations
47                );
48            }
49            let (best_idx, worst_idx) = self.find_best_worst(&population);
50            let best_vars = population[best_idx].variables.clone();
51            let worst_vars = population[worst_idx].variables.clone();
52            let best_fitness = population[best_idx].fitness;
53
54            history.push(best_fitness);
55
56            // Update population
57            population = population
58                .into_par_iter()
59                .map(|mut ind| {
60                    let mut local_rng = thread_rng();
61                    let mut new_vars = Array1::zeros(dim);
62
63                    let r1: f64 = local_rng.gen();
64                    let r2: f64 = local_rng.gen();
65
66                    // For Rao-2 and Rao-3, we need a random individual.
67                    // In parallel map, we approximate by generating a random point in bounds.
68                    // This is "Best-Worst-Random" style.
69                    let mut rand_vars = Array1::zeros(dim);
70                    if matches!(self.variant, RaoVariant::Rao2 | RaoVariant::Rao3) {
71                        for j in 0..dim {
72                            rand_vars[j] = local_rng.gen_range(lower[j]..upper[j]);
73                        }
74                    }
75                    // We assume fitness of random point is worse? Or compare?
76                    // Standard Rao compares fitness. We'll compute it if needed.
77                    let rand_fitness =
78                        if matches!(self.variant, RaoVariant::Rao2 | RaoVariant::Rao3) {
79                            problem.fitness(&rand_vars)
80                        } else {
81                            0.0
82                        };
83
84                    for j in 0..dim {
85                        let term1 = best_vars[j] - worst_vars[j];
86
87                        let delta = match self.variant {
88                            RaoVariant::Rao1 => r1 * term1,
89                            RaoVariant::Rao2 => {
90                                let term2 = if ind.fitness < rand_fitness {
91                                    ind.variables[j] - rand_vars[j]
92                                } else {
93                                    rand_vars[j] - ind.variables[j]
94                                };
95                                r1 * term1 + r2 * term2
96                            }
97                            RaoVariant::Rao3 => {
98                                // Rao-3 often uses |worst| and |best|
99                                let term1_abs = best_vars[j] - worst_vars[j].abs();
100                                let term2_abs = if ind.fitness < rand_fitness {
101                                    ind.variables[j] - rand_vars[j]
102                                } else {
103                                    rand_vars[j] - ind.variables[j]
104                                };
105                                r1 * term1_abs + r2 * term2_abs
106                            }
107                        };
108
109                        new_vars[j] = (ind.variables[j] + delta).clamp(lower[j], upper[j]);
110                    }
111
112                    let new_fitness = problem.fitness(&new_vars);
113                    if new_fitness < ind.fitness {
114                        ind.variables = new_vars;
115                        ind.fitness = new_fitness;
116                    }
117                    ind
118                })
119                .collect();
120        }
121
122        let (final_best_idx, _) = self.find_best_worst(&population);
123        let final_best = &population[final_best_idx];
124
125        OptimizationResult {
126            best_variables: final_best.variables.clone(),
127            best_fitness: final_best.fitness,
128            history,
129        }
130    }
131
132    fn find_best_worst(&self, population: &[Individual]) -> (usize, usize) {
133        let mut best_idx = 0;
134        let mut worst_idx = 0;
135        for (i, ind) in population.iter().enumerate() {
136            if ind.fitness < population[best_idx].fitness {
137                best_idx = i;
138            }
139            if ind.fitness > population[worst_idx].fitness {
140                worst_idx = i;
141            }
142        }
143        (best_idx, worst_idx)
144    }
145}