Skip to main content

scirs2_optimize/multi_objective/algorithms/
nsga3.rs

1//! NSGA-III (Non-dominated Sorting Genetic Algorithm III)
2//!
3//! Reference-point-based many-objective evolutionary algorithm.
4//! NSGA-III extends NSGA-II for many-objective optimization (>3 objectives)
5//! by replacing crowding distance with reference-point-based selection.
6//!
7//! Key features:
8//! - Das-Dennis reference point generation
9//! - Adaptive normalization of objective values
10//! - Reference-point association using perpendicular distance
11//! - Niching-based selection for population diversity
12//!
13//! # References
14//!
15//! - Deb & Jain, "An Evolutionary Many-Objective Optimization Algorithm Using
16//!   Reference-Point-Based Nondominated Sorting Approach, Part I", IEEE TEC 2014
17
18use super::{utils, MultiObjectiveConfig, MultiObjectiveOptimizer};
19use crate::error::OptimizeError;
20use crate::multi_objective::crossover::{CrossoverOperator, SimulatedBinaryCrossover};
21use crate::multi_objective::mutation::{MutationOperator, PolynomialMutation};
22use crate::multi_objective::selection::{SelectionOperator, TournamentSelection};
23use crate::multi_objective::solutions::{MultiObjectiveResult, MultiObjectiveSolution, Population};
24use scirs2_core::ndarray::{Array1, ArrayView1};
25use scirs2_core::random::rngs::StdRng;
26use scirs2_core::random::{Rng, SeedableRng};
27use std::cmp::Ordering;
28
29/// NSGA-III optimizer for many-objective optimization
30pub struct NSGAIII {
31    config: MultiObjectiveConfig,
32    n_objectives: usize,
33    n_variables: usize,
34    /// Structured reference points (Das-Dennis)
35    reference_points: Vec<Array1<f64>>,
36    population: Population,
37    generation: usize,
38    n_evaluations: usize,
39    rng: StdRng,
40    crossover: SimulatedBinaryCrossover,
41    mutation: PolynomialMutation,
42    selection: TournamentSelection,
43    /// Ideal point (best value per objective)
44    ideal_point: Array1<f64>,
45    /// Niche count for each reference point
46    niche_count: Vec<usize>,
47    convergence_history: Vec<f64>,
48}
49
50impl NSGAIII {
51    /// Create new NSGA-III optimizer with default reference points
52    pub fn new(population_size: usize, n_objectives: usize, n_variables: usize) -> Self {
53        let config = MultiObjectiveConfig {
54            population_size,
55            ..Default::default()
56        };
57        Self::with_config(config, n_objectives, n_variables, None)
58    }
59
60    /// Create NSGA-III with full configuration and optional custom reference points
61    pub fn with_config(
62        config: MultiObjectiveConfig,
63        n_objectives: usize,
64        n_variables: usize,
65        custom_reference_points: Option<Vec<Array1<f64>>>,
66    ) -> Self {
67        let seed = config.random_seed.unwrap_or_else(|| {
68            use std::time::{SystemTime, UNIX_EPOCH};
69            SystemTime::now()
70                .duration_since(UNIX_EPOCH)
71                .map(|d| d.as_secs())
72                .unwrap_or(42)
73        });
74
75        let rng = StdRng::seed_from_u64(seed);
76
77        let crossover =
78            SimulatedBinaryCrossover::new(config.crossover_eta, config.crossover_probability);
79        let mutation = PolynomialMutation::new(config.mutation_probability, config.mutation_eta);
80        let selection = TournamentSelection::new(2);
81
82        // Generate reference points
83        let reference_points = custom_reference_points.unwrap_or_else(|| {
84            let n_partitions = Self::auto_partitions(n_objectives, config.population_size);
85            utils::generate_das_dennis_points(n_objectives, n_partitions)
86        });
87
88        let niche_count = vec![0; reference_points.len()];
89
90        Self {
91            config,
92            n_objectives,
93            n_variables,
94            reference_points,
95            population: Population::new(),
96            generation: 0,
97            n_evaluations: 0,
98            rng,
99            crossover,
100            mutation,
101            selection,
102            ideal_point: Array1::from_elem(n_objectives, f64::INFINITY),
103            niche_count,
104            convergence_history: Vec::new(),
105        }
106    }
107
108    /// Automatically determine the number of partitions based on objectives and pop size
109    fn auto_partitions(n_objectives: usize, pop_size: usize) -> usize {
110        // Heuristic: choose partitions so Das-Dennis points ~ pop_size
111        // Number of points = C(H + M - 1, M - 1) where H = partitions, M = objectives
112        if n_objectives <= 2 {
113            return pop_size.max(4);
114        }
115        // For higher dimensions, use fewer partitions
116        let mut p = 1;
117        loop {
118            let n_points = binomial_coefficient(p + n_objectives - 1, n_objectives - 1);
119            if n_points >= pop_size / 2 {
120                return p;
121            }
122            p += 1;
123            if p > 50 {
124                return p;
125            }
126        }
127    }
128
129    /// Evaluate a single individual
130    fn evaluate_individual<F>(
131        &mut self,
132        variables: &Array1<f64>,
133        objective_function: &F,
134    ) -> Result<Array1<f64>, OptimizeError>
135    where
136        F: Fn(&ArrayView1<f64>) -> Array1<f64> + Send + Sync,
137    {
138        self.n_evaluations += 1;
139
140        if let Some(max_evals) = self.config.max_evaluations {
141            if self.n_evaluations > max_evals {
142                return Err(OptimizeError::MaxEvaluationsReached);
143            }
144        }
145
146        let objectives = objective_function(&variables.view());
147        if objectives.len() != self.n_objectives {
148            return Err(OptimizeError::InvalidInput(format!(
149                "Expected {} objectives, got {}",
150                self.n_objectives,
151                objectives.len()
152            )));
153        }
154
155        Ok(objectives)
156    }
157
158    /// Update the ideal point from a set of solutions
159    fn update_ideal_point(&mut self, solutions: &[MultiObjectiveSolution]) {
160        for sol in solutions {
161            for (i, &obj) in sol.objectives.iter().enumerate() {
162                if obj < self.ideal_point[i] {
163                    self.ideal_point[i] = obj;
164                }
165            }
166        }
167    }
168
169    /// Normalize objectives using ideal and extreme points (ASF-based)
170    fn normalize_objectives(&self, solutions: &[MultiObjectiveSolution]) -> Vec<Array1<f64>> {
171        let n = solutions.len();
172        if n == 0 {
173            return vec![];
174        }
175
176        // Translate by ideal point
177        let translated: Vec<Array1<f64>> = solutions
178            .iter()
179            .map(|sol| &sol.objectives - &self.ideal_point)
180            .collect();
181
182        // Find extreme points using Achievement Scalarizing Function (ASF)
183        let mut extreme_points = Vec::with_capacity(self.n_objectives);
184        for j in 0..self.n_objectives {
185            let mut best_asf = f64::INFINITY;
186            let mut best_idx = 0;
187            for (i, t) in translated.iter().enumerate() {
188                // ASF: max(t_k / w_k) where w is the j-th axis weight
189                let asf = (0..self.n_objectives)
190                    .map(|k| {
191                        if k == j {
192                            t[k] // weight = 1 for target objective
193                        } else {
194                            t[k] * 1e6 // large weight for others
195                        }
196                    })
197                    .fold(f64::NEG_INFINITY, f64::max);
198                if asf < best_asf {
199                    best_asf = asf;
200                    best_idx = i;
201                }
202            }
203            extreme_points.push(translated[best_idx].clone());
204        }
205
206        // Compute intercepts from the extreme points
207        let mut intercepts = Array1::from_elem(self.n_objectives, 1.0);
208        if extreme_points.len() == self.n_objectives {
209            // Build the matrix of extreme points and solve for intercepts
210            // Simplified: use the diagonal elements as intercepts
211            for (j, ep) in extreme_points.iter().enumerate() {
212                let val = ep[j];
213                if val > 1e-10 {
214                    intercepts[j] = val;
215                }
216            }
217        }
218
219        // Normalize each translated objective
220        translated
221            .iter()
222            .map(|t| {
223                let mut normalized = Array1::zeros(self.n_objectives);
224                for j in 0..self.n_objectives {
225                    normalized[j] = if intercepts[j] > 1e-10 {
226                        t[j] / intercepts[j]
227                    } else {
228                        t[j]
229                    };
230                }
231                normalized
232            })
233            .collect()
234    }
235
236    /// Associate each solution to its nearest reference point
237    /// Returns (reference_point_index, perpendicular_distance) for each solution
238    fn associate_to_reference_points(
239        &self,
240        normalized_objectives: &[Array1<f64>],
241    ) -> Vec<(usize, f64)> {
242        normalized_objectives
243            .iter()
244            .map(|obj| {
245                let mut min_dist = f64::INFINITY;
246                let mut min_idx = 0;
247
248                for (rp_idx, rp) in self.reference_points.iter().enumerate() {
249                    let dist = perpendicular_distance(obj, rp);
250                    if dist < min_dist {
251                        min_dist = dist;
252                        min_idx = rp_idx;
253                    }
254                }
255
256                (min_idx, min_dist)
257            })
258            .collect()
259    }
260
261    /// Niching-based selection from the last front
262    fn niching_selection(
263        &mut self,
264        last_front_indices: &[usize],
265        all_solutions: &[MultiObjectiveSolution],
266        associations: &[(usize, f64)],
267        selected_so_far: &[usize],
268        remaining: usize,
269    ) -> Vec<usize> {
270        // Compute niche counts for already-selected solutions
271        let mut niche_count = vec![0usize; self.reference_points.len()];
272        for &idx in selected_so_far {
273            let (rp_idx, _) = associations[idx];
274            niche_count[rp_idx] += 1;
275        }
276
277        let mut selected = Vec::with_capacity(remaining);
278        let mut available: Vec<usize> = last_front_indices.to_vec();
279
280        for _ in 0..remaining {
281            if available.is_empty() {
282                break;
283            }
284
285            // Find the reference point with minimum niche count
286            // among those that have at least one member in the last front
287            let relevant_rps: Vec<usize> = available
288                .iter()
289                .map(|&idx| associations[idx].0)
290                .collect::<std::collections::HashSet<_>>()
291                .into_iter()
292                .collect();
293
294            if relevant_rps.is_empty() {
295                break;
296            }
297
298            let min_niche = relevant_rps
299                .iter()
300                .map(|&rp| niche_count[rp])
301                .min()
302                .unwrap_or(0);
303
304            // Collect reference points with minimum niche count
305            let min_niche_rps: Vec<usize> = relevant_rps
306                .iter()
307                .filter(|&&rp| niche_count[rp] == min_niche)
308                .copied()
309                .collect();
310
311            // Randomly pick one
312            let chosen_rp_idx = self.rng.random_range(0..min_niche_rps.len());
313            let chosen_rp = min_niche_rps[chosen_rp_idx];
314
315            // Find members of the last front associated with this reference point
316            let rp_members: Vec<usize> = available
317                .iter()
318                .filter(|&&idx| associations[idx].0 == chosen_rp)
319                .copied()
320                .collect();
321
322            if rp_members.is_empty() {
323                continue;
324            }
325
326            let chosen_member = if min_niche == 0 {
327                // Pick the one with smallest perpendicular distance
328                *rp_members
329                    .iter()
330                    .min_by(|&&a, &&b| {
331                        associations[a]
332                            .1
333                            .partial_cmp(&associations[b].1)
334                            .unwrap_or(Ordering::Equal)
335                    })
336                    .unwrap_or(&rp_members[0])
337            } else {
338                // Randomly pick
339                rp_members[self.rng.random_range(0..rp_members.len())]
340            };
341
342            selected.push(chosen_member);
343            niche_count[chosen_rp] += 1;
344            available.retain(|&idx| idx != chosen_member);
345        }
346
347        selected
348    }
349
350    /// Create offspring through crossover and mutation
351    fn create_offspring<F>(
352        &mut self,
353        objective_function: &F,
354    ) -> Result<Vec<MultiObjectiveSolution>, OptimizeError>
355    where
356        F: Fn(&ArrayView1<f64>) -> Array1<f64> + Send + Sync,
357    {
358        let mut offspring = Vec::new();
359        let pop_solutions = self.population.solutions().to_vec();
360
361        while offspring.len() < self.config.population_size {
362            let selected = self.selection.select(&pop_solutions, 2);
363            if selected.len() < 2 {
364                break;
365            }
366
367            let p1_vars = selected[0].variables.as_slice().unwrap_or(&[]);
368            let p2_vars = selected[1].variables.as_slice().unwrap_or(&[]);
369
370            let (mut c1_vars, mut c2_vars) = self.crossover.crossover(p1_vars, p2_vars);
371
372            let bounds: Vec<(f64, f64)> = if let Some((lower, upper)) = &self.config.bounds {
373                lower
374                    .iter()
375                    .zip(upper.iter())
376                    .map(|(&l, &u)| (l, u))
377                    .collect()
378            } else {
379                vec![(-1.0, 1.0); self.n_variables]
380            };
381
382            self.mutation.mutate(&mut c1_vars, &bounds);
383            self.mutation.mutate(&mut c2_vars, &bounds);
384
385            let c1_arr = Array1::from_vec(c1_vars);
386            let c1_obj = self.evaluate_individual(&c1_arr, objective_function)?;
387            offspring.push(MultiObjectiveSolution::new(c1_arr, c1_obj));
388
389            if offspring.len() < self.config.population_size {
390                let c2_arr = Array1::from_vec(c2_vars);
391                let c2_obj = self.evaluate_individual(&c2_arr, objective_function)?;
392                offspring.push(MultiObjectiveSolution::new(c2_arr, c2_obj));
393            }
394        }
395
396        Ok(offspring)
397    }
398
399    /// NSGA-III environmental selection using reference-point-based niching
400    fn environmental_selection(
401        &mut self,
402        combined: Vec<MultiObjectiveSolution>,
403    ) -> Vec<MultiObjectiveSolution> {
404        let target_size = self.config.population_size;
405        if combined.len() <= target_size {
406            return combined;
407        }
408
409        // Non-dominated sorting
410        let mut temp_pop = Population::from_solutions(combined.clone());
411        let fronts = temp_pop.non_dominated_sort();
412
413        // Fill until the critical front
414        let mut selected_indices: Vec<usize> = Vec::new();
415        let mut last_front_idx = 0;
416
417        for (fi, front) in fronts.iter().enumerate() {
418            if selected_indices.len() + front.len() <= target_size {
419                selected_indices.extend(front);
420            } else {
421                last_front_idx = fi;
422                break;
423            }
424        }
425
426        if selected_indices.len() >= target_size {
427            return selected_indices
428                .iter()
429                .take(target_size)
430                .map(|&i| combined[i].clone())
431                .collect();
432        }
433
434        // Need to select from the last (critical) front using niching
435        let remaining = target_size - selected_indices.len();
436        let last_front = &fronts[last_front_idx];
437
438        // Update ideal point
439        self.update_ideal_point(&combined);
440
441        // Normalize objectives for ALL solutions considered
442        let all_considered: Vec<usize> = selected_indices
443            .iter()
444            .chain(last_front.iter())
445            .copied()
446            .collect();
447        let all_solutions: Vec<MultiObjectiveSolution> = all_considered
448            .iter()
449            .map(|&i| combined[i].clone())
450            .collect();
451
452        let normalized = self.normalize_objectives(&all_solutions);
453        let associations = self.associate_to_reference_points(&normalized);
454
455        // Map indices back
456        let n_selected = selected_indices.len();
457        let last_front_local: Vec<usize> = (n_selected..all_solutions.len()).collect();
458        let selected_local: Vec<usize> = (0..n_selected).collect();
459
460        let niching_result = self.niching_selection(
461            &last_front_local,
462            &all_solutions,
463            &associations,
464            &selected_local,
465            remaining,
466        );
467
468        // Build final selection
469        let mut result: Vec<MultiObjectiveSolution> = selected_indices
470            .iter()
471            .map(|&i| combined[i].clone())
472            .collect();
473
474        for local_idx in niching_result {
475            let global_idx = all_considered[local_idx];
476            result.push(combined[global_idx].clone());
477        }
478
479        result
480    }
481
482    /// Calculate metrics for the current generation
483    fn calculate_metrics(&mut self) {
484        if let Some(ref_point) = &self.config.reference_point {
485            let pareto_front = self.population.extract_pareto_front();
486            let hv = utils::calculate_hypervolume(&pareto_front, ref_point);
487            self.convergence_history.push(hv);
488        }
489    }
490}
491
492impl MultiObjectiveOptimizer for NSGAIII {
493    fn optimize<F>(&mut self, objective_function: F) -> Result<MultiObjectiveResult, OptimizeError>
494    where
495        F: Fn(&ArrayView1<f64>) -> Array1<f64> + Send + Sync,
496    {
497        self.initialize_population()?;
498
499        // Generate and evaluate initial population
500        let initial_vars = utils::generate_random_population(
501            self.config.population_size,
502            self.n_variables,
503            &self.config.bounds,
504        );
505
506        let mut initial_solutions = Vec::new();
507        for vars in initial_vars {
508            let objs = self.evaluate_individual(&vars, &objective_function)?;
509            initial_solutions.push(MultiObjectiveSolution::new(vars, objs));
510        }
511
512        self.update_ideal_point(&initial_solutions);
513        self.population = Population::from_solutions(initial_solutions);
514
515        // Main evolution loop
516        while self.generation < self.config.max_generations {
517            if self.check_convergence() {
518                break;
519            }
520            self.evolve_generation(&objective_function)?;
521        }
522
523        // Extract results
524        let pareto_front = self.population.extract_pareto_front();
525        let hypervolume = self
526            .config
527            .reference_point
528            .as_ref()
529            .map(|rp| utils::calculate_hypervolume(&pareto_front, rp));
530
531        let mut result = MultiObjectiveResult::new(
532            pareto_front,
533            self.population.solutions().to_vec(),
534            self.n_evaluations,
535            self.generation,
536        );
537        result.hypervolume = hypervolume;
538        result.metrics.convergence_history = self.convergence_history.clone();
539        result.metrics.population_stats = self.population.calculate_statistics();
540
541        Ok(result)
542    }
543
544    fn evolve_generation<F>(&mut self, objective_function: &F) -> Result<(), OptimizeError>
545    where
546        F: Fn(&ArrayView1<f64>) -> Array1<f64> + Send + Sync,
547    {
548        let offspring = self.create_offspring(objective_function)?;
549
550        let mut combined = self.population.solutions().to_vec();
551        combined.extend(offspring);
552
553        let next_pop = self.environmental_selection(combined);
554        self.population = Population::from_solutions(next_pop);
555
556        self.generation += 1;
557        self.calculate_metrics();
558
559        Ok(())
560    }
561
562    fn initialize_population(&mut self) -> Result<(), OptimizeError> {
563        self.population.clear();
564        self.generation = 0;
565        self.n_evaluations = 0;
566        self.ideal_point = Array1::from_elem(self.n_objectives, f64::INFINITY);
567        self.niche_count = vec![0; self.reference_points.len()];
568        self.convergence_history.clear();
569        Ok(())
570    }
571
572    fn check_convergence(&self) -> bool {
573        if let Some(max_evals) = self.config.max_evaluations {
574            if self.n_evaluations >= max_evals {
575                return true;
576            }
577        }
578
579        if self.convergence_history.len() >= 10 {
580            let recent = &self.convergence_history[self.convergence_history.len() - 10..];
581            let max_hv = recent.iter().fold(f64::NEG_INFINITY, |a, &b| a.max(b));
582            let min_hv = recent.iter().fold(f64::INFINITY, |a, &b| a.min(b));
583            if (max_hv - min_hv) < self.config.tolerance {
584                return true;
585            }
586        }
587
588        false
589    }
590
591    fn get_population(&self) -> &Population {
592        &self.population
593    }
594
595    fn get_generation(&self) -> usize {
596        self.generation
597    }
598
599    fn get_evaluations(&self) -> usize {
600        self.n_evaluations
601    }
602
603    fn name(&self) -> &str {
604        "NSGA-III"
605    }
606}
607
608/// Compute perpendicular distance from a point to a reference line (through origin)
609fn perpendicular_distance(point: &Array1<f64>, direction: &Array1<f64>) -> f64 {
610    let dir_norm_sq = direction.dot(direction);
611    if dir_norm_sq < 1e-30 {
612        return point.dot(point).sqrt();
613    }
614
615    // Projection of point onto direction
616    let proj_scalar = point.dot(direction) / dir_norm_sq;
617    let proj = proj_scalar * direction;
618
619    // Distance = |point - projection|
620    let diff = point - &proj;
621    diff.dot(&diff).sqrt()
622}
623
624/// Compute binomial coefficient C(n, k)
625fn binomial_coefficient(n: usize, k: usize) -> usize {
626    if k > n {
627        return 0;
628    }
629    if k == 0 || k == n {
630        return 1;
631    }
632    let k = k.min(n - k);
633    let mut result = 1usize;
634    for i in 0..k {
635        result = result.saturating_mul(n - i);
636        result /= i + 1;
637    }
638    result
639}
640
641#[cfg(test)]
642mod tests {
643    use super::*;
644    use scirs2_core::ndarray::{array, s};
645
646    fn zdt1(x: &ArrayView1<f64>) -> Array1<f64> {
647        let f1 = x[0];
648        let g = 1.0 + 9.0 * x.slice(s![1..]).sum() / (x.len() - 1) as f64;
649        let f2 = g * (1.0 - (f1 / g).sqrt());
650        array![f1, f2]
651    }
652
653    fn dtlz1(x: &ArrayView1<f64>) -> Array1<f64> {
654        // 3-objective DTLZ1
655        let n = x.len();
656        let k = n - 2; // for 3 objectives, k = n - M + 1
657        let g: f64 = (0..k)
658            .map(|i| {
659                let xi = x[n - k + i];
660                (xi - 0.5).powi(2) - (20.0 * std::f64::consts::PI * (xi - 0.5)).cos()
661            })
662            .sum::<f64>();
663        let g = 100.0 * (k as f64 + g);
664
665        let f1 = 0.5 * x[0] * x[1] * (1.0 + g);
666        let f2 = 0.5 * x[0] * (1.0 - x[1]) * (1.0 + g);
667        let f3 = 0.5 * (1.0 - x[0]) * (1.0 + g);
668        array![f1, f2, f3]
669    }
670
671    #[test]
672    fn test_nsga3_creation() {
673        let nsga3 = NSGAIII::new(100, 3, 5);
674        assert_eq!(nsga3.n_objectives, 3);
675        assert_eq!(nsga3.n_variables, 5);
676        assert!(!nsga3.reference_points.is_empty());
677        assert_eq!(nsga3.generation, 0);
678    }
679
680    #[test]
681    fn test_nsga3_with_config() {
682        let mut config = MultiObjectiveConfig::default();
683        config.population_size = 50;
684        config.max_generations = 10;
685        config.random_seed = Some(42);
686
687        let nsga3 = NSGAIII::with_config(config, 2, 3, None);
688        assert_eq!(nsga3.n_objectives, 2);
689        assert!(!nsga3.reference_points.is_empty());
690    }
691
692    #[test]
693    fn test_nsga3_custom_reference_points() {
694        let rps = vec![
695            array![1.0, 0.0, 0.0],
696            array![0.0, 1.0, 0.0],
697            array![0.0, 0.0, 1.0],
698            array![0.5, 0.5, 0.0],
699            array![0.5, 0.0, 0.5],
700            array![0.0, 0.5, 0.5],
701            array![0.333, 0.333, 0.334],
702        ];
703
704        let config = MultiObjectiveConfig {
705            population_size: 20,
706            max_generations: 5,
707            random_seed: Some(42),
708            ..Default::default()
709        };
710
711        let nsga3 = NSGAIII::with_config(config, 3, 5, Some(rps.clone()));
712        assert_eq!(nsga3.reference_points.len(), 7);
713    }
714
715    #[test]
716    fn test_nsga3_optimize_zdt1() {
717        let mut config = MultiObjectiveConfig::default();
718        config.max_generations = 10;
719        config.population_size = 20;
720        config.bounds = Some((Array1::zeros(3), Array1::ones(3)));
721        config.random_seed = Some(42);
722
723        let mut nsga3 = NSGAIII::with_config(config, 2, 3, None);
724        let result = nsga3.optimize(zdt1);
725
726        assert!(result.is_ok());
727        let res = result.expect("should succeed");
728        assert!(res.success);
729        assert!(!res.pareto_front.is_empty());
730        assert!(res.n_evaluations > 0);
731    }
732
733    #[test]
734    fn test_nsga3_optimize_dtlz1() {
735        let mut config = MultiObjectiveConfig::default();
736        config.max_generations = 10;
737        config.population_size = 20;
738        config.bounds = Some((Array1::zeros(5), Array1::ones(5)));
739        config.random_seed = Some(42);
740
741        let mut nsga3 = NSGAIII::with_config(config, 3, 5, None);
742        let result = nsga3.optimize(dtlz1);
743
744        assert!(result.is_ok());
745        let res = result.expect("should succeed");
746        assert!(res.success);
747        assert!(!res.pareto_front.is_empty());
748    }
749
750    #[test]
751    fn test_nsga3_max_evaluations() {
752        let mut config = MultiObjectiveConfig::default();
753        config.max_generations = 1000;
754        config.max_evaluations = Some(50);
755        config.population_size = 10;
756        config.bounds = Some((Array1::zeros(3), Array1::ones(3)));
757        config.random_seed = Some(42);
758
759        let mut nsga3 = NSGAIII::with_config(config, 2, 3, None);
760        let result = nsga3.optimize(zdt1);
761        assert!(result.is_ok());
762        let res = result.expect("should succeed");
763        assert!(res.n_evaluations <= 60); // Allow slight overshoot
764    }
765
766    #[test]
767    fn test_perpendicular_distance() {
768        let point = array![1.0, 1.0];
769        let direction = array![1.0, 0.0];
770
771        let dist = perpendicular_distance(&point, &direction);
772        assert!(
773            (dist - 1.0).abs() < 1e-10,
774            "Distance should be 1.0, got {}",
775            dist
776        );
777
778        // Point on the line
779        let point2 = array![2.0, 0.0];
780        let dist2 = perpendicular_distance(&point2, &direction);
781        assert!(dist2.abs() < 1e-10, "Distance should be 0.0, got {}", dist2);
782    }
783
784    #[test]
785    fn test_binomial_coefficient() {
786        assert_eq!(binomial_coefficient(5, 2), 10);
787        assert_eq!(binomial_coefficient(10, 3), 120);
788        assert_eq!(binomial_coefficient(4, 0), 1);
789        assert_eq!(binomial_coefficient(4, 4), 1);
790        assert_eq!(binomial_coefficient(0, 0), 1);
791    }
792
793    #[test]
794    fn test_nsga3_name() {
795        let nsga3 = NSGAIII::new(50, 2, 3);
796        assert_eq!(nsga3.name(), "NSGA-III");
797    }
798
799    #[test]
800    fn test_nsga3_convergence_check() {
801        let config = MultiObjectiveConfig {
802            tolerance: 1e-10,
803            max_generations: 2,
804            population_size: 10,
805            ..Default::default()
806        };
807        let nsga3 = NSGAIII::with_config(config, 2, 2, None);
808        assert!(!nsga3.check_convergence());
809    }
810}