#![no_std]
#![deny(unsafe_code)]
#![warn(missing_docs)]
#![deny(missing_debug_implementations)]
#![cfg_attr(docs_rs, feature(doc_cfg))]

//! Simple and minimalist randomization library.
//!
//! NOT FOR CRYPTOGRAPHIC PURPOSES.
//!
//! ## Usage
//!
//! You should make a [`PCG32`] value, and then generally you'll call methods
//! from the [`Gen32`] trait:
//!
//! ```rust
//! use randomize::{Gen32, PCG32};
//! let mut g = PCG32::seed(5, 6);
//! println!("rolling 1d6: {}", g.dice(1, 6));
//! ```
//!
//! ## Cargo Features
//! * `getrandom`: adds the [`from_getrandom`](PCG32::from_getrandom) method to
//!   the [`PCG32`] type, which makes a new generator from OS randomness. This
//!   depends on the [getrandom](https://docs.rs/getrandom) crate.

use core::convert::{TryFrom, TryInto};

/// A default seed for any PCG. Downcast to fit, as necessary.
const DEFAULT_PCG_SEED: u128 = 201526561274146932589719779721328219291;

/// A default `inc` for any PCG. Downcast to fit, as necessary.
const DEFAULT_PCG_INC: u128 = 34172814569070222299;

/// The PCG multiplier for 64 bits of state.
const PCG_MULTIPLIER_64: u64 = 6364136223846793005;

/// For random number generators with a primary output of `u32` per use.
///
/// All other methods then default in terms of the `next_u32` method.
pub trait Gen32 {
  /// Produce a `u32`
  fn next_u32(&mut self) -> u32;

  /// Produce a `bool`
  #[inline(always)]
  fn next_bool(&mut self) -> bool {
    (self.next_u32() as i32) < 0
  }

  /// Produce a `u8`
  #[inline(always)]
  fn next_u8(&mut self) -> u8 {
    (self.next_u32() >> 24) as u8
  }

  /// Produce a `u16`
  #[inline(always)]
  fn next_u16(&mut self) -> u16 {
    (self.next_u32() >> 16) as u16
  }

  /// Produce a `u64`
  #[inline(always)]
  fn next_u64(&mut self) -> u64 {
    let l = self.next_u32() as u64;
    let h = self.next_u32() as u64;
    h << 32 | l
  }

  /// Returns an `f32` in the unsigned unit range, `[0, 1]`
  #[inline]
  fn next_f32_unit(&mut self) -> f32 {
    ieee754_random_f32(self, true)
  }

  /// Returns an `f32` in the signed unit range, `[-1, 1]`
  #[inline]
  fn next_f32_signed_unit(&mut self) -> f32 {
    ieee754_random_f32(self, false)
  }

  /// Gives a value within `0 .. B`
  ///
  /// This is often more efficient than making a [`BoundedRandU32`] if you don't
  /// need to use a specific bound value more than once.
  ///
  /// ## Panics
  /// * If the input is 0.
  #[inline]
  fn next_bounded(&mut self, b: u32) -> u32 {
    assert!(b != 0, "Gen32::next_bounded> Bound must be non-zero.");
    let mut x = self.next_u32() as u64;
    let mut mul = (b as u64).wrapping_mul(x);
    let mut low = mul as u32;
    if low < b {
      let threshold = b.wrapping_neg() % b;
      while low < threshold {
        x = self.next_u32() as u64;
        mul = (b as u64).wrapping_mul(x);
        low = mul as u32;
      }
    }
    let high = (mul >> 32) as u32;
    high
  }

  /// Performs an `XdY` style dice roll.
  ///
  /// * If `count` or `sides` are less than 0, the output is 0.
  /// * Requires linear time to compute based on `count`. Expected inputs are 20
  ///   or less.
  #[inline]
  fn dice(&mut self, mut count: i32, sides: i32) -> i32 {
    use core::cmp::Ordering;
    let range = match sides.cmp(&1) {
      Ordering::Less => return 0,
      Ordering::Equal => return count.max(0),
      Ordering::Greater => match sides {
        4 => D4,
        6 => D6,
        8 => D8,
        10 => D10,
        12 => D12,
        20 => D20,
        _ => StandardDie::new(sides as u32),
      },
    };
    let mut t = 0_i32;
    while count > 0 {
      t = t.wrapping_add(range.sample(self));
      count -= 1;
    }
    t
  }

  /// Performs a "step" roll according to the 4e chart.
  ///
  /// This relates to a particular paper and pencil RPG. If you're not familiar
  /// with the game that's fine.
  /// * The average output of any positive value is approximately equal to the
  ///   input, with no hard upper bound.
  /// * The output of any non-positive value is 1.
  /// * Requires linear time to compute. Expected inputs are 30 or less.
  #[inline]
  fn step_ed4(&mut self, mut step: i32) -> i32 {
    if step < 1 {
      1
    } else {
      let mut total: i32 = 0;
      while step > 13 {
        total = total.wrapping_add(X12.sample(self));
        step -= 7;
      }
      total.wrapping_add(match step {
        1 => X4.sample(self).wrapping_sub(2).max(1),
        2 => X4.sample(self).wrapping_sub(1).max(1),
        3 => X4.sample(self),
        4 => X6.sample(self),
        5 => X8.sample(self),
        6 => X10.sample(self),
        7 => X12.sample(self),
        8 => X6.sample(self).wrapping_add(X6.sample(self)),
        9 => X8.sample(self).wrapping_add(X6.sample(self)),
        10 => X8.sample(self).wrapping_add(X8.sample(self)),
        11 => X10.sample(self).wrapping_add(X8.sample(self)),
        12 => X10.sample(self).wrapping_add(X10.sample(self)),
        13 => X12.sample(self).wrapping_add(X10.sample(self)),
        _ => unreachable!(),
      })
    }
  }

  /// Rolls an After Sundown style dice pool.
  ///
  /// This relates to a particular paper and pencil RPG. If you're not familiar
  /// with the game that's fine.
  /// * `size` D6s are rolled. This returns the number of them that are a 5 or
  ///   6.
  #[inline]
  fn sundown_pool(&mut self, mut size: u32) -> u32 {
    let mut hits = 0;
    while size > 0 {
      if D6.sample(self) >= 5 {
        hits += 1
      }
      size -= 1;
    }
    hits
  }

  /// Returns a value in `0..x` with the odds modified by `luck`.
  ///
  /// This pertains to a particular video game. If you're not familiar with the
  /// game that's fine.
  /// * This is a constant time operation.
  /// * higher luck pushes the output toward zero.
  /// * lower luck pushes the output towards the upper value.
  /// * `luck` is expected to be +/-30
  ///
  /// ## Panics
  /// * If `x` is 0 or less.
  #[inline]
  fn rn_bounded_luck(&mut self, x: i32, luck: i32) -> i32 {
    assert!(x > 0);
    let adjustment =
      if x <= 15 { (luck.abs() + 1) / 3 * luck.signum() } else { luck };
    let mut i = self.next_bounded(x as u32) as i32;
    if adjustment != 0 && self.next_bounded(37 + adjustment.abs() as u32) != 0 {
      i -= adjustment;
      i = i.max(0).min(x - 1);
    }
    i
  }

  /// Returns a value of 1 or more.
  ///
  /// * The output starts at 1, then has a repeated `1/x` chance of getting +1.
  /// * As soon as the value doesn't get a +1, it is returned.
  ///
  /// ## Panics
  /// * If `x` is less than 2.
  #[inline]
  fn rn_exponential_decay(&mut self, x: i32) -> i32 {
    assert!(x > 1);
    let mut temp = 1;
    while self.next_bounded(x as u32) == 0 {
      temp += 1;
    }
    temp
  }
  
  /// Returns a value.
  ///
  /// This pertains to a particular video game. If you're not familiar with the
  /// game that's fine.
  /// * The input value affects the output.
  /// * This runs in constant time.
  #[inline]
  fn rn_z(&mut self, i: i32) -> i32 {
    let mut x = i as i64;
    let mut temp = 1000_i64;
    temp += self.next_bounded(1000) as i64;
    temp *= self.rn_exponential_decay(4).min(5) as i64;
    if self.next_bool() {
      x *= temp;
      x /= 1000;
    } else {
      x *= 1000;
      x /= temp;
    }
    x as i32
  }

  /// Gets a value out of the slice given (by copy).
  ///
  /// * The default impl will not pick past index `u32::MAX`.
  #[inline(always)]
  fn pick<T>(&mut self, buf: &[T]) -> T
  where
    Self: Sized,
    T: Copy,
  {
    let end: u32 = saturating_usize_as_u32(buf.len());
    buf[usize::try_from(self.next_bounded(end)).unwrap()]
  }

  /// Gets a value out of the slice given (by shared ref).
  ///
  /// * The default impl will not pick past index `u32::MAX`.
  #[inline(always)]
  fn pick_ref<'b, T>(&mut self, buf: &'b [T]) -> &'b T
  where
    Self: Sized,
  {
    let end: u32 = saturating_usize_as_u32(buf.len());
    &buf[usize::try_from(self.next_bounded(end)).unwrap()]
  }

  /// Gets a value out of the slice given (by unique ref).
  ///
  /// * The default impl will not pick past index `u32::MAX`.
  #[inline(always)]
  fn pick_mut<'b, T>(&mut self, buf: &'b mut [T]) -> &'b mut T
  where
    Self: Sized,
  {
    let end: u32 = saturating_usize_as_u32(buf.len());
    &mut buf[usize::try_from(self.next_bounded(end)).unwrap()]
  }

  /// Shuffles a slice in `O(len)` time.
  ///
  /// * The default impl shuffles only the first `u32::MAX` elements.
  #[inline]
  fn shuffle<T>(&mut self, buf: &mut [T])
  where
    Self: Sized,
  {
    // Note(Lokathor): The "standard" Fisher-Yates shuffle goes backward from
    // the end of the slice, but this version allows us to access memory forwrd
    // from the start to the end, so that we play more nicely with the
    // fetch-ahead of most modern CPUs.
    let mut possibility_count: u32 =
      buf.len().try_into().unwrap_or(u32::max_value());
    let mut this_index: usize = 0;
    let end = buf.len() - 1;
    while this_index < end {
      let offset = self.next_bounded(possibility_count) as usize;
      buf.swap(this_index, this_index + offset);
      possibility_count -= 1;
      this_index += 1;
    }
  }
}

// Asserts that `Gen32` is an object-safe trait.
const _: [&mut dyn Gen32; 0] = [];

/// A [permuted congruential generator](https://en.wikipedia.org/wiki/Permuted_congruential_generator)
/// with 32 bits of output per step.
#[derive(Debug, Clone, PartialEq, Eq, Hash)]
#[repr(C)]
pub struct PCG32 {
  state: u64,
  inc: u64,
}
impl Default for PCG32 {
  #[inline(always)]
  fn default() -> Self {
    PCG32::seed(DEFAULT_PCG_SEED as u64, DEFAULT_PCG_INC as u64)
  }
}
impl PCG32 {
  /// Seed a new generator.
  #[inline]
  pub const fn seed(seed: u64, inc: u64) -> Self {
    let inc = (inc << 1) | 1;
    let mut state = pcg_core_state64(0, inc);
    state = state.wrapping_add(seed);
    state = pcg_core_state64(state, inc);
    Self { state, inc }
  }

  /// Generates a new PCG32 by calling [`getrandom`](getrandom::getrandom).
  ///
  /// If `getrandom` fails then you get the default generator.
  #[inline]
  #[cfg(feature = "getrandom")]
  #[cfg_attr(docs_rs, doc(cfg(feature = "getrandom")))]
  pub fn from_getrandom() -> Self {
    const SIZE_U64: usize = core::mem::size_of::<u64>();
    let mut buf = [0_u8; SIZE_U64 * 2];
    match getrandom::getrandom(&mut buf) {
      Ok(_) => {
        let (seed_slice, inc_slice) = buf.split_at(SIZE_U64);
        let seed = u64::from_ne_bytes(seed_slice.try_into().unwrap());
        let inc = u64::from_ne_bytes(inc_slice.try_into().unwrap());
        Self::seed(seed, inc)
      }
      Err(_) => Self::default(),
    }
  }

  /// Runs the generator once and gets a `u32` as output.
  #[inline]
  pub fn next_u32(&mut self) -> u32 {
    let out = xsh_rr_64_32(self.state);
    self.state = pcg_core_state64(self.state, self.inc);
    out
  }

  /// Advances the generator `delta` steps in `log(delta)` time.
  #[inline]
  pub fn jump(&mut self, delta: u64) {
    self.state = jump_lcg64(delta, self.state, PCG_MULTIPLIER_64, self.inc)
  }
}
impl From<[u64; 2]> for PCG32 {
  #[must_use]
  #[inline(always)]
  fn from([state, inc]: [u64; 2]) -> Self {
    Self { state, inc }
  }
}
impl Gen32 for PCG32 {
  #[inline(always)]
  fn next_u32(&mut self) -> u32 {
    PCG32::next_u32(self)
  }
  #[inline(always)]
  fn next_u64(&mut self) -> u64 {
    let out = xsl_rr_rr_64_64(self.state);
    self.state = pcg_core_state64(self.state, self.inc);
    out
  }
}

/// Permutation: XSH RR `u64` to `u32`.
#[must_use]
#[inline(always)]
const fn xsh_rr_64_32(state: u64) -> u32 {
  ((((state >> 18) ^ state) >> 27) as u32).rotate_right((state >> 59) as u32)
}

/// Permutation: XSL RR RR `u64`
#[must_use]
#[inline(always)]
const fn xsl_rr_rr_64_64(state: u64) -> u64 {
  let rot1: u32 = (state >> 59) as u32;
  let high: u32 = (state >> 32) as u32;
  let low: u32 = state as u32;
  let xor_d: u32 = high ^ low;
  let new_low: u32 = xor_d.rotate_right(rot1);
  let new_high: u32 = high.rotate_right(new_low & 31);
  ((new_high as u64) << 32) | new_low as u64
}

/// Advances a PCG with 64 bits of state.
#[must_use]
#[inline(always)]
const fn pcg_core_state64(state: u64, inc: u64) -> u64 {
  lcg64(state, PCG_MULTIPLIER_64, inc)
}

/// The `u64` LCG, not suitably random on its own.
#[must_use]
#[inline(always)]
const fn lcg64(state: u64, mult: u64, inc: u64) -> u64 {
  state.wrapping_mul(mult).wrapping_add(inc)
}

macro_rules! make_jump_lcgX {
  ($(#[$attr:meta])* $f:ident, $u:ty) => {
    $(#[$attr])*
    /// Gives the state `delta` steps from now in `log(delta)` time.
    #[must_use]
    #[inline(always)]
    const fn $f(mut delta: $u, state: $u, mult: $u, inc: $u) -> $u {
      let mut cur_mult: $u = mult;
      let mut cur_plus: $u = inc;
      let mut acc_mult: $u = 1;
      let mut acc_plus: $u = 0;
      while delta > 0 {
        if (delta & 1) > 0 {
          acc_mult = acc_mult.wrapping_mul(cur_mult);
          acc_plus = acc_plus.wrapping_mul(cur_mult).wrapping_add(cur_plus);
        }
        cur_plus = cur_mult.wrapping_add(1).wrapping_mul(cur_plus);
        cur_mult = cur_mult.wrapping_mul(cur_mult);
        delta /= 2;
      }
      acc_mult.wrapping_mul(state).wrapping_add(acc_plus)
    }
  };
}
make_jump_lcgX!(jump_lcg64, u64);

/// Converts the `usize` into a `u32`, or gives `u32::MAX` if that wouldn't fit.
#[inline(always)]
const fn saturating_usize_as_u32(val: usize) -> u32 {
  #[cfg(target_pointer_width = "16")]
  {
    val as u32
  }
  #[cfg(target_pointer_width = "32")]
  {
    val as u32
  }
  #[cfg(target_pointer_width = "64")]
  {
    if val <= core::u32::MAX as usize {
      val as u32
    } else {
      core::u32::MAX
    }
  }
}

/// Stores the values to sample a number in `0 .. N`
///
/// Making one of these performs a division operation. In comparison,
/// [`Gen32::next_bounded`] will avoid needing to do a division much of the
/// time. Thus, unless you need to sample repeatedly from a specific bounded
/// range, simply calling `next_bounded` directly might be more efficient.
#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
pub struct BoundedRandU32 {
  /// number of possible outputs. outputs will be in `0 .. count`
  count: u32,
  /// Multiplication threshold thing. https://arxiv.org/abs/1805.10941
  threshold: u32,
}
impl BoundedRandU32 {
  /// Constructs a new `BoundedRandU32`.
  ///
  /// ## Panics
  /// If the count is 0.
  #[inline]
  pub const fn new(count: u32) -> Self {
    let threshold = count.wrapping_neg() % count;
    Self { count, threshold }
  }

  /// Constructs a new `BoundedRandU32`, or `None` on failure.
  ///
  /// ## Failure
  /// If the count is 0.
  #[inline]
  pub const fn try_new(count: u32) -> Option<Self> {
    if count > 0 {
      Some(Self::new(count))
    } else {
      None
    }
  }

  /// The number of possible outputs.
  #[inline]
  pub const fn count(self) -> u32 {
    self.count
  }

  /// Given a `u32`, place it into this bounded range.
  ///
  /// ## Failure
  /// * If the value is such that it doesn't fit evenly it is rejected.
  #[inline]
  pub const fn place_in_range(self, val: u32) -> Option<u32> {
    let mul: u64 = (val as u64).wrapping_mul(self.count as u64);
    let low_part: u32 = mul as u32;
    if low_part < self.threshold {
      None
    } else {
      //debug_assert!(((mul >> 32) as u32) < self.count());
      Some((mul >> 32) as u32)
    }
  }

  /// Given a gen, sample from the gen until `place_in_range` succeeds.
  #[inline]
  pub fn sample<G: Gen32 + ?Sized>(self, gen: &mut G) -> u32 {
    loop {
      if let Some(output) = self.place_in_range(gen.next_u32()) {
        return output;
      }
    }
  }
}

/// Stores data for a standard 1 through `N` sided die.
///
/// Produces values in `1 ..= N`
#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
#[repr(transparent)]
pub struct StandardDie(BoundedRandU32);
impl StandardDie {
  /// Constructs a new die.
  ///
  /// ## Panics
  /// If the count is 0.
  #[inline]
  pub const fn new(sides: u32) -> Self {
    Self(BoundedRandU32::new(sides))
  }

  /// The number of sides of this die.
  #[inline]
  pub const fn sides(self) -> i32 {
    self.0.count() as i32
  }

  /// Sample from the generator to get a die roll.
  #[inline]
  pub fn sample<G: Gen32 + ?Sized>(self, gen: &mut G) -> i32 {
    1 + self.0.sample(gen) as i32
  }
}

#[doc(hidden)]
pub const D4: StandardDie = StandardDie::new(4);
#[doc(hidden)]
pub const D6: StandardDie = StandardDie::new(6);
#[doc(hidden)]
pub const D8: StandardDie = StandardDie::new(8);
#[doc(hidden)]
pub const D10: StandardDie = StandardDie::new(10);
#[doc(hidden)]
pub const D12: StandardDie = StandardDie::new(12);
#[doc(hidden)]
pub const D20: StandardDie = StandardDie::new(20);

/// Stores data for an "exploding" 1 through `N` sided die.
///
/// When rolled, if a maximum value is rolled, then the die is rolled again and
/// added to the total. Successive rolls can also trigger additional rolls on a
/// maximum value.
#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
#[repr(transparent)]
pub struct ExplodingDie(StandardDie);
impl ExplodingDie {
  /// Constructs an exploding die.
  ///
  /// ## Panics
  /// If the count is 0.
  #[inline]
  pub const fn new(sides: u32) -> Self {
    Self(StandardDie::new(sides))
  }

  /// The number of sides of this die.
  #[inline]
  pub const fn sides(self) -> i32 {
    self.0.sides()
  }

  /// Sample from the generator to perform an exploding roll.
  #[inline]
  pub fn sample<G: Gen32 + ?Sized>(self, gen: &mut G) -> i32 {
    let mut t: i32 = 0;
    while self.0.sample(gen) == self.0.sides() {
      t = t.wrapping_add(1).wrapping_add(self.sides());
    }
    t.wrapping_add(self.0.sample(gen) as i32)
  }
}

#[doc(hidden)]
pub const X4: ExplodingDie = ExplodingDie::new(4);
#[doc(hidden)]
pub const X6: ExplodingDie = ExplodingDie::new(6);
#[doc(hidden)]
pub const X8: ExplodingDie = ExplodingDie::new(8);
#[doc(hidden)]
pub const X10: ExplodingDie = ExplodingDie::new(10);
#[doc(hidden)]
pub const X12: ExplodingDie = ExplodingDie::new(12);
#[doc(hidden)]
pub const X20: ExplodingDie = ExplodingDie::new(20);

/// Returns `k` with probability `2^(-k-1)`, a "binary exponential
/// distribution".
#[inline]
fn next_binary_exp_distr<G: Gen32 + ?Sized>(g: &mut G) -> u32 {
  let r: u32 = g.next_u32();
  if r > 0 {
    r.trailing_zeros()
  } else {
    32 + next_binary_exp_distr(g)
  }
}

/// Gives an `f32` output, in the unsigned (`[0,1]`) or signed (`[-1, 1]`)
/// range.
fn ieee754_random_f32<G: Gen32 + ?Sized>(g: &mut G, signed: bool) -> f32 {
  // Returns random number in [0, 1] or [-1, 1] depending on signed.
  let bit_width = 32;
  let exponent_bias = 127;
  let num_mantissa_bits = 23;
  let num_rest_bits = bit_width - num_mantissa_bits - 1 - signed as i32;
  let r: u32 = g.next_u32();

  debug_assert!(num_rest_bits >= 0);
  debug_assert!(core::mem::size_of::<u32>() * 8 == bit_width as _);

  let mantissa = r >> (bit_width - num_mantissa_bits);
  let (sign_mask, rand_bit, rest_bits);
  if signed {
    sign_mask = r << (bit_width - 1);
    rand_bit = (r & 2) != 0;
    rest_bits = (r >> 2) & ((1 << num_rest_bits) - 1);
  } else {
    sign_mask = 0;
    rand_bit = (r & 1) != 0;
    rest_bits = (r >> 1) & ((1 << num_rest_bits) - 1);
  }

  // If our mantissa is zero, half of the time we must increase our exponent.
  let increment_exponent = (mantissa == 0 && rand_bit) as i32;

  // We can use rest_bits to save more calls to the rng.
  let mut exponent: i32 = -1 + increment_exponent
    - if rest_bits > 0 {
      rest_bits.trailing_zeros() as i32
    } else {
      num_rest_bits + next_binary_exp_distr(g) as i32
    };

  // It is very unlikely our exponent is invalid at this point, but keep
  // regenerating it until it is valid.
  while exponent < -exponent_bias || exponent > 0 {
    exponent = -1 + increment_exponent - next_binary_exp_distr(g) as i32;
  }

  f32::from_bits(
    sign_mask
      | (((exponent + exponent_bias) as u32) << num_mantissa_bits)
      | mantissa,
  )
}

/// Gives an `f64` output, in the unsigned (`[0,1]`) or signed (`[-1, 1]`)
/// range.
#[allow(dead_code)]
fn ieee754_random_f64<G: Gen32 + ?Sized>(g: &mut G, signed: bool) -> f64 {
  // Returns random number in [0, 1] or [-1, 1] depending on signed.
  let bit_width = 64;
  let exponent_bias = 1023;
  let num_mantissa_bits = 52;
  let num_rest_bits = bit_width - num_mantissa_bits - 1 - signed as i32;
  let r: u64 = g.next_u64();

  debug_assert!(num_rest_bits >= 0);
  debug_assert!(core::mem::size_of::<u64>() * 8 == bit_width as _);

  let mantissa = r >> (bit_width - num_mantissa_bits);
  let (sign_mask, rand_bit, rest_bits);
  if signed {
    sign_mask = r << (bit_width - 1);
    rand_bit = (r & 2) != 0;
    rest_bits = (r >> 2) & ((1 << num_rest_bits) - 1);
  } else {
    sign_mask = 0;
    rand_bit = (r & 1) != 0;
    rest_bits = (r >> 1) & ((1 << num_rest_bits) - 1);
  }

  // If our mantissa is zero, half of the time we must increase our exponent.
  let increment_exponent = (mantissa == 0 && rand_bit) as i32;

  // We can use rest_bits to save more calls to the rng.
  let mut exponent: i32 = -1 + increment_exponent
    - if rest_bits > 0 {
      rest_bits.trailing_zeros() as i32
    } else {
      num_rest_bits + next_binary_exp_distr(g) as i32
    };

  // It is very unlikely our exponent is invalid at this point, but keep
  // regenerating it until it is valid.
  while exponent < -exponent_bias || exponent > 0 {
    exponent = -1 + increment_exponent - next_binary_exp_distr(g) as i32;
  }

  f64::from_bits(
    sign_mask
      | (((exponent + exponent_bias) as u64) << num_mantissa_bits)
      | mantissa,
  )
}