oxigrid 0.1.2

Pure Rust Energy Systems Simulation & Optimization Library
Documentation
/// Battery price arbitrage optimisation.
///
/// Schedules a battery energy storage system (BESS) to buy cheap energy
/// (charge) and sell/use stored energy when prices are high (discharge),
/// maximising economic value subject to SoC, power, and ramp constraints.
///
/// Uses a greedy priority-based dispatch (equivalent to LP solution for
/// simple convex price profiles).
use serde::{Deserialize, Serialize};

/// Battery parameters for arbitrage optimisation.
#[derive(Debug, Clone, Serialize, Deserialize)]
pub struct ArbitrageBattery {
    /// Energy capacity `kWh`
    pub capacity_kwh: f64,
    /// Maximum charge power `kW`
    pub p_charge_max_kw: f64,
    /// Maximum discharge power `kW`
    pub p_discharge_max_kw: f64,
    /// Round-trip efficiency (0–1)
    pub efficiency_rt: f64,
    /// Minimum SoC limit (0–1)
    pub soc_min: f64,
    /// Maximum SoC limit (0–1)
    pub soc_max: f64,
    /// Initial SoC (0–1)
    pub soc_init: f64,
}

impl ArbitrageBattery {
    /// Standard lithium-ion BESS with 90% round-trip efficiency.
    pub fn lithium_ion(capacity_kwh: f64, power_kw: f64) -> Self {
        Self {
            capacity_kwh,
            p_charge_max_kw: power_kw,
            p_discharge_max_kw: power_kw,
            efficiency_rt: 0.90,
            soc_min: 0.10,
            soc_max: 0.90,
            soc_init: 0.50,
        }
    }

    /// Charge efficiency (one-way): √η_rt
    pub fn charge_efficiency(&self) -> f64 {
        self.efficiency_rt.sqrt()
    }

    /// Discharge efficiency (one-way): √η_rt
    pub fn discharge_efficiency(&self) -> f64 {
        self.efficiency_rt.sqrt()
    }
}

/// Schedule for one interval.
#[derive(Debug, Clone, Copy, Serialize, Deserialize)]
pub struct DispatchInterval {
    /// Power `kW`: positive = charge, negative = discharge
    pub power_kw: f64,
    /// SoC at end of interval (0–1)
    pub soc_end: f64,
    /// Revenue/cost for this interval [$/interval] (positive = revenue)
    pub revenue: f64,
}

/// Result of an arbitrage optimisation.
#[derive(Debug, Clone, Serialize, Deserialize)]
pub struct ArbitrageResult {
    /// Per-interval schedule
    pub schedule: Vec<DispatchInterval>,
    /// Total profit [currency units]
    pub total_profit: f64,
    /// Total energy charged `kWh`
    pub energy_charged_kwh: f64,
    /// Total energy discharged `kWh`
    pub energy_discharged_kwh: f64,
    /// Number of equivalent full cycles
    pub cycles: f64,
}

/// Optimise battery dispatch for day-ahead price arbitrage.
///
/// # Arguments
/// - `battery`  — battery specification
/// - `prices`   — electricity prices for each interval [$/kWh]
/// - `dt_h`     — time step duration `hours`
///
/// # Algorithm
/// Greedy: rank intervals by price; charge at cheapest N intervals,
/// discharge at most expensive N intervals, subject to energy balance.
pub fn optimise_arbitrage(
    battery: &ArbitrageBattery,
    prices: &[f64],
    dt_h: f64,
) -> ArbitrageResult {
    let n = prices.len();
    let mut power = vec![0.0_f64; n]; // kW per interval

    // Compute available energy for charging/discharging
    let e_available = (battery.soc_max - battery.soc_min) * battery.capacity_kwh;
    let e_chargeable = e_available;
    let e_dischargeable = e_available;

    // Sort intervals by price for greedy assignment
    let mut idx_sorted: Vec<usize> = (0..n).collect();
    idx_sorted.sort_by(|&a, &b| {
        prices[a]
            .partial_cmp(&prices[b])
            .unwrap_or(std::cmp::Ordering::Equal)
    });

    let mut e_charged = 0.0_f64;
    let mut e_discharged = 0.0_f64;
    let eta_c = battery.charge_efficiency();
    let eta_d = battery.discharge_efficiency();

    // Assign charging to cheapest intervals
    for &i in &idx_sorted {
        let e_step_max = battery.p_charge_max_kw * dt_h; // kWh per step
        let e_remaining = e_chargeable - e_charged;
        if e_remaining <= 1e-6 {
            break;
        }
        let e_this = e_step_max.min(e_remaining);
        power[i] = e_this / dt_h; // kW
        e_charged += e_this;
    }

    // Assign discharging to most expensive intervals (reverse order)
    for &i in idx_sorted.iter().rev() {
        // Don't discharge if we assigned charging here
        if power[i] > 1e-6 {
            continue;
        }
        let e_step_max = battery.p_discharge_max_kw * dt_h;
        let e_remaining = e_dischargeable - e_discharged;
        if e_remaining <= 1e-6 {
            break;
        }
        let e_this = e_step_max.min(e_remaining);
        power[i] = -e_this / dt_h; // negative = discharge
        e_discharged += e_this;
    }

    // Simulate SoC and compute revenue
    let mut soc = battery.soc_init;
    let mut schedule = Vec::with_capacity(n);
    let mut total_profit = 0.0_f64;

    for (i, &p) in power.iter().enumerate() {
        let revenue;
        if p > 0.0 {
            // Charging: pay for electricity, store η_c * E
            let e_input = p * dt_h;
            let e_stored = e_input * eta_c;
            soc = (soc + e_stored / battery.capacity_kwh).min(battery.soc_max);
            revenue = -e_input * prices[i]; // pay
        } else if p < 0.0 {
            // Discharging: deliver η_d * E, earn revenue
            let e_delivered = p.abs() * dt_h;
            let e_from_battery = e_delivered / eta_d;
            soc = (soc - e_from_battery / battery.capacity_kwh).max(battery.soc_min);
            revenue = e_delivered * prices[i]; // earn
        } else {
            revenue = 0.0;
        }
        total_profit += revenue;
        schedule.push(DispatchInterval {
            power_kw: p,
            soc_end: soc,
            revenue,
        });
    }

    ArbitrageResult {
        schedule,
        total_profit,
        energy_charged_kwh: e_charged,
        energy_discharged_kwh: e_discharged,
        cycles: e_discharged / battery.capacity_kwh,
    }
}

/// Compute the spread (max − min price) that makes arbitrage profitable.
///
/// Accounts for round-trip efficiency losses.
/// Minimum viable spread = price_low / η_rt to break even.
pub fn minimum_viable_spread(efficiency_rt: f64) -> f64 {
    // Need: (revenue from discharge) > (cost of charge)
    // p_high * η_rt > p_low → p_high - p_low > p_low * (1/η_rt - 1)
    // Min spread ≈ p_low * (1 - η_rt) / η_rt
    // Normalised (for p_low = 1 $/kWh):
    (1.0 - efficiency_rt) / efficiency_rt
}

#[cfg(test)]
mod tests {
    use super::*;

    fn sample_battery() -> ArbitrageBattery {
        ArbitrageBattery::lithium_ion(100.0, 50.0)
    }

    #[test]
    fn test_arbitrage_charges_at_low_price() {
        let bat = sample_battery();
        // First 12h cheap, last 12h expensive
        let mut prices = vec![0.05_f64; 12]; // cheap
        prices.extend(vec![0.20_f64; 12]); // expensive
        let result = optimise_arbitrage(&bat, &prices, 1.0);
        // Should charge in first half
        let charge_in_cheap: f64 = result.schedule[..12]
            .iter()
            .filter(|d| d.power_kw > 0.0)
            .map(|d| d.power_kw)
            .sum();
        assert!(charge_in_cheap > 0.0, "Should charge during cheap hours");
    }

    #[test]
    fn test_arbitrage_discharges_at_high_price() {
        let bat = sample_battery();
        let mut prices = vec![0.05_f64; 12];
        prices.extend(vec![0.20_f64; 12]);
        let result = optimise_arbitrage(&bat, &prices, 1.0);
        let discharge_in_peak: f64 = result.schedule[12..]
            .iter()
            .filter(|d| d.power_kw < 0.0)
            .map(|d| d.power_kw.abs())
            .sum();
        assert!(
            discharge_in_peak > 0.0,
            "Should discharge during expensive hours"
        );
    }

    #[test]
    fn test_soc_stays_in_bounds() {
        let bat = sample_battery();
        let prices: Vec<f64> = (0..24).map(|h| if h < 12 { 0.05 } else { 0.25 }).collect();
        let result = optimise_arbitrage(&bat, &prices, 1.0);
        for d in &result.schedule {
            assert!(
                d.soc_end >= bat.soc_min - 1e-6,
                "SoC below min: {:.4}",
                d.soc_end
            );
            assert!(
                d.soc_end <= bat.soc_max + 1e-6,
                "SoC above max: {:.4}",
                d.soc_end
            );
        }
    }

    #[test]
    fn test_profit_positive_with_spread() {
        let bat = sample_battery();
        let mut prices = vec![0.05_f64; 12];
        prices.extend(vec![0.20_f64; 12]);
        let result = optimise_arbitrage(&bat, &prices, 1.0);
        assert!(
            result.total_profit > 0.0,
            "Should be profitable: ${:.2}",
            result.total_profit
        );
    }

    #[test]
    fn test_minimum_viable_spread() {
        let spread = minimum_viable_spread(0.90);
        assert!(
            spread > 0.0 && spread < 1.0,
            "Min spread={:.4} $/kWh",
            spread
        );
    }

    #[test]
    fn test_flat_price_no_profit() {
        let bat = sample_battery();
        let prices = vec![0.10_f64; 24];
        let result = optimise_arbitrage(&bat, &prices, 1.0);
        // With flat prices, any cycling loses money due to round-trip losses
        // The optimiser should ideally not cycle, but since we use greedy,
        // just verify cycles are small or profit is ≤ 0
        let net = result.total_profit;
        assert!(
            net <= 1e-6,
            "Flat price should give ≤ 0 profit: ${:.4}",
            net
        );
    }

    #[test]
    fn test_schedule_length_matches_prices() {
        let bat = sample_battery();
        let prices = vec![0.05_f64, 0.10, 0.15, 0.20, 0.25, 0.30];
        let result = optimise_arbitrage(&bat, &prices, 1.0);
        assert_eq!(
            result.schedule.len(),
            prices.len(),
            "Schedule length must equal number of price intervals"
        );
    }

    #[test]
    fn test_charge_efficiency_sqrt() {
        let bat = ArbitrageBattery::lithium_ion(100.0, 50.0);
        let expected = 0.90_f64.sqrt();
        assert!(
            (bat.charge_efficiency() - expected).abs() < 1e-12,
            "charge_efficiency should be sqrt(0.90), got {:.12}",
            bat.charge_efficiency()
        );
    }

    #[test]
    fn test_discharge_efficiency_sqrt() {
        let bat = ArbitrageBattery::lithium_ion(100.0, 50.0);
        let expected = 0.90_f64.sqrt();
        assert!(
            (bat.discharge_efficiency() - expected).abs() < 1e-12,
            "discharge_efficiency should be sqrt(0.90), got {:.12}",
            bat.discharge_efficiency()
        );
    }

    #[test]
    fn test_minimum_viable_spread_perfect_efficiency() {
        let spread = minimum_viable_spread(1.0);
        assert!(
            spread.abs() < 1e-12,
            "With perfect efficiency, min viable spread should be 0.0, got {:.12}",
            spread
        );
    }

    #[test]
    fn test_cycles_nonnegative() {
        let bat = sample_battery();
        let mut prices = vec![0.05_f64; 12];
        prices.extend(vec![0.20_f64; 12]);
        let result = optimise_arbitrage(&bat, &prices, 1.0);
        assert!(
            result.cycles >= 0.0,
            "Cycles must be non-negative, got {:.6}",
            result.cycles
        );
    }

    #[test]
    fn test_empty_prices() {
        let bat = sample_battery();
        let result = optimise_arbitrage(&bat, &[], 1.0);
        assert!(
            result.schedule.is_empty(),
            "Empty price list should yield empty schedule"
        );
    }

    #[test]
    fn test_energy_charged_equals_discharged_symmetric() {
        // Symmetric profile: 12h cheap then 12h expensive, equal durations.
        // The optimizer charges as much as it can and then discharges as much as it can.
        // Both are bounded by the same e_available limit, so they should be equal.
        let bat = sample_battery();
        let mut prices = vec![0.05_f64; 12];
        prices.extend(vec![0.20_f64; 12]);
        let result = optimise_arbitrage(&bat, &prices, 1.0);
        let diff = (result.energy_charged_kwh - result.energy_discharged_kwh).abs();
        // Both are capped at (soc_max - soc_min) * capacity = 0.8 * 100 = 80 kWh
        assert!(
            diff < 1e-6,
            "For symmetric profile, charged ({:.4} kWh) should equal discharged ({:.4} kWh)",
            result.energy_charged_kwh,
            result.energy_discharged_kwh
        );
    }
}