drand48 0.2.0

/*
DRAND48 Linear congruential generator
implementation by Radim Kolar <hsn@sendmail.cz> 2025
https://gitlab.com/hsn10/drand48

This is free and unencumbered software released into the public domain.
SPDX-License-Identifier: Unlicense OR CC0-1.0

For more information, please refer to <http://unlicense.org/>
*/

/// Default seed is 0x1234_abcd_330e
pub const SEED: i64 = 0x_1234_abcd_330e;

/// Multiplier a = 0x5DEECE66D (25214903917_i64)
///
/// 25214903917 = 7 × 443 × 739 × 11003
pub const A: i64 = 25214903917;

/// Increment c = 0xB (11_i64)
pub const C: i64 = 11;

/// Modulo m = 2^48 (281_474_976_710_656_i64)
pub const M: i64 = 2i64.pow(48);

/// Period is 2^48 = m (full period LCG)
pub const PERIOD: i64 = M;

/// Minimum returned value 0
pub const MIN: i64 = 0;

/// Maximum returned value is m - 1
pub const MAX: i64 = M - 1;

/// Structure for the 48-bit Linear Congruential Generator (LCG) drand48.
///
/// Generator have 48-bit internal state.
///
/// Parameters:
/// Modulus (m): 2^48
/// Multiplier (a): 0x5DEECE66D (25214903917_i64)
/// Increment (c): 0xB (11_i64)
pub struct DRAND48(i64);

/// Creates a generator seeded by srand48() POSIX.1 method
///
/// 0x_330e is used for lowest seed bits.
/// seed can be negative.
pub fn srand48(seed: i32) -> DRAND48 {
    DRAND48::seed( ( (seed as i64) << 16 | 0x330e ) & 0xFFFF_FFFF_FFFF_i64 )
}

impl DRAND48 {
    /// Creates a new instance of the drand48 generator with the default seed.
    ///
    /// The default seed for drand48 is 0x1234_abcd_330e_i64.
    pub fn new() -> Self {
        Self(SEED)
    }

    /// Creates a new instance of the drand48 generator with a custom seed.
    ///
    /// Function will panic if called with a negative or out of range seed.
    /// Maximum seed value is 2^48-1.
    /// ### Note:
    /// Use [`validate_seed`] or [`clamp_seed`] to check or sanitize input.
    pub fn seed(seed: i64) -> Self {
        if !validate_seed(seed) {
           panic!("Invalid seed. Must be non negative and lower than 2^48.");
        }
        Self(seed & 0xFFFF_FFFF_FFFF_i64) // Trimmed to 48 bits
    }

    /// Generates the next 48-bit pseudorandom integer in the sequence.
    ///
    /// The internal state is updated according to the LCG formula:
    /// X_n+1 = (a * X_n + c) mod 2^48
    ///
    /// The 48 lowest bits of the previous state are used for the calculation.
    ///
    /// # Return value
    /// `i64` representing the 48-bit pseudorandom integer.
    ///
    /// ### Note:
    /// This method returns the internal 48-bit state.
    /// To obtain a `double` in the range [0.0, 1.0), use the [`drand48()`] method.
    pub fn next(&mut self) -> i64 {
        // Apply the LCG formula: X_n+1 = (a * X_n + c) mod 2^48
        // The modulo 2^48 operation is implicitly performed by ensuring that
        // our state doesn't exceed 48 bits (we do this via bitmasking).
        // But in this case, it's simpler to just perform the calculation on u64,
        // and then trim it to 48 bits because the multiplication can exceed 48 bits,
        // but the result mod 2^48 is still valid.

        // Calculate the new state
        // In Rust, to achieve the modulo 2^48 effect, we just ensure that
        // the result of the computation doesn’t exceed 48 bits.
        // Since both the multiplier and increment are 48-bit (or smaller),
        // and the state is 48-bit, the result fits into a u64.
        // The `wrapping_mul` and `wrapping_add` operations in Rust automatically
        // handle the "modulo 2^64" behavior, which is fine for us if we then
        // use only the lowest 48 bits.
        // drand48 internally uses exactly these 48 bits.
        self.0 = (A.wrapping_mul(self.0).wrapping_add(C)) & M-1;
        self.0

        // Return the current state (before converting it to double or a 31-bit number),
        // because this is what "next" typically returns in an internal LCG implementation.
        // From this, various types of output (double, long, etc.) are generated.
    }

    /// Generates a pseudorandom number of type double in the range [0.0, 1.0).
    ///
    /// This method internally calls `next()` and converts the 48-bit integer
    /// to a float by discarding lower 24 bits.
    ///
    /// Result is normalized to the interval [0.0, 1.0).
    pub fn frand48(&mut self) -> f32 {
        let next_val = self.next();

        // We want top 24 bits from 48-bit output.
        let top_24_bits = next_val >> 24;
        // Divide by 2^24, to normalize into [0.0, 1.0).
        // (1u32 << 24) is 2^24.
        (top_24_bits as f32) / (1u32 << 24) as f32
    }

    /// Generates a pseudorandom number of type double in the range [0.0, 1.0).
    ///
    /// This method internally calls `next()` and converts the 48-bit integer
    /// to a double, normalized to the interval [0.0, 1.0).
    pub fn drand48(&mut self) -> f64 {
        let next_val = self.next();
        // Divide by 2^48 (which is 0x1_0000_0000_0000_f64)
        // or simply 2.0_f64.powi(48)
        (next_val as f64) / M as f64
    }

    /// Generates a pseudorandom integer of type i32 in the range [0, 2^31 - 1].
    ///
    /// This corresponds to the behavior of the standard `lrand48()` function,
    /// which returns the highest 31 bits (or the 31 most significant bits from the 48-bit state).
    pub fn lrand48(&mut self) -> i32 {
        // From drand48, the 31 most significant bits are taken,
        // which corresponds to bits 16 to 46 (counting from 0).
        // Or it's easier to take the 48-bit state and shift it right by 17 bits.
        // This moves bits 17-47 into positions 0-30, which corresponds to i32.
        ((self.next() >> 17) & 0x7FFF_FFFF_i64) as i32
    }

/// Generates a pseudorandom signed 32-bit integer (`i32`) in the range [-2^31, 2^31 - 1].
///
/// This corresponds to the behavior of the standard POSIX `mrand48()` function,
/// which returns a signed 32-bit integer derived from the top 32 bits of the internal 48-bit state.
///
/// The function shifts the 48-bit internal state right by 16 bits, yielding a 32-bit value.
/// This includes both the sign bit and 31 data bits, matching the expected behavior of `mrand48()`.
///
/// # Return value
/// A signed 32-bit integer in the full `i32` range.
pub fn mrand48(&mut self) -> i32 {
    // Shift the 48-bit state right by 16 bits to get the upper 32 bits
    // This matches the expected output of mrand48 (signed 32-bit integer)
    ((self.next() >> 16) & 0xFFFF_FFFF_i64) as i32
}
}

/**
  Check if seed is valid for drand48 generator.

  Seed must be:
  1. non negative. Can be zero.
  1. less than 2^48
*/
pub fn validate_seed(seed_val: i64) -> bool {
    if seed_val < 0 {
        false
    } else {
       seed_val < M
    }
}

/**
  clamp seed into valid drand48 range
*/
pub fn clamp_seed(seed_val: i64) -> i64 {
   seed_val.abs() % M
}

/** extracting values from DRAND48.next() output */
pub mod extract;

#[cfg(test)]
#[path = "const_test.rs"]
mod consts;

#[cfg(test)]
#[path = "mrand48_test.rs"]
mod mrand;

#[cfg(test)]
#[path = "lrand48_test.rs"]
mod lrand;

#[cfg(test)]
#[path = "drand48_test.rs"]
mod drand;

#[cfg(test)]
#[path = "validate_test.rs"]
mod validate;

#[cfg(test)]
#[path = "clamp_test.rs"]
mod clamp;