///|
let gas_constant_j_per_mol_k : Double = 8.31446261815324

///|
fn abs_double(x : Double) -> Double {
  if x < 0.0 {
    0.0 - x
  } else {
    x
  }
}

///|
fn clamp_positive(x : Double) -> Double {
  if x < 1.0E-12 {
    1.0E-12
  } else {
    x
  }
}

///|
fn sum(xs : Array[Double]) -> Double {
  for i = 0, acc = 0.0; i < xs.length(); {
    continue i + 1, acc + xs[i]
  } nobreak {
    acc
  }
}

///|
pub fn normalize(xs : Array[Double]) -> Array[Double] raise VleError {
  if xs.length() == 0 {
    raise VleError::EmptyMixture
  }
  for i = 0; i < xs.length(); i = i + 1 {
    if xs[i] < 0.0 {
      raise VleError::InvalidFraction(index=i, value=xs[i])
    }
  }
  let total = sum(xs)
  if total <= 0.0 {
    raise VleError::InvalidFraction(index=0, value=total)
  }
  [
    for x in xs => x / total
  ]
}

///|
fn assert_same_length(a : Int, b : Int) -> Unit raise VleError {
  if a != b {
    raise VleError::LengthMismatch(expected=a, actual=b)
  }
}

///|
fn assert_temperature(t : Double) -> Unit raise VleError {
  if t <= 0.0 {
    raise VleError::NonPositiveTemperature(t)
  }
}

///|
fn assert_pressure(p : Double) -> Unit raise VleError {
  if p <= 0.0 {
    raise VleError::NonPositivePressure(p)
  }
}

///|
fn bisect(
  low~ : Double,
  high~ : Double,
  tolerance~ : Double,
  max_iter~ : Int,
  f : (Double) -> Double raise VleError,
) -> (Double, Int) raise VleError {
  let f_low = f(low)
  let f_high = f(high)
  if f_low * f_high > 0.0 {
    raise VleError::SolverDidNotBracket(low~, high~)
  }
  for i = 0, lo = low, hi = high, flo = f_low; i < max_iter; {
    let mid = (lo + hi) / 2.0
    let fm = f(mid)
    if abs_double(fm) <= tolerance || abs_double(hi - lo) <= tolerance {
      break (mid, i + 1)
    }
    if flo * fm <= 0.0 {
      continue i + 1, lo, mid, flo
    } else {
      continue i + 1, mid, hi, fm
    }
  } nobreak {
    raise VleError::SolverDidNotConverge(iterations=max_iter)
  }
}