~ruther/ctu-fee-eoa

84ad01848b0004b524d31948eb05066bf3e5d9c9 — Rutherther 6 days ago 458589f
feat(constr_hw02): Add nsga_constr algorithm
1 files changed, 254 insertions(+), 5 deletions(-)

M codes/constr_hw02/src/main.rs
M codes/constr_hw02/src/main.rs => codes/constr_hw02/src/main.rs +254 -5
@@ 1,7 1,7 @@
use std::{convert::Infallible, env, fs, io::Write, rc::Rc};

use eoa_lib::{
    bounded::BoundedOVector, comparison::MinimizingOperator, constraints::{stochastic_ranking_evolution_algorithm, ConstrainedEvalFitness, ConstrainedFitnessFunction, ConstraintFunction, LowerThanConstraintFunction}, crossover::{ArithmeticCrossover, BoundedCrossover, BoundedCrossoverStrategy}, evolution::{EvolutionCandidate, EvolutionResult, EvolutionStats}, fitness::FitnessFunction, initializer::{Initializer, RandomInitializer}, multi_objective_evolution::nsga_2, pairing::AdjacentPairing, perturbation::{BoundedPerturbation, BoundedPerturbationStrategy, MutationPerturbation, RandomDistributionPerturbation}, population::{EvaluatedChromosome, EvaluatedPopulation, Population}, replacement::{BestReplacement, GenerationalReplacement}, selection::{BestSelection, TournamentSelection}
    bounded::BoundedOVector, comparison::MinimizingOperator, constraints::{stochastic_ranking_evolution_algorithm, ConstrainedEvalFitness, ConstrainedFitnessFunction, ConstraintFunction, LowerThanConstraintFunction}, crossover::{ArithmeticCrossover, BoundedCrossover, BoundedCrossoverStrategy}, evolution::{EvolutionCandidate, EvolutionResult, EvolutionStats}, fitness::FitnessFunction, initializer::{Initializer, RandomInitializer}, multi_objective_evolution::{nsga_2, constrained_nsga_2}, pairing::AdjacentPairing, perturbation::{BoundedPerturbation, BoundedPerturbationStrategy, MutationPerturbation, RandomDistributionPerturbation}, population::{EvaluatedChromosome, EvaluatedPopulation, Population}, replacement::{BestReplacement, GenerationalReplacement}, selection::{BestSelection, TournamentSelection}
};
use nalgebra::{ArrayStorage, Const, Matrix, SVector};
use rand::RngCore;


@@ 460,11 460,11 @@ pub fn solve_with_stochastic_ranking<const DIM: usize, const CONSTRAINTS: usize>

    // Extract feasible fractions from the result
    let (evolution_result, feasible_fractions) = result;
    

    // For now, create placeholder constraint violations data
    // TODO: This needs library-level changes to properly track constraint violations
    let avg_constraint_violations = vec![0.0; feasible_fractions.len()];
    

    Ok((evolution_result, feasible_fractions, avg_constraint_violations))
}



@@ 710,6 710,108 @@ pub fn solve_with_nsga_multi<const DIM: usize, const CONSTRAINTS: usize, const C
    Ok((result, feasible_fractions, avg_constraint_violations))
}

/// Solve a constrained optimization problem using constrained NSGA-II
pub fn solve_with_nsga_constr<const DIM: usize, const CONSTRAINTS: usize>(
    problem: &ConstrainedProblem<DIM, CONSTRAINTS>,
    population_size: usize,
    parents_count: usize,
    iterations: usize,
    mutation_std_dev: f64,
    rng: &mut dyn RngCore,
) -> Result<(EvolutionResult<SVector<f64, DIM>, [f64; 1]>, Vec<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::<DIM>, population_size, rng)
    );

    // Setup components
    let mut pairing = AdjacentPairing::new();
    let mut crossover = ArithmeticCrossover::new();

    // Setup bounded random distribution perturbation with Normal distribution
    let normal_perturbation = RandomDistributionPerturbation::<DIM, Normal<f64>>::normal(mutation_std_dev)?;
    let mut mutation = MutationPerturbation::new(Box::new(normal_perturbation), 0.1);

    let better_than = MinimizingOperator::new();

    // Clone the problem to access its fields
    let cloned_problem = problem.clone();

    // Create objectives array
    let objectives: [Box<dyn FitnessFunction<In = SVector<f64, DIM>, Out = f64, Err = Infallible>>; 1] = [
        Box::new(cloned_problem.objective),
    ];

    // Convert constraints to boxed trait objects
    let constraints: [Box<dyn ConstraintFunction<Chromosome = SVector<f64, DIM>, Out = f64, Err = Infallible>>; CONSTRAINTS] =
        cloned_problem.constraints.into_iter().map(|c| Box::new(c) as Box<dyn ConstraintFunction<Chromosome = SVector<f64, DIM>, Out = f64, Err = Infallible>>).collect::<Vec<_>>().try_into()
            .map_err(|_| "Failed to convert constraint references")?;

    let mut feasible_fractions = Vec::with_capacity(iterations);
    let mut avg_constraint_violations = Vec::with_capacity(iterations);

    let result = constrained_nsga_2::<1, CONSTRAINTS, SVector<f64, DIM>, f64, Infallible, 2, _, _, _>(
        initial_population,
        parents_count,
        objectives,
        constraints,
        &mut pairing,
        &mut crossover,
        &mut mutation,
        &better_than,
        iterations,
        rng,
        |_iteration: usize, _stats: &EvolutionStats<SVector<f64, DIM>, _>, population: &EvaluatedPopulation<SVector<f64, DIM>, _>| {
            // Calculate feasible fraction and average constraint violation
            let feasible_count = population.population
                .iter()
                .filter(|individual| {
                    // Check if the individual satisfies all constraints
                    problem.constraints.iter().all(|constraint| {
                        constraint.is_feasible(&individual.chromosome).unwrap_or(false)
                    })
                })
                .count();

            // Calculate average constraint violation
            let total_violation: f64 = population.population
                .iter()
                .map(|individual| {
                    problem.constraints
                        .iter()
                        .map(|constraint| constraint.evaluate(&individual.chromosome).unwrap_or(0.0).max(0.0))
                        .sum::<f64>()
                })
                .sum();
            let avg_constraint_violation = total_violation / population.population.len() as f64;

            let feasible_fraction = feasible_count as f64 / population.population.len() as f64;
            feasible_fractions.push(feasible_fraction);
            avg_constraint_violations.push(avg_constraint_violation);
        },
        |chromosome, evaluation, best_candidate| {
            // Only save feasible solutions (all constraints satisfied)
            if !problem.constraints.iter().all(|constraint| {
                constraint.is_feasible(chromosome).unwrap_or(false)
            }) {
                return false;
            }

            if best_candidate.is_none() {
                return true;
            }

            // Compare based on fitness (first objective)
            evaluation.evaluations[0] < best_candidate.as_ref().unwrap().evaluated_chromosome.evaluation.evaluations[0]
        }
    )?;

    Ok((result, feasible_fractions, avg_constraint_violations))
}

const ITERATIONS: usize = 1000;
const POPULATION: usize = 500;
const PARENTS_COUNT: usize = 500;


@@ 1121,6 1223,145 @@ fn handle_g24_nsga() -> Result<(), Box<dyn std::error::Error>> {
    run_nsga_ii(problem, config)
}

fn run_nsga_constr<const DIM: usize, const CONSTRAINTS: usize>(
    problem: ConstrainedProblem<DIM, CONSTRAINTS>,
    config: NsgaConfig,
) -> Result<(), Box<dyn std::error::Error>> {
    let mut rng = rand::rng();

    let result = solve_with_nsga_constr(
        &problem,
        config.population_size,
        config.parents_count,
        config.iterations,
        config.mutation_std_dev,
        &mut rng,
    )?;

    let (evolution_result, feasible_fractions, avg_constraint_violations) = result;

    // Save results to CSV files
    save_evolution_results("nsga_constr", &problem, &evolution_result.clone().map(|x| x[0]), &feasible_fractions, &avg_constraint_violations)?;

    if let Some(best_candidate) = evolution_result.best_candidate {
        println!("Best solution found:");
        println!("  Chromosome: {:?}", best_candidate.chromosome);
        println!("  Objectives: {:?} ({} %)", best_candidate.evaluation, ((best_candidate.evaluation[0] - problem.optimal_value) / problem.optimal_value).abs() * 100.0);
        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!")
    }

    Ok(())
}

// NSGA-Constr handler functions
fn handle_g04_nsga_constr() -> Result<(), Box<dyn std::error::Error>> {
    let problem = problem_g04();
    let config = NsgaConfig {
        population_size: POPULATION,
        parents_count: PARENTS_COUNT,
        iterations: ITERATIONS,
        mutation_std_dev: 1.0,
    };
    run_nsga_constr(problem, config)
}

fn handle_g05_nsga_constr() -> Result<(), Box<dyn std::error::Error>> {
    let problem = problem_g05();
    let config = NsgaConfig {
        population_size: POPULATION,
        parents_count: PARENTS_COUNT,
        iterations: ITERATIONS,
        mutation_std_dev: 10.0,
    };
    run_nsga_constr(problem, config)
}

fn handle_g06_nsga_constr() -> Result<(), Box<dyn std::error::Error>> {
    let problem = problem_g06();
    let config = NsgaConfig {
        population_size: POPULATION,
        parents_count: PARENTS_COUNT,
        iterations: ITERATIONS,
        mutation_std_dev: 1.0,
    };
    run_nsga_constr(problem, config)
}

fn handle_g08_nsga_constr() -> Result<(), Box<dyn std::error::Error>> {
    let problem = problem_g08();
    let config = NsgaConfig {
        population_size: POPULATION,
        parents_count: PARENTS_COUNT,
        iterations: ITERATIONS,
        mutation_std_dev: 0.5,
    };
    run_nsga_constr(problem, config)
}

fn handle_g09_nsga_constr() -> Result<(), Box<dyn std::error::Error>> {
    let problem = problem_g09();
    let config = NsgaConfig {
        population_size: POPULATION,
        parents_count: PARENTS_COUNT,
        iterations: ITERATIONS,
        mutation_std_dev: 1.0,
    };
    run_nsga_constr(problem, config)
}

fn handle_g11_nsga_constr() -> Result<(), Box<dyn std::error::Error>> {
    let problem = problem_g11(G11_EPS);
    let config = NsgaConfig {
        population_size: POPULATION,
        parents_count: PARENTS_COUNT,
        iterations: ITERATIONS,
        mutation_std_dev: 0.01 / 50.0,
    };
    run_nsga_constr(problem, config)
}

fn handle_g21_nsga_constr() -> Result<(), Box<dyn std::error::Error>> {
    let problem = problem_g21();
    let config = NsgaConfig {
        population_size: POPULATION,
        parents_count: PARENTS_COUNT,
        iterations: ITERATIONS,
        mutation_std_dev: 10.0,
    };
    run_nsga_constr(problem, config)
}

fn handle_g24_nsga_constr() -> Result<(), Box<dyn std::error::Error>> {
    let problem = problem_g24();
    let config = NsgaConfig {
        population_size: POPULATION,
        parents_count: PARENTS_COUNT,
        iterations: ITERATIONS,
        mutation_std_dev: 0.1,
    };
    run_nsga_constr(problem, config)
}

// NSGA-Multi handler functions for individual problems
fn handle_nsga_multi_g06(capped: bool) -> Result<(), Box<dyn std::error::Error>> {
    let mut rng = rand::rng();


@@ 1319,7 1560,7 @@ fn main() {

    if args.len() != 3 {
        eprintln!("Usage: {} <method> <problem>", args[0]);
        eprintln!("Methods: srank, nsga, nsga_multi, nsga_multi_noncapped");
        eprintln!("Methods: srank, nsga, nsga_constr, nsga_multi, nsga_multi_noncapped");
        eprintln!("Problems: g04, g05, g06, g08, g09, g11, g21, g24");
        std::process::exit(1);
    }


@@ 1344,6 1585,14 @@ fn main() {
        ("nsga", "g11") => handle_g11_nsga(),
        ("nsga", "g21") => handle_g21_nsga(),
        ("nsga", "g24") => handle_g24_nsga(),
        ("nsga_constr", "g04") => handle_g04_nsga_constr(),
        ("nsga_constr", "g05") => handle_g05_nsga_constr(),
        ("nsga_constr", "g06") => handle_g06_nsga_constr(),
        ("nsga_constr", "g08") => handle_g08_nsga_constr(),
        ("nsga_constr", "g09") => handle_g09_nsga_constr(),
        ("nsga_constr", "g11") => handle_g11_nsga_constr(),
        ("nsga_constr", "g21") => handle_g21_nsga_constr(),
        ("nsga_constr", "g24") => handle_g24_nsga_constr(),
        ("nsga_multi", "g04") => handle_nsga_multi_g04(true),
        ("nsga_multi", "g05") => handle_nsga_multi_g05(true),
        ("nsga_multi", "g06") => handle_nsga_multi_g06(true),


@@ 1362,7 1611,7 @@ fn main() {
        ("nsga_multi_noncapped", "g24") => handle_nsga_multi_g24(false),
        (_, _) => {
            eprintln!("Invalid method '{}' or problem '{}'", method, problem);
            eprintln!("Methods: srank, nsga, nsga_multi, nsga_multi_noncapped");
            eprintln!("Methods: srank, nsga, nsga_constr, nsga_multi, nsga_multi_noncapped");
            eprintln!("Problems: g04, g05, g06, g08, g09, g11, g21, g24");
            std::process::exit(1);
        }