///|
pub fn dot(a : Array[Double], b : Array[Double]) -> Double {
  if a.length() != b.length() {
    return 0.0
  }
  let mut result = 0.0
  for i in 0.. Array[Double] {
  let n = if a.length() < b.length() { a.length() } else { b.length() }
  let result = Array::new(capacity=n)
  for i in 0.. Array[Double] {
  let n = if a.length() < b.length() { a.length() } else { b.length() }
  let result = Array::new(capacity=n)
  for i in 0.. Array[Double] {
  let result = Array::new(capacity=a.length())
  for value in a {
    result.push(value * factor)
  }
  result
}

///|
pub fn l2_norm(a : Array[Double]) -> Double {
  dot(a, a).sqrt()
}

///|
pub fn l1_norm(a : Array[Double]) -> Double {
  let mut result = 0.0
  for value in a {
    result += value.abs()
  }
  result
}

///|
pub fn matrix_zeros(rows : Int, columns : Int) -> Array[Array[Double]] {
  let result : Array[Array[Double]] = Array::new(
    capacity=if rows > 0 { rows } else { 0 },
  )
  for _ in 0.. Array[Array[Double]] {
  let result = matrix_zeros(size, size)
  for i in 0.. Array[Array[Double]] {
  if matrix.length() == 0 {
    return []
  }
  let columns = matrix[0].length()
  let result = matrix_zeros(columns, matrix.length())
  for i in 0.. Array[Array[Double]] {
  if left.length() == 0 ||
    right.length() == 0 ||
    left[0].length() != right.length() {
    return []
  }
  let result = matrix_zeros(left.length(), right[0].length())
  for i in 0.. Array[Double] {
  let result = Array::new(capacity=matrix.length())
  for row in matrix {
    result.push(dot(row, vector))
  }
  result
}

///|
/// Solves a dense linear system using Gaussian elimination with partial pivoting.
pub fn solve_linear_system(
  matrix : Array[Array[Double]],
  right_hand_side : Array[Double],
) -> Array[Double] {
  let n = matrix.length()
  if n == 0 || right_hand_side.length() != n {
    return []
  }
  let augmented = matrix_zeros(n, n + 1)
  for i in 0.. augmented[best][pivot].abs() {
        best = row
      }
    }
    if augmented[best][pivot].abs() < 1.0e-12 {
      return []
    }
    if best != pivot {
      let temporary = augmented[pivot]
      augmented[pivot] = augmented[best]
      augmented[best] = temporary
    }
    let divisor = augmented[pivot][pivot]
    for column in pivot..<(n + 1) {
      augmented[pivot][column] /= divisor
    }
    for row in 0.. Array[Array[Double]] {
  let n = matrix.length()
  if n == 0 {
    return []
  }
  let result = matrix_zeros(n, n)
  for column in 0.. Array[Array[Double]] {
  let result = matrix.copy()
  for i in 0.. Array[Array[Double]] {
  if rows.length() == 0 {
    return []
  }
  let columns = rows[0].length()
  let means = Array::make(columns, 0.0)
  for row in rows {
    for column in 0.. 1 {
    rows.length() - 1
  } else {
    rows.length()
  }
  if denominator > 0 {
    for i in 0.. Array[Array[Double]] {
  if rows.length() == 0 || rows.length() != weights.length() {
    return []
  }
  let columns = rows[0].length()
  let means = Array::make(columns, 0.0)
  let total = sum(weights)
  if total == 0.0 {
    return matrix_zeros(columns, columns)
  }
  for i in 0..