spg_example/
spg_example.rs

1use nalgebra::DVector;
2use optimization_solvers::{
3    BackTracking, FuncEvalMultivariate, LineSearchSolver, SpectralProjectedGradient, Tracer,
4};
5
6fn main() {
7    // Setting up logging
8    std::env::set_var("RUST_LOG", "info");
9    let _ = Tracer::default().with_normal_stdout_layer().build();
10
11    // Convex function: f(x,y) = x^2 + y^2 + exp(x^2 + y^2)
12    // This function is convex and has a minimum at (0, 0)
13    let f_and_g = |x: &DVector<f64>| -> FuncEvalMultivariate {
14        let x1 = x[0];
15        let x2 = x[1];
16
17        // Function value
18        let f = x1.powi(2) + x2.powi(2) + (x1.powi(2) + x2.powi(2)).exp();
19
20        // Gradient
21        let exp_term = (x1.powi(2) + x2.powi(2)).exp();
22        let g1 = 2.0 * x1 * (1.0 + exp_term);
23        let g2 = 2.0 * x2 * (1.0 + exp_term);
24        let g = DVector::from_vec(vec![g1, g2]);
25
26        FuncEvalMultivariate::new(f, g)
27    };
28
29    // Setting up the line search (backtracking)
30    let armijo_factor = 1e-4;
31    let beta = 0.5;
32    let mut ls = BackTracking::new(armijo_factor, beta);
33
34    // Setting up the solver with box constraints
35    let tol = 1e-6;
36    let x0 = DVector::from_vec(vec![0.5, 0.5]); // Starting point
37    let lower_bound = DVector::from_vec(vec![-1.0, -1.0]); // -1 <= x <= 1, -1 <= y <= 1
38    let upper_bound = DVector::from_vec(vec![1.0, 1.0]);
39
40    // Create a mutable oracle for SPG initialization
41    let mut oracle_for_init = f_and_g;
42    let mut solver = SpectralProjectedGradient::new(
43        tol,
44        x0.clone(),
45        &mut oracle_for_init,
46        lower_bound.clone(),
47        upper_bound.clone(),
48    );
49
50    // Running the solver
51    let max_iter_solver = 100;
52    let max_iter_line_search = 20;
53
54    println!("=== Spectral Projected Gradient (SPG) Example ===");
55    println!("Objective: f(x,y) = x^2 + y^2 + exp(x^2 + y^2) (convex)");
56    println!("Global minimum: (0, 0) with f(0,0) = 1");
57    println!("Constraints: -1 <= x <= 1, -1 <= y <= 1");
58    println!("Starting point: {:?}", x0);
59    println!("Lower bounds: {:?}", lower_bound);
60    println!("Upper bounds: {:?}", upper_bound);
61    println!("Tolerance: {}", tol);
62    println!();
63
64    match solver.minimize(
65        &mut ls,
66        f_and_g,
67        max_iter_solver,
68        max_iter_line_search,
69        None,
70    ) {
71        Ok(()) => {
72            let x = solver.x();
73            let eval = f_and_g(x);
74            println!("✅ Optimization completed successfully!");
75            println!("Final iterate: {:?}", x);
76            println!("Function value: {:.6}", eval.f());
77            println!("Gradient norm: {:.6}", eval.g().norm());
78            println!("Iterations: {}", solver.k());
79
80            // Check constraint satisfaction
81            println!("Constraint satisfaction:");
82            for i in 0..x.len() {
83                println!(
84                    "  x[{}] = {:.6} (bounds: [{:.1}, {:.1}])",
85                    i, x[i], lower_bound[i], upper_bound[i]
86                );
87            }
88
89            // Check if we're close to the known minimum
90            let true_min = DVector::from_vec(vec![0.0, 0.0]);
91            let distance_to_min = (x - true_min).norm();
92            println!("Distance to true minimum: {:.6}", distance_to_min);
93            println!("Expected function value: 1.0");
94
95            // Show some properties of SPG
96            println!("SPG properties:");
97            println!("  - Uses spectral step length estimation");
98            println!("  - Handles box constraints efficiently");
99            println!("  - Often faster than standard projected gradient");
100        }
101        Err(e) => {
102            println!("❌ Optimization failed: {:?}", e);
103        }
104    }
105}