math_audio_optimisation/
mutant_current_to_pbest1.rs

1//! Current-to-pbest/1 mutation strategy
2//!
3//! This is the key mutation strategy used by SHADE, L-SHADE and their variants.
4//! It blends the current individual with one of the top-p% best individuals,
5//! plus a difference vector that can include archive solutions.
6
7use ndarray::{Array1, Array2, Zip};
8use rand::Rng;
9
10use crate::distinct_indices::distinct_indices;
11use crate::external_archive::ExternalArchive;
12
13/// Current-to-pbest/1 mutation with optional archive
14///
15/// v_i = x_i + F * (x_pbest - x_i) + F * (x_r1 - x_r2)
16///
17/// Where:
18/// - x_pbest is randomly selected from top p% of population
19/// - x_r1 is a random individual from population
20/// - x_r2 can be from population OR from external archive
21pub(crate) fn mutant_current_to_pbest1<R: Rng + ?Sized>(
22    i: usize,
23    pop: &Array2<f64>,
24    sorted_indices: &[usize],
25    p_best_size: usize,
26    archive: Option<&ExternalArchive>,
27    f: f64,
28    rng: &mut R,
29) -> Array1<f64> {
30    let npop = pop.nrows();
31
32    let pbest_idx = if p_best_size >= npop {
33        sorted_indices[0]
34    } else {
35        sorted_indices[rng.random_range(0..p_best_size)]
36    };
37
38    let r1 = distinct_indices(i, 1, npop, rng)[0];
39
40    let r2 = if let Some(arch) = archive {
41        if !arch.is_empty() && rng.random::<f64>() < 0.5 {
42            if let Some(sol) = arch.random_select(rng) {
43                sol.clone()
44            } else {
45                let available: Vec<usize> =
46                    (0..npop).filter(|&idx| idx != i && idx != r1).collect();
47                if available.is_empty() {
48                    pop.row(r1).to_owned()
49                } else {
50                    pop.row(available[rng.random_range(0..available.len())])
51                        .to_owned()
52                }
53            }
54        } else {
55            let available: Vec<usize> = (0..npop).filter(|&idx| idx != i && idx != r1).collect();
56            if available.is_empty() {
57                pop.row(r1).to_owned()
58            } else {
59                pop.row(available[rng.random_range(0..available.len())])
60                    .to_owned()
61            }
62        }
63    } else {
64        let available: Vec<usize> = (0..npop).filter(|&idx| idx != i && idx != r1).collect();
65        if available.is_empty() {
66            pop.row(r1).to_owned()
67        } else {
68            pop.row(available[rng.random_range(0..available.len())])
69                .to_owned()
70        }
71    };
72
73    Zip::from(pop.row(i))
74        .and(pop.row(pbest_idx))
75        .and(pop.row(r1))
76        .and(r2.view())
77        .map_collect(|&curr, &pbest, &x1, &x2| curr + f * (pbest - curr) + f * (x1 - x2))
78}
79
80/// Compute p_best_size from p parameter (0 < p <= 1)
81/// p=0.1 means top 10% of population
82#[inline]
83pub(crate) fn compute_pbest_size(p: f64, npop: usize) -> usize {
84    ((p * npop as f64).ceil() as usize).max(1).min(npop)
85}
86
87#[cfg(test)]
88mod tests {
89    use super::*;
90    use ndarray::array;
91
92    #[test]
93    fn test_pbest_size() {
94        assert_eq!(compute_pbest_size(0.1, 100), 10);
95        assert_eq!(compute_pbest_size(0.1, 50), 5);
96        assert_eq!(compute_pbest_size(0.1, 10), 1);
97        assert_eq!(compute_pbest_size(0.5, 10), 5);
98    }
99
100    #[test]
101    fn test_mutant_basic() {
102        let pop = array![[0.0, 0.0], [1.0, 1.0], [2.0, 2.0], [3.0, 3.0]];
103        let sorted_indices = vec![3, 2, 1, 0];
104
105        let mut rng = rand::rng();
106        let mutant = mutant_current_to_pbest1(0, &pop, &sorted_indices, 2, None, 0.5, &mut rng);
107
108        assert_eq!(mutant.len(), 2);
109    }
110}
math_audio_optimisation/mutant_current_to_pbest1.rs

math_audio_optimisation/
mutant_current_to_pbest1.rs
