@@ 1,7 1,11 @@
use std::convert::Infallible;
-use eoa_lib::{constraints::LowerThanConstraintFunction, fitness::FitnessFunction};
-use nalgebra::{OVector, SVector};
+use eoa_lib::{
+ bounded::BoundedOVector, comparison::MinimizingOperator, constraints::{stochastic_ranking_evolution_algorithm, ConstraintFunction, LowerThanConstraintFunction}, crossover::{ArithmeticCrossover, BoundedCrossover, BoundedCrossoverStrategy}, evolution::EvolutionResult, fitness::FitnessFunction, initializer::{Initializer, RandomInitializer}, pairing::AdjacentPairing, perturbation::{BoundedPerturbation, BoundedPerturbationStrategy, MutationPerturbation, RandomDistributionPerturbation}, population::Population, replacement::{BestReplacement, GenerationalReplacement}, selection::{BestSelection, TournamentSelection}
+};
+use nalgebra::SVector;
+use rand::RngCore;
+use rand_distr::Normal;
pub struct ArbitraryFitness<const SIZE: usize> {
fun: Box<dyn Fn(SVector<f64, SIZE>) -> f64>
@@ 25,12 29,20 @@ impl<const SIZE: usize> FitnessFunction for ArbitraryFitness<SIZE> {
}
}
-fn problem_g06() -> (ArbitraryFitness<2>, [LowerThanConstraintFunction<SVector<f64, 2>, f64>; 2], f64) {
- (
- ArbitraryFitness::new(
+/// A constrained optimization problem with clear field definitions
+pub struct ConstrainedProblem<const CONSTRAINTS: usize> {
+ pub objective: ArbitraryFitness<2>,
+ pub constraints: [LowerThanConstraintFunction<SVector<f64, 2>, f64>; CONSTRAINTS],
+ pub bounds: (SVector<f64, 2>, SVector<f64, 2>), // (min, max)
+ pub optimal_value: f64,
+}
+
+fn problem_g06() -> ConstrainedProblem<2> {
+ ConstrainedProblem {
+ objective: ArbitraryFitness::new(
Box::new(|vec| (vec[0] - 10.0).powi(3) + (vec[1] - 20.0).powi(3))
),
- [
+ constraints: [
LowerThanConstraintFunction::new(
Box::new(|vec| -(vec[0] - 5.0).powi(2) - (vec[1] - 5.0).powi(2) + 100.0)
),
@@ 38,13 50,17 @@ fn problem_g06() -> (ArbitraryFitness<2>, [LowerThanConstraintFunction<SVector<f
Box::new(|vec| (vec[0] - 6.0).powi(2) + (vec[1] - 5.0).powi(2) - 82.81)
),
],
- -6961.8137558015
- )
+ bounds: (
+ SVector::<f64, 2>::new(0.0, 0.0), // min bounds
+ SVector::<f64, 2>::new(100.0, 100.0) // max bounds
+ ),
+ optimal_value: -6961.8137558015,
+ }
}
-fn problem_g08(eps: f64) -> (ArbitraryFitness<2>, [LowerThanConstraintFunction<SVector<f64, 2>, f64>; 2], f64) {
- (
- ArbitraryFitness::new(
+fn problem_g08() -> ConstrainedProblem<2> {
+ ConstrainedProblem {
+ objective: ArbitraryFitness::new(
Box::new(|vec| {
let num = (2.0 * std::f64::consts::PI * vec[0]).sin().powi(3)
* (2.0 * std::f64::consts::PI * vec[1]).sin();
@@ 52,7 68,7 @@ fn problem_g08(eps: f64) -> (ArbitraryFitness<2>, [LowerThanConstraintFunction<S
-num / den
})
),
- [
+ constraints: [
LowerThanConstraintFunction::new(
Box::new(move |vec| {
let x1 = vec[0];
@@ 68,19 84,23 @@ fn problem_g08(eps: f64) -> (ArbitraryFitness<2>, [LowerThanConstraintFunction<S
})
),
],
- -0.0958250414180359
- )
+ bounds: (
+ SVector::<f64, 2>::new(0.0, 0.0), // min bounds
+ SVector::<f64, 2>::new(10.0, 10.0) // max bounds
+ ),
+ optimal_value: -0.0958250414180359,
+ }
}
-pub fn problem_g11(eps: f64) -> (ArbitraryFitness<2>, [LowerThanConstraintFunction<SVector<f64, 2>, f64>; 2], f64) {
- (
- ArbitraryFitness::new(
+pub fn problem_g11(eps: f64) -> ConstrainedProblem<2> {
+ ConstrainedProblem {
+ objective: ArbitraryFitness::new(
Box::new(|vec| {
// Minimize f(x) = x1^2 + (x2 - 1)^2
vec[0].powi(2) + (vec[1] - 1.0).powi(2)
})
),
- [
+ constraints: [
// Equality h(x) = x2 - x1^2 = 0
// Transformed 1: h - eps >= 0 => -(h - eps) <= 0
LowerThanConstraintFunction::new(
@@ 97,19 117,23 @@ pub fn problem_g11(eps: f64) -> (ArbitraryFitness<2>, [LowerThanConstraintFuncti
})
),
],
- 0.7499 // Best known optimum
- )
+ bounds: (
+ SVector::<f64, 2>::new(-1.0, -1.0), // min bounds
+ SVector::<f64, 2>::new(1.0, 1.0) // max bounds
+ ),
+ optimal_value: 0.7499, // Best known optimum
+ }
}
-pub fn problem_g24() -> (ArbitraryFitness<2>, [LowerThanConstraintFunction<SVector<f64, 2>, f64>; 2], f64) {
- (
- ArbitraryFitness::new(
+pub fn problem_g24() -> ConstrainedProblem<2> {
+ ConstrainedProblem {
+ objective: ArbitraryFitness::new(
Box::new(|vec| {
// Minimize f(x) = -x1 - x2
-vec[0] - vec[1]
})
),
- [
+ constraints: [
// g1(x) = -2x1^4 + 8x1^3 - 8x1^2 + x2 - 2 <= 0
LowerThanConstraintFunction::new(
Box::new(|vec| {
@@ 123,10 147,134 @@ pub fn problem_g24() -> (ArbitraryFitness<2>, [LowerThanConstraintFunction<SVect
})
),
],
- -5.50801327159536 // Best known optimum
- )
+ bounds: (
+ SVector::<f64, 2>::new(0.0, 0.0), // min bounds
+ SVector::<f64, 2>::new(3.0, 4.0) // max bounds
+ ),
+ optimal_value: -5.50801327159536, // Best known optimum
+ }
+}
+
+/// Solve a constrained optimization problem using stochastic ranking
+///
+/// Returns the evolution result with feasible fractions for each iteration
+pub fn solve_with_stochastic_ranking<const CONSTRAINTS: usize>(
+ mut problem: ConstrainedProblem<CONSTRAINTS>,
+ population_size: usize,
+ parents_count: usize,
+ iterations: usize,
+ stochastic_params: (usize, f64), // (N, p) for stochastic ranking
+ mutation_std_dev: f64,
+ rng: &mut dyn RngCore,
+) -> Result<(ConstrainedProblem<CONSTRAINTS>, EvolutionResult<SVector<f64, 2>, f64>, Vec<f64>), Box<dyn std::error::Error>> {
+ // Create initial population
+ let initializer = RandomInitializer::new(
+ Box::new(BoundedOVector::new(problem.bounds.0, problem.bounds.1))
+ );
+ let initial_population = Population::from_vec(
+ initializer.initialize(nalgebra::Const::<2>, population_size, rng)
+ );
+
+ // Setup components as specified
+ // let mut selection = TournamentSelection::new(5, 0.9);
+ let mut selection = TournamentSelection::new(5, 0.95);
+ let mut replacement = GenerationalReplacement;
+ let mut pairing = AdjacentPairing::new();
+ let mut crossover = BoundedCrossover::<nalgebra::Const<2>, 2, _>::new(
+ ArithmeticCrossover::new(),
+ problem.bounds.0,
+ problem.bounds.1,
+ BoundedCrossoverStrategy::Trim
+ );
+
+ // Setup bounded random distribution perturbation with Normal distribution
+ let normal_perturbation = RandomDistributionPerturbation::<2, Normal<f64>>::normal(mutation_std_dev)?;
+ let perturbation = BoundedPerturbation::new(
+ normal_perturbation,
+ problem.bounds.0,
+ problem.bounds.1,
+ BoundedPerturbationStrategy::Trim
+ );
+ let mut mutation = MutationPerturbation::new(Box::new(perturbation), 0.1);
+
+ // Constraint weights (initially set to 1.0)
+ let constraint_weights = [1.0; CONSTRAINTS];
+
+ let (N, p) = stochastic_params;
+ let better_than = MinimizingOperator::new();
+
+ // Convert constraint array references
+ let constraint_refs = problem.constraints.iter().collect::<Vec<_>>().try_into()
+ .map_err(|_| "Failed to convert constraint references")?;
+
+ stochastic_ranking_evolution_algorithm(
+ initial_population,
+ parents_count,
+ N,
+ p,
+ &mut problem.objective,
+ constraint_refs,
+ constraint_weights,
+ &mut pairing,
+ &mut selection,
+ &mut crossover,
+ &mut mutation,
+ &mut replacement,
+ &better_than,
+ iterations,
+ rng,
+ ).map(|(x, y)| (problem, x, y))
}
fn main() {
- println!("Hello, world!");
+ let mut rng = rand::rng();
+
+ // Example usage with problem G06
+ let problem = problem_g11(0.01);
+
+ let population_size = 500;
+ let result = solve_with_stochastic_ranking(
+ problem,
+ population_size, // population_size
+ population_size / 2, // parents_count
+ 1000, // iterations
+ (population_size * 2, 0.2), // (N, p) for stochastic ranking
+ 0.001, // mutation_std_dev
+ &mut rng,
+ );
+
+ match result {
+ Ok((problem, evolution_result, feasible_fractions)) => {
+ if let Some(best_candidate) = evolution_result.best_candidate {
+ println!("{:?}", feasible_fractions);
+ println!("Best solution found:");
+ println!(" Chromosome: {:?}", best_candidate.chromosome);
+ println!(" Fitness: {}", best_candidate.evaluation);
+ println!(" Iterations: {}", evolution_result.iterations);
+ println!(" Evaluations: {}", evolution_result.evaluations);
+ println!(" Final feasible fraction: {:.2}%", feasible_fractions.last().unwrap_or(&0.0) * 100.0);
+
+ // Sanity check: verify the best candidate is feasible
+ let best_chromosome = &best_candidate.chromosome;
+
+ println!("\nFeasibility check for best solution:");
+ for (i, constraint) in problem.constraints.iter().enumerate() {
+ match constraint.evaluate(best_chromosome) {
+ Ok(value) => {
+ let is_feasible = value <= 0.0;
+ println!(" Constraint {}: {} ({})", i+1, value, if is_feasible { "FEASIBLE" } else { "INFEASIBLE" });
+ }
+ Err(e) => {
+ println!(" Constraint {}: Error evaluating - {}", i+1, e);
+ }
+ }
+ }
+ } else {
+ println!("Could not find any feasible solution!")
+ }
+ }
+ Err(e) => {
+ println!("Error: {}", e);
+ }
+ }
}