Matrix Operations in JS

matrix.js

withAs = (obj, cb) => cb(obj)
sum = arr => arr.reduce((a, b) => a + b)
mul = arr => arr.reduce((a, b) => a * b)
sub = arr => arr.splice(1).reduce((a, b) => a - b, arr[0])
deepClone = obj => JSON.parse(JSON.stringify(obj))

shifter = (arr, step) => [
  ...arr.splice(step),
  ...arr.splice(arr.length - step)
]

makeArray = (n, cb = i => i) =>
  [...Array(n).keys()].map(cb)

makeMatrix = (len, wid, fill) =>
  makeArray(len).map(i => makeArray(
    wid, j => fill ? fill(i, j) : 0
  ))

matrixRandom =
  (len, wid, min=0, max=100) =>
  makeMatrix(len, wid, x => Math.round(
    Math.random() * (max - min)
  ) + min)

matrixSize = matrix =>
  [matrix.length, matrix[0].length]

arr2mat = (len, wid, arr) =>
  makeArray(len).map(i => arr.slice(
    i * wid, i * wid + wid
  ))

matrixMap = (matrix, cb) =>
  deepClone(matrix).map((i, ix) => i.map(
    (j, jx) => cb({i, ix, j, jx, matrix})
  ))

matrixScalar = (n, matrix) =>
  matrixMap(matrix, ({j}) => n * j)

matrixAdd = matrices =>
  matrices.reduce((acc, inc) => matrixMap(
    acc, ({j, ix, jx}) => j + inc[ix][jx]
  ), makeMatrix(...matrixSize(matrices[0])))

matrixSub = matrices => matrices.splice(1)
  .reduce((acc, inc) => matrixMap(
    acc, ({j, ix, jx}) => j - inc[ix][jx]
  ), matrices[0])

matrixMul = (m1, m2) => makeMatrix(
  m1.length, m2[0].length, (i, j) => sum(
    m1[i].map((k, kx) => k * m2[kx][j])
  )
)

matrixMuls = matrices =>
  deepClone(matrices).splice(1)
  .reduce((acc, inc) => makeMatrix(
    acc.length, inc[0].length,
    (ix, jx) => sum(acc[ix].map(
      (k, kx) => k * inc[kx][jx]
    ))
  ), deepClone(matrices[0]))

matrixMinor = (matrix, row, col) =>
  matrix.length < 3 ? matrix :
  matrix.filter((i, ix) => ix !== row - 1)
  .map(i => i.filter((j, jx) => jx !== col - 1))

matrixTrans = matrix => makeMatrix(
  ...shifter(matrixSize(matrix), 1),
  (i, j) => matrix[j][i]
)

matrixDet = matrix => withAs(
  deepClone(matrix), clone => matrix.length < 3
  ? sub(matrixTrans(clone.map(shifter)).map(mul))
  : sum(clone[0].map((i, ix) => matrixDet(
    matrixMinor(matrix, 1, ix+1)
  ) * Math.pow(-1, ix+2) * i))
)

matrixCofactor = matrix => matrixMap(
  matrix, ({i, ix, j, jx}) =>
  matrix[0].length > 2 ? 
  Math.pow(-1, ix + jx + 2) * matrixDet(
    matrixMinor(matrix, ix+1, jx+1)
  ) : (ix != jx ? -matrix[jx][ix]
    : matrix[+!ix][+!jx]
  )
)

matrixInverse = matrix => matrixMap(
  matrixTrans(matrixCofactor(matrix)),
  ({j}) => j / matrixDet(matrix)
)

matrixLinSolve = (conds, res) => matrixMul(
  matrixInverse(conds), res.map(i => [i])
).map(i => i[0])

matrixPolySolve = (xs, ys) =>
  matrixLinSolve(xs.map(i => makeArray(
    ys.length, j => Math.pow(i, j)
  ).reverse()), ys)

makeIdentity = (len, wid) => makeMatrix(
  len, wid, (i, j) => +(i === j)
)

matrixPow = (pow, matrix) =>
  matrixMuls(makeArray(pow, x => matrix))

facto = n => n === 1
  ? 1 : n * facto(n - 1)

matrixExpo = (depth, matrix) => matrixAdd([
  makeIdentity(...matrixSize(matrix)),
  ...makeArray(depth, i => matrixMap(
    matrixPow(i + 1, matrix),
    ({j}) => j / facto(i + 1)
  ))
])