// Prevent roundoff at unit-vector boundaries from leaking NaN into API results.

///|
fn clamp_unit(x : Double) -> Double {
  x.max(-1.0).min(1.0)
}

///|
fn finite(x : Double) -> Bool {
  x == x && x.abs() < 1.0 / 0.0
}

// Roots in [-1,1] of the parabola through (-1,h0), (0,h1), (1,h2).
// Preserve upstream arithmetic for normal quadratics; handle the degenerate linear case.

///|
fn interval_roots(
  h0 : Double,
  h1 : Double,
  h2 : Double,
) -> (Int, Double, Double, Double) {
  let a = (h0 + h2) / 2.0 - h1
  let b = (h2 - h0) / 2.0
  if a.abs() < 1.0e-14 {
    if b.abs() < 1.0e-14 {
      return (0, 0.0, 0.0, h1)
    }
    let x = -h1 / b
    return (if x.abs() <= 1.0 { 1 } else { 0 }, x, x, h1)
  }
  let xe = -b / (2.0 * a)
  let disc = b * b - 4.0 * a * h1
  let ye = (a * xe + b) * xe + h1
  if disc < 0.0 {
    return (0, 0.0, 0.0, ye)
  }
  let dx = disc.sqrt() / (a.abs() * 2.0)
  let mut x1 = xe - dx
  let x2 = xe + dx
  let mut roots = 0
  if x1.abs() <= 1.0 {
    roots += 1
  }
  if x2.abs() <= 1.0 {
    roots += 1
  }
  if x1 < -1.0 {
    x1 = x2
  }
  (roots, x1, x2, ye)
}