librebound-sys 4.6.0

//! Safe RAII wrappers around REBOUND C objects.

use crate::ffi;
use crate::{Error, Result};

// ---------------------------------------------------------------------------
// IAS15 adaptive-timestep mode (mirrors `enum reb_ias15.adaptive_mode`)
// ---------------------------------------------------------------------------

/// IAS15 timestep-adaptation rule. REBOUND default since 2024-01 is
/// [`Self::Prs23`]; older code used [`Self::Global`].
#[derive(Debug, Clone, Copy, PartialEq, Eq)]
#[repr(i32)]
pub enum Ias15AdaptiveMode {
    /// Per-particle fractional error.
    Individual = 0,
    /// Single global fractional-error estimate (default before 2024-01).
    Global = 1,
    /// Pham, Rein & Spiegel (2023) criterion (default since 2024-01).
    Prs23 = 2,
    /// Aarseth (1985) timestep criterion.
    Aarseth85 = 3,
}

impl Ias15AdaptiveMode {
    fn from_raw(raw: i32) -> Self {
        match raw {
            0 => Self::Individual,
            1 => Self::Global,
            2 => Self::Prs23,
            3 => Self::Aarseth85,
            // Forward-compatibility: unknown values reported as the current
            // REBOUND default. The setter only accepts the four enum values,
            // so this branch is reachable only if REBOUND adds a new mode.
            _ => Self::Prs23,
        }
    }
}

// ---------------------------------------------------------------------------
// IntegratorConfig
// ---------------------------------------------------------------------------

/// Per-call IAS15 integrator knobs. `None` for any field leaves the REBOUND
/// default in place.
///
/// Applied to a freshly-created [`Simulation`] before any particles are added.
#[derive(Debug, Clone, Copy, Default)]
pub struct IntegratorConfig {
    /// REBOUND `r->dt` (initial timestep). The integrator picks the sign
    /// automatically once the first step is taken.
    pub initial_dt: Option<f64>,
    /// REBOUND `r->ri_ias15.min_dt`. When the adaptive step would shrink
    /// below this, it clamps. `None` (or 0) = no floor.
    pub min_dt: Option<f64>,
    /// REBOUND `r->ri_ias15.epsilon` (precision). REBOUND default 1e-9.
    pub epsilon: Option<f64>,
    /// REBOUND `r->ri_ias15.adaptive_mode`. Default `Prs23` since 2024-01.
    pub adaptive_mode: Option<Ias15AdaptiveMode>,
}

impl IntegratorConfig {
    /// Apply each `Some` field to `sim`. Should be called immediately after
    /// `Simulation::new()` and before adding particles.
    pub fn apply(&self, sim: &mut Simulation) {
        if let Some(dt) = self.initial_dt {
            sim.set_dt(dt);
        }
        if let Some(eps) = self.epsilon {
            sim.set_ias15_epsilon(eps);
        }
        if let Some(min_dt) = self.min_dt {
            sim.set_ias15_min_dt(min_dt);
        }
        if let Some(mode) = self.adaptive_mode {
            sim.set_ias15_adaptive_mode(mode);
        }
    }
}

// ---------------------------------------------------------------------------
// Simulation
// ---------------------------------------------------------------------------

/// Owned REBOUND simulation. Freed on drop.
///
/// Not `Send`/`Sync` — each thread must create its own simulation.
pub struct Simulation {
    ptr: *mut ffi::reb_simulation,
}

impl Simulation {
    /// Create a new, empty REBOUND simulation.
    pub fn new() -> Result<Self> {
        let ptr = unsafe { ffi::reb_simulation_create() };
        if ptr.is_null() {
            return Err(Error::Other("reb_simulation_create returned null".into()));
        }
        Ok(Self { ptr })
    }

    /// Raw const pointer to the underlying REBOUND simulation. Downstream
    /// `-sys` crates layered on top of REBOUND use this to call ASSIST /
    /// other C APIs that take `*const reb_simulation`.
    pub fn as_ptr(&self) -> *const ffi::reb_simulation {
        self.ptr
    }

    /// Raw mutable pointer to the underlying REBOUND simulation. Required by
    /// C APIs that take `*mut reb_simulation` (assist_attach, assist_detach,
    /// etc.) and by benchmark probes.
    pub fn as_mut_ptr(&mut self) -> *mut ffi::reb_simulation {
        self.ptr
    }

    /// Alias for [`as_mut_ptr`] kept for backward compatibility.
    #[doc(hidden)]
    pub fn raw_ptr_mut(&mut self) -> *mut ffi::reb_simulation {
        self.ptr
    }

    pub fn t(&self) -> f64 {
        unsafe { ffi::assist_rs_sim_get_t(self.ptr) }
    }
    pub fn set_t(&mut self, t: f64) {
        unsafe { ffi::assist_rs_sim_set_t(self.ptr, t) }
    }

    pub fn dt(&self) -> f64 {
        unsafe { ffi::assist_rs_sim_get_dt(self.ptr) }
    }
    pub fn set_dt(&mut self, dt: f64) {
        unsafe { ffi::assist_rs_sim_set_dt(self.ptr, dt) }
    }

    pub fn n_particles(&self) -> usize {
        unsafe { ffi::assist_rs_sim_get_N(self.ptr) as usize }
    }

    /// Total IAS15 steps (accepted + rejected) since this simulation was
    /// created. Useful for diagnosing adaptive-timestep behavior.
    pub fn steps_done(&self) -> u64 {
        unsafe { ffi::assist_rs_sim_get_steps_done(self.ptr) }
    }

    pub fn n_var(&self) -> i32 {
        unsafe { ffi::assist_rs_sim_get_N_var(self.ptr) }
    }

    pub fn status(&self) -> i32 {
        unsafe { ffi::assist_rs_sim_get_status(self.ptr) }
    }

    pub fn set_exact_finish_time(&mut self, v: bool) {
        unsafe { ffi::assist_rs_sim_set_exact_finish_time(self.ptr, v as i32) }
    }

    /// IAS15 precision parameter (REBOUND `r->ri_ias15.epsilon`). Default 1e-9.
    /// Larger values are looser but faster.
    pub fn ias15_epsilon(&self) -> f64 {
        unsafe { ffi::assist_rs_sim_get_ias15_epsilon(self.ptr) }
    }
    pub fn set_ias15_epsilon(&mut self, eps: f64) {
        unsafe { ffi::assist_rs_sim_set_ias15_epsilon(self.ptr, eps) }
    }

    /// IAS15 minimum timestep floor (REBOUND `r->ri_ias15.min_dt`). When the
    /// adaptive step would shrink below this, it clamps instead of grinding.
    /// Default 0 = no floor.
    pub fn ias15_min_dt(&self) -> f64 {
        unsafe { ffi::assist_rs_sim_get_ias15_min_dt(self.ptr) }
    }
    pub fn set_ias15_min_dt(&mut self, min_dt: f64) {
        unsafe { ffi::assist_rs_sim_set_ias15_min_dt(self.ptr, min_dt) }
    }

    /// IAS15 adaptive-timestep selector. Default `Prs23` (REBOUND default
    /// since 2024-01).
    pub fn ias15_adaptive_mode(&self) -> Ias15AdaptiveMode {
        let raw = unsafe { ffi::assist_rs_sim_get_ias15_adaptive_mode(self.ptr) };
        Ias15AdaptiveMode::from_raw(raw)
    }
    pub fn set_ias15_adaptive_mode(&mut self, mode: Ias15AdaptiveMode) {
        unsafe { ffi::assist_rs_sim_set_ias15_adaptive_mode(self.ptr, mode as i32) }
    }

    /// Diagnostic counter: how many IAS15 steps hit the predictor-corrector
    /// iteration cap without converging. Monotone-increasing across the
    /// simulation lifetime; nonzero indicates the integrator was working
    /// at the edge of convergence (typically a hint to tighten `epsilon` or
    /// raise `min_dt`).
    pub fn ias15_iterations_max_exceeded(&self) -> u64 {
        unsafe { ffi::assist_rs_sim_get_ias15_iterations_max_exceeded(self.ptr) }
    }

    /// Add a particle to the simulation.
    pub fn add_particle(&mut self, p: ffi::reb_particle) {
        unsafe { ffi::reb_simulation_add(self.ptr, p) }
    }

    /// Add a test particle with given position and velocity (mass = 0).
    pub fn add_test_particle(&mut self, x: f64, y: f64, z: f64, vx: f64, vy: f64, vz: f64) {
        let p = ffi::reb_particle {
            x,
            y,
            z,
            vx,
            vy,
            vz,
            m: 0.0,
            ..Default::default()
        };
        self.add_particle(p);
    }

    /// Read-only access to the particle array.
    pub fn particles(&self) -> &[ffi::reb_particle] {
        let n = self.n_particles();
        if n == 0 {
            return &[];
        }
        let ptr = unsafe { ffi::assist_rs_sim_get_particles(self.ptr) };
        if ptr.is_null() {
            return &[];
        }
        unsafe { std::slice::from_raw_parts(ptr, n) }
    }

    /// Integrate to target time. Returns the status code.
    pub fn integrate(&mut self, tmax: f64) -> Result<()> {
        let status = unsafe { ffi::reb_simulation_integrate(self.ptr, tmax) };
        match status {
            ffi::REB_STATUS_SUCCESS | ffi::REB_STATUS_RUNNING => Ok(()),
            ffi::REB_STATUS_NO_PARTICLES => Err(Error::NoParticles),
            ffi::REB_STATUS_ENCOUNTER => Err(Error::CloseEncounter),
            ffi::REB_STATUS_ESCAPE => Err(Error::Escape),
            ffi::REB_STATUS_COLLISION => Err(Error::Collision),
            other => Err(Error::IntegrationFailed(other)),
        }
    }

    /// Add first-order variational particles for a test particle.
    /// Returns the index of the first variational particle.
    pub fn add_variation_1st_order(&mut self, testparticle: i32) -> i32 {
        unsafe { ffi::reb_simulation_add_variation_1st_order(self.ptr, testparticle) }
    }
}

impl Drop for Simulation {
    fn drop(&mut self) {
        if !self.ptr.is_null() {
            unsafe { ffi::reb_simulation_free(self.ptr) }
        }
    }
}