1use rand::prelude::*;
19use std::sync::atomic::{AtomicBool, Ordering};
20use std::sync::Arc;
21use std::time::Duration;
22
23use crate::error::{check_above, check_at_least, check_range, Result};
24use crate::result::StopReason;
25use crate::timing::{expired, Timer};
26
27#[cfg(feature = "serde")]
28use serde::{Deserialize, Serialize};
29
30#[derive(Debug, Clone, Copy, Default)]
32#[cfg_attr(feature = "serde", derive(Serialize, Deserialize))]
33pub enum CoolingSchedule {
34 #[default]
36 Geometric,
37 Linear,
39 Adaptive,
41 LundyMees,
43}
44
45#[derive(Debug, Clone)]
47#[cfg_attr(feature = "serde", derive(Serialize, Deserialize))]
48pub struct SaConfig {
49 pub initial_temp: f64,
51 pub final_temp: f64,
53 pub cooling_rate: f64,
55 pub iterations_per_temp: usize,
57 pub max_iterations: Option<u64>,
59 pub cooling_schedule: CoolingSchedule,
61 pub time_limit: Option<Duration>,
63 pub target_fitness: Option<f64>,
65 pub enable_reheating: bool,
67 pub reheat_threshold: u64,
69 pub reheat_factor: f64,
71}
72
73impl Default for SaConfig {
74 fn default() -> Self {
75 Self {
76 initial_temp: 1000.0,
77 final_temp: 0.001,
78 cooling_rate: 0.95,
79 iterations_per_temp: 100,
80 max_iterations: Some(100_000),
81 cooling_schedule: CoolingSchedule::Geometric,
82 time_limit: None,
83 target_fitness: None,
84 enable_reheating: false,
85 reheat_threshold: 1000,
86 reheat_factor: 2.0,
87 }
88 }
89}
90
91impl SaConfig {
92 pub fn new() -> Self {
94 Self::default()
95 }
96
97 pub fn with_initial_temp(mut self, temp: f64) -> Self {
99 self.initial_temp = temp;
100 self
101 }
102
103 pub fn with_final_temp(mut self, temp: f64) -> Self {
105 self.final_temp = temp;
106 self
107 }
108
109 pub fn with_cooling_rate(mut self, rate: f64) -> Self {
111 self.cooling_rate = rate;
112 self
113 }
114
115 pub fn validate(&self) -> Result<()> {
124 check_range("cooling_rate", self.cooling_rate, 0.001, 0.9999)?;
125 check_above("initial_temp", self.initial_temp, 0.0)?;
126 check_above("final_temp", self.final_temp, 0.0)?;
127 check_range("final_temp", self.final_temp, 0.0, self.initial_temp)?;
128 check_at_least("iterations_per_temp", self.iterations_per_temp as f64, 1.0)?;
129 check_at_least("reheat_factor", self.reheat_factor, 1.0)?;
130 Ok(())
131 }
132
133 pub fn with_iterations_per_temp(mut self, iterations: usize) -> Self {
135 self.iterations_per_temp = iterations;
136 self
137 }
138
139 pub fn with_max_iterations(mut self, iterations: u64) -> Self {
141 self.max_iterations = Some(iterations);
142 self
143 }
144
145 pub fn with_cooling_schedule(mut self, schedule: CoolingSchedule) -> Self {
147 self.cooling_schedule = schedule;
148 self
149 }
150
151 pub fn with_time_limit(mut self, duration: Duration) -> Self {
153 self.time_limit = Some(duration);
154 self
155 }
156
157 pub fn with_target_fitness(mut self, fitness: f64) -> Self {
159 self.target_fitness = Some(fitness);
160 self
161 }
162
163 pub fn with_reheating(mut self, threshold: u64, factor: f64) -> Self {
165 self.enable_reheating = true;
166 self.reheat_threshold = threshold;
167 self.reheat_factor = factor;
168 self
169 }
170}
171
172pub trait SaSolution: Clone + Send + Sync {
174 fn objective(&self) -> f64;
177
178 fn set_objective(&mut self, value: f64);
180}
181
182#[derive(Debug, Clone, Copy, PartialEq, Eq)]
184pub enum NeighborhoodOperator {
185 Swap,
187 Relocate,
189 Inversion,
191 Rotation,
193 Chain,
195 MirrorFlip,
200}
201
202pub trait SaProblem: Send + Sync {
204 type Solution: SaSolution;
206
207 fn initial_solution<R: Rng>(&self, rng: &mut R) -> Self::Solution;
209
210 fn neighbor<R: Rng>(
212 &self,
213 solution: &Self::Solution,
214 operator: NeighborhoodOperator,
215 rng: &mut R,
216 ) -> Self::Solution;
217
218 fn evaluate(&self, solution: &mut Self::Solution);
220
221 fn available_operators(&self) -> Vec<NeighborhoodOperator> {
223 vec![
224 NeighborhoodOperator::Swap,
225 NeighborhoodOperator::Relocate,
226 NeighborhoodOperator::Inversion,
227 ]
228 }
229
230 fn on_temperature_change(
232 &self,
233 _temperature: f64,
234 _iteration: u64,
235 _best: &Self::Solution,
236 _current: &Self::Solution,
237 ) {
238 }
240}
241
242#[derive(Debug, Clone)]
244pub struct SaProgress {
245 pub temperature: f64,
247 pub iteration: u64,
249 pub best_fitness: f64,
251 pub current_fitness: f64,
253 pub acceptance_rate: f64,
255 pub elapsed: Duration,
257 pub running: bool,
259}
260
261#[derive(Debug, Clone)]
263pub struct SaResult<S: SaSolution> {
264 pub best: S,
266 pub final_temperature: f64,
268 pub iterations: u64,
270 pub elapsed: Duration,
272 pub stop_reason: StopReason,
275 pub reheat_count: u32,
277 pub history: Vec<f64>,
279}
280
281pub struct SaRunner<P: SaProblem> {
283 config: SaConfig,
284 problem: P,
285 cancelled: Arc<AtomicBool>,
286}
287
288impl<P: SaProblem> SaRunner<P> {
289 pub fn new(config: SaConfig, problem: P) -> Result<Self> {
294 Self::with_cancellation(config, problem, Arc::new(AtomicBool::new(false)))
295 }
296
297 pub fn with_cancellation(
303 config: SaConfig,
304 problem: P,
305 cancelled: Arc<AtomicBool>,
306 ) -> Result<Self> {
307 config.validate()?;
308 Ok(Self {
309 config,
310 problem,
311 cancelled,
312 })
313 }
314
315 pub fn cancel_handle(&self) -> Arc<AtomicBool> {
317 self.cancelled.clone()
318 }
319
320 pub fn run(&self) -> SaResult<P::Solution> {
322 self.run_with_rng(&mut rand::rng())
323 }
324
325 pub fn run_with_rng<R: Rng>(&self, rng: &mut R) -> SaResult<P::Solution> {
327 self.run_with_rng_and_progress::<R, fn(SaProgress)>(rng, None)
328 }
329
330 pub fn run_with_rng_and_progress<R: Rng, F>(
334 &self,
335 rng: &mut R,
336 progress_callback: Option<F>,
337 ) -> SaResult<P::Solution>
338 where
339 F: Fn(SaProgress),
340 {
341 let start = Timer::now();
342 let report = |progress: SaProgress| {
343 if let Some(ref callback) = progress_callback {
344 callback(progress);
345 }
346 };
347 let mut history = Vec::new();
348
349 let mut current = self.problem.initial_solution(rng);
351 self.problem.evaluate(&mut current);
352 let mut best = current.clone();
353 let mut best_fitness = best.objective();
354
355 let mut temperature = self.config.initial_temp;
356 let mut iteration = 0u64;
357 let reached = |fitness: f64| self.config.target_fitness.is_some_and(|t| fitness >= t);
358 let mut reheat_count = 0u32;
359 let mut stagnation_count = 0u64;
360
361 let operators = self.problem.available_operators();
362 let temp_delta = if matches!(self.config.cooling_schedule, CoolingSchedule::Linear) {
363 (self.config.initial_temp - self.config.final_temp)
364 / (self.config.max_iterations.unwrap_or(10000) as f64
365 / self.config.iterations_per_temp as f64)
366 } else {
367 0.0
368 };
369
370 let mut accepted_count = 0usize;
372 let mut total_count = 0usize;
373
374 let stop_reason = loop {
377 if reached(best_fitness) {
378 break StopReason::TargetReached;
379 }
380 if temperature <= self.config.final_temp {
381 break StopReason::Completed;
382 }
383 if self
384 .config
385 .max_iterations
386 .is_some_and(|max| iteration >= max)
387 {
388 break StopReason::IterationLimit;
389 }
390 if self.cancelled.load(Ordering::Relaxed) {
391 break StopReason::Cancelled;
392 }
393 if expired(&start, self.config.time_limit) {
394 break StopReason::TimeLimit;
395 }
396
397 for _ in 0..self.config.iterations_per_temp {
400 if expired(&start, self.config.time_limit) {
401 break;
402 }
403 iteration += 1;
404 total_count += 1;
405
406 let operator = operators[rng.random_range(0..operators.len())];
408
409 let mut neighbor = self.problem.neighbor(¤t, operator, rng);
411 self.problem.evaluate(&mut neighbor);
412
413 let current_obj = current.objective();
414 let neighbor_obj = neighbor.objective();
415 let delta = neighbor_obj - current_obj;
416
417 let accept = if delta >= 0.0 {
419 true
421 } else {
422 let probability = (delta / temperature).exp();
424 rng.random::<f64>() < probability
425 };
426
427 if accept {
428 accepted_count += 1;
429 current = neighbor;
430
431 if current.objective() > best_fitness {
433 best = current.clone();
434 best_fitness = best.objective();
435 stagnation_count = 0;
436 } else {
437 stagnation_count += 1;
438 }
439 } else {
440 stagnation_count += 1;
441 }
442
443 report(SaProgress {
447 temperature,
448 iteration,
449 best_fitness,
450 current_fitness: current.objective(),
451 acceptance_rate: accepted_count as f64 / total_count.max(1) as f64,
452 elapsed: start.elapsed(),
453 running: true,
454 });
455
456 if let Some(max) = self.config.max_iterations {
458 if iteration >= max {
459 break;
460 }
461 }
462 }
463
464 history.push(best_fitness);
466
467 self.problem
469 .on_temperature_change(temperature, iteration, &best, ¤t);
470
471 if self.config.enable_reheating && stagnation_count >= self.config.reheat_threshold {
473 temperature *= self.config.reheat_factor;
474 temperature = temperature.min(self.config.initial_temp);
475 stagnation_count = 0;
476 reheat_count += 1;
477 }
478
479 temperature = self.cool_down(temperature, temp_delta, accepted_count, total_count);
481
482 accepted_count = 0;
484 total_count = 0;
485 };
486
487 history.push(best_fitness);
489 report(SaProgress {
490 temperature,
491 iteration,
492 best_fitness,
493 current_fitness: current.objective(),
494 acceptance_rate: 0.0,
495 elapsed: start.elapsed(),
496 running: false,
497 });
498
499 SaResult {
500 best,
501 final_temperature: temperature,
502 iterations: iteration,
503 elapsed: start.elapsed(),
504 stop_reason,
505 reheat_count,
506 history,
507 }
508 }
509
510 fn cool_down(&self, current_temp: f64, delta: f64, accepted: usize, total: usize) -> f64 {
512 match self.config.cooling_schedule {
513 CoolingSchedule::Geometric => current_temp * self.config.cooling_rate,
514 CoolingSchedule::Linear => (current_temp - delta).max(self.config.final_temp),
515 CoolingSchedule::Adaptive => {
516 let acceptance_rate = if total > 0 {
518 accepted as f64 / total as f64
519 } else {
520 0.5
521 };
522
523 let adjusted_rate = if acceptance_rate > 0.5 {
525 self.config.cooling_rate * 0.95 } else if acceptance_rate < 0.1 {
527 self.config.cooling_rate.powf(0.5) } else {
529 self.config.cooling_rate
530 };
531
532 current_temp * adjusted_rate
533 }
534 CoolingSchedule::LundyMees => {
535 current_temp / (1.0 + self.config.cooling_rate * current_temp)
537 }
538 }
539 }
540}
541
542#[derive(Debug, Clone)]
544pub struct PermutationSolution {
545 pub sequence: Vec<usize>,
547 pub rotations: Vec<usize>,
549 pub rotation_options: usize,
551 pub mirrors: Vec<bool>,
554 objective: f64,
556}
557
558impl PermutationSolution {
559 pub fn new(size: usize, rotation_options: usize) -> Self {
561 Self {
562 sequence: (0..size).collect(),
563 rotations: vec![0; size],
564 rotation_options,
565 mirrors: vec![false; size],
566 objective: f64::NEG_INFINITY,
567 }
568 }
569
570 pub fn random<R: Rng>(size: usize, rotation_options: usize, rng: &mut R) -> Self {
572 let mut sequence: Vec<usize> = (0..size).collect();
573 sequence.shuffle(rng);
574
575 let rotations: Vec<usize> = (0..size)
576 .map(|_| rng.random_range(0..rotation_options.max(1)))
577 .collect();
578
579 let mirrors: Vec<bool> = (0..size).map(|_| rng.random()).collect();
580
581 Self {
582 sequence,
583 rotations,
584 rotation_options,
585 mirrors,
586 objective: f64::NEG_INFINITY,
587 }
588 }
589
590 pub fn len(&self) -> usize {
592 self.sequence.len()
593 }
594
595 pub fn is_empty(&self) -> bool {
597 self.sequence.is_empty()
598 }
599
600 pub fn apply_swap<R: Rng>(&self, rng: &mut R) -> Self {
602 let mut result = self.clone();
603 if result.sequence.len() < 2 {
604 return result;
605 }
606
607 let i = rng.random_range(0..result.sequence.len());
608 let j = rng.random_range(0..result.sequence.len());
609 result.sequence.swap(i, j);
610 result.objective = f64::NEG_INFINITY;
611 result
612 }
613
614 pub fn apply_relocate<R: Rng>(&self, rng: &mut R) -> Self {
616 let mut result = self.clone();
617 if result.sequence.len() < 2 {
618 return result;
619 }
620
621 let from = rng.random_range(0..result.sequence.len());
622 let to = rng.random_range(0..result.sequence.len());
623
624 if from != to {
625 let elem = result.sequence.remove(from);
626 let insert_pos = if to > from { to - 1 } else { to };
627 result
628 .sequence
629 .insert(insert_pos.min(result.sequence.len()), elem);
630 }
631
632 result.objective = f64::NEG_INFINITY;
633 result
634 }
635
636 pub fn apply_inversion<R: Rng>(&self, rng: &mut R) -> Self {
638 let mut result = self.clone();
639 let n = result.sequence.len();
640 if n < 2 {
641 return result;
642 }
643
644 let (mut p1, mut p2) = (rng.random_range(0..n), rng.random_range(0..n));
645 if p1 > p2 {
646 std::mem::swap(&mut p1, &mut p2);
647 }
648
649 result.sequence[p1..=p2].reverse();
650 result.objective = f64::NEG_INFINITY;
651 result
652 }
653
654 pub fn apply_rotation<R: Rng>(&self, rng: &mut R) -> Self {
656 let mut result = self.clone();
657 if result.rotations.is_empty() || result.rotation_options <= 1 {
658 return result;
659 }
660
661 let idx = rng.random_range(0..result.rotations.len());
662 result.rotations[idx] = rng.random_range(0..result.rotation_options);
663 result.objective = f64::NEG_INFINITY;
664 result
665 }
666
667 pub fn apply_mirror_flip<R: Rng>(&self, rng: &mut R) -> Self {
672 let mut result = self.clone();
673 if result.mirrors.is_empty() {
674 return result;
675 }
676
677 let idx = rng.random_range(0..result.mirrors.len());
678 result.mirrors[idx] = !result.mirrors[idx];
679 result.objective = f64::NEG_INFINITY;
680 result
681 }
682
683 pub fn apply_chain<R: Rng>(&self, rng: &mut R) -> Self {
685 let mut result = self.clone();
686 let n = result.sequence.len();
687 if n < 4 {
688 return self.apply_swap(rng);
690 }
691
692 let mut positions: Vec<usize> = (0..n).collect();
694 positions.shuffle(rng);
695 let mut selected: Vec<usize> = positions.into_iter().take(3).collect();
696 selected.sort();
697
698 let (p1, p2, p3) = (selected[0], selected[1], selected[2]);
699
700 let seg1: Vec<usize> = result.sequence[..p1].to_vec();
703 let seg2: Vec<usize> = result.sequence[p1..p2].to_vec();
704 let seg3: Vec<usize> = result.sequence[p2..p3].to_vec();
705 let seg4: Vec<usize> = result.sequence[p3..].to_vec();
706
707 result.sequence = [seg1, seg3, seg2, seg4].concat();
708 result.objective = f64::NEG_INFINITY;
709 result
710 }
711}
712
713impl SaSolution for PermutationSolution {
714 fn objective(&self) -> f64 {
715 self.objective
716 }
717
718 fn set_objective(&mut self, value: f64) {
719 self.objective = value;
720 }
721}
722
723#[cfg(test)]
724mod tests {
725 use super::*;
726
727 struct SimpleMaxProblem {
728 size: usize,
729 }
730
731 impl SaProblem for SimpleMaxProblem {
732 type Solution = PermutationSolution;
733
734 fn initial_solution<R: Rng>(&self, rng: &mut R) -> Self::Solution {
735 PermutationSolution::random(self.size, 1, rng)
736 }
737
738 fn neighbor<R: Rng>(
739 &self,
740 solution: &Self::Solution,
741 operator: NeighborhoodOperator,
742 rng: &mut R,
743 ) -> Self::Solution {
744 match operator {
745 NeighborhoodOperator::Swap => solution.apply_swap(rng),
746 NeighborhoodOperator::Relocate => solution.apply_relocate(rng),
747 NeighborhoodOperator::Inversion => solution.apply_inversion(rng),
748 NeighborhoodOperator::Rotation => solution.apply_rotation(rng),
749 NeighborhoodOperator::Chain => solution.apply_chain(rng),
750 NeighborhoodOperator::MirrorFlip => solution.apply_mirror_flip(rng),
751 }
752 }
753
754 fn evaluate(&self, solution: &mut Self::Solution) {
755 let mut inversions = 0i64;
758 for i in 0..solution.sequence.len() {
759 for j in (i + 1)..solution.sequence.len() {
760 if solution.sequence[i] > solution.sequence[j] {
761 inversions += 1;
762 }
763 }
764 }
765 solution.set_objective(-inversions as f64);
766 }
767 }
768
769 #[test]
770 fn test_sa_basic() {
771 let config = SaConfig::default()
772 .with_initial_temp(100.0)
773 .with_final_temp(0.1)
774 .with_cooling_rate(0.9)
775 .with_iterations_per_temp(50)
776 .with_max_iterations(5000);
777
778 let problem = SimpleMaxProblem { size: 10 };
779 let runner = SaRunner::new(config, problem).expect("valid config");
780 let result = runner.run();
781
782 assert!(result.best.objective() > -20.0);
784 assert!(result.iterations > 0);
785 }
786
787 #[test]
788 fn test_cooling_schedules() {
789 let problem = SimpleMaxProblem { size: 5 };
790
791 for schedule in [
792 CoolingSchedule::Geometric,
793 CoolingSchedule::Linear,
794 CoolingSchedule::Adaptive,
795 CoolingSchedule::LundyMees,
796 ] {
797 let config = SaConfig::default()
798 .with_cooling_schedule(schedule)
799 .with_max_iterations(1000);
800
801 let runner = SaRunner::new(config, problem.clone()).expect("valid config");
802 let result = runner.run();
803
804 assert!(result.iterations > 0);
806 }
807 }
808
809 #[test]
810 fn test_neighborhood_operators() {
811 let mut rng = rand::rng();
812 let solution = PermutationSolution::random(10, 4, &mut rng);
813
814 let swap = solution.apply_swap(&mut rng);
816 let relocate = solution.apply_relocate(&mut rng);
817 let inversion = solution.apply_inversion(&mut rng);
818 let rotation = solution.apply_rotation(&mut rng);
819 let chain = solution.apply_chain(&mut rng);
820
821 for sol in [&swap, &relocate, &inversion, &rotation, &chain] {
822 let mut sorted = sol.sequence.clone();
823 sorted.sort();
824 assert_eq!(sorted, (0..10).collect::<Vec<_>>());
825 }
826 }
827
828 #[test]
829 fn test_reheating() {
830 let config = SaConfig::default()
831 .with_initial_temp(10.0)
832 .with_final_temp(0.1)
833 .with_max_iterations(500)
834 .with_reheating(50, 1.5);
835
836 let problem = SimpleMaxProblem { size: 8 };
837 let runner = SaRunner::new(config, problem).expect("valid config");
838 let result = runner.run();
839
840 assert!(result.iterations > 0);
842 }
843
844 impl Clone for SimpleMaxProblem {
845 fn clone(&self) -> Self {
846 Self { size: self.size }
847 }
848 }
849}