lbfgsbrs 0.1.2

Rust port of L-BFGS-B-C
Documentation

use lbfgsbrs::lbfgsb::{LbfgsbMinimizer, LbfgsbParameters};
use env_logger;
use log;

fn main() {
    env_logger::init();
    
    log::info!("Solving sample problem (Rosenbrock test fcn).");
    log::info!("      (f = 0.0 at the optimal solution.)\n\0");

    // Problem dimensions
    let n = 25;
    
    // Initialize the solution vector
    let mut x: Vec<f64> = vec![3.0; n];
    
    // Define the objective function (Rosenbrock function)
    let f = |x: &Vec<f64>| -> f64 {
        let mut f = 0.25 * (x[0] - 1.0).powi(2);
        
        for i in 1..n {
            let t = x[i] - x[i-1].powi(2);
            f += t.powi(2);
        }
        
        f * 4.0
    };
    
    // Define the gradient function
    let g = |x: &Vec<f64>| -> Vec<f64> {
        let mut grad = vec![0.0; n];
        
        let mut t1 = x[1] - x[0].powi(2);
        grad[0] = 2.0 * (x[0] - 1.0) - 16.0 * x[0] * t1;
        
        for i in 1..n-1 {
            let t2 = t1;
            t1 = x[i+1] - x[i].powi(2);
            grad[i] = 8.0 * t2 - 16.0 * x[i] * t1;
        }
        
        grad[n-1] = 8.0 * t1;
        
        grad
    };
    
    // Create parameters for the optimizer
    let params = LbfgsbParameters {
        m: 5,
        factr: 0.0, 
        pgtol: 0.0, 
        time_limit: 0.2,
        max_iter: 1000,
    };
    
    // Create the optimizer
    let mut optimizer = LbfgsbMinimizer::new(&mut x, &f, &g, Some(params));
    
    // Set bounds for even indices (1-indexed in the original code)
    for i in (0..n).step_by(2) {
        optimizer.set_lower_bound(i, 1.0);
        optimizer.set_upper_bound(i, 100.0);
    }
    
    // Set bounds for odd indices (1-indexed in the original code)
    for i in (1..n).step_by(2) {
        optimizer.set_lower_bound(i, -100.0);
        optimizer.set_upper_bound(i, 100.0);
    }
    
    // Run the optimization
    match optimizer.minimize() {
        Ok(_) => {
            log::trace!("Optimization successful");
        },
        Err(e) => {
            log::warn!("Error during optimization: {:?}", e);
        },
    }
    
    // Get the optimized solution
    let solution = optimizer.get_x();
    
    // Print the solution
    log::info!("Optimization complete.");
    log::info!("Final X = {:?}", solution);
    
    // Add test code to verify the results
    let ref_final_x = vec![
        1.0000000732135395, 1.0000001555788756, 1.0000003157333517, 1.0000006337587277, 
        1.0000012686655855, 1.0000025379074822, 1.0000050761076995, 1.0000101523819898, 
        1.0000203049371699, 1.0000406103219177, 1.0000812223109525, 1.0001624512291802, 
        1.0003249288533613, 1.0006499632878765, 1.0013003490268968, 1.0026023889590963, 
        1.0052115503421177, 1.0104502609362256, 1.0210097298211367, 1.0424608683836587, 
        1.0867246621088944, 1.1809704912345709, 1.3946913011665660, 1.9451638255508485, 
        3.7836623082330907
    ];

    // Compare results with reference solution
    let mut all_match = true;
    for i in 0..25 {
        if (solution[i] - ref_final_x[i]).abs() >= 1e-12 {
            log::info!("Mismatch at index {}: expected {}, got {}", 
                    i, ref_final_x[i], solution[i]);
            all_match = false;
        }
    }
    
    if all_match {
        log::info!("All values match the reference solution within tolerance of 1e-12");
    } else {
        log::info!("Some values don't match the reference solution!");
    }
    
    log::info!("Objective value: {}", f(&solution));
}