From 27e7e9b2f6e80c930d090dd203d8243f179f6e6e Mon Sep 17 00:00:00 2001 From: jobo322 Date: Fri, 13 Jun 2025 15:04:54 -0500 Subject: [PATCH 01/33] feat: add new matrix addition method and improve performance of mmul and kroneckerProduct method --- src/__tests__/test.js | 20 +++++++++++++++- src/index.js | 54 ++++++++++++++++++++++++++++++++----------- 2 files changed, 60 insertions(+), 14 deletions(-) diff --git a/src/__tests__/test.js b/src/__tests__/test.js index 3e29942..56076b1 100644 --- a/src/__tests__/test.js +++ b/src/__tests__/test.js @@ -1,6 +1,24 @@ -import { SparseMatrix } from '..'; +import { SparseMatrix } from '../index'; describe('Sparse Matrix', () => { + it('add', () => { + let m1 = new SparseMatrix([ + [2, 0, 1], + [0, 0, 3], + [2, 0, 1], + ]); + let m2 = new SparseMatrix([ + [0, 1, 5], + [2, 0, 0], + [-2, 0, -1], + ]); + let m3 = m1.add(m2).to2DArray(); + expect(m3).toStrictEqual([ + [2, 1, 6], + [2, 0, 3], + [0, 0, 0], + ]); + }); it('mmul', () => { let m1 = new SparseMatrix([ [2, 0, 1], diff --git a/src/index.js b/src/index.js index 1c2907f..eff2077 100644 --- a/src/index.js +++ b/src/index.js @@ -17,11 +17,12 @@ export class SparseMatrix { if (Array.isArray(rows)) { const matrix = rows; - rows = matrix.length; - options = columns || {}; + const nbRows = matrix.length; + const nbColmuns = matrix[0].length; + options = columns || { initialCapacity: nbRows * nbColmuns }; columns = matrix[0].length; - this._init(rows, columns, new HashTable(options), options.threshold); - for (let i = 0; i < rows; i++) { + this._init(nbRows, columns, new HashTable(options), options.threshold); + for (let i = 0; i < nbRows; i++) { for (let j = 0; j < columns; j++) { let value = matrix[i][j]; if (this.threshold && Math.abs(value) < this.threshold) value = 0; @@ -30,6 +31,7 @@ export class SparseMatrix { } } } + this.elements.maybeShrinkCapacity(); } else { this._init(rows, columns, new HashTable(options), options.threshold); } @@ -132,6 +134,19 @@ export class SparseMatrix { return this; } + add(other) { + if (typeof other === 'number') return this.addS(other); + if (this.rows !== other.rows || this.columns !== other.columns) { + throw new RangeError('Matrices dimensions must be equal'); + } + + other.withEachNonZero((i, j, v) => { + this.set(i, j, v + this.get(i, j)); + }); + + return this; + } + mmul(other) { if (this.columns !== other.rows) { // eslint-disable-next-line no-console @@ -144,14 +159,12 @@ export class SparseMatrix { const p = other.columns; const result = new SparseMatrix(m, p); - this.forEachNonZero((i, j, v1) => { - other.forEachNonZero((k, l, v2) => { + this.withEachNonZero((i, j, v1) => { + other.withEachNonZero((k, l, v2) => { if (j === k) { result.set(i, l, result.get(i, l) + v1 * v2); } - return v2; }); - return v1; }); return result; } @@ -165,16 +178,31 @@ export class SparseMatrix { const result = new SparseMatrix(m * p, n * q, { initialCapacity: this.cardinality * other.cardinality, }); - this.forEachNonZero((i, j, v1) => { - other.forEachNonZero((k, l, v2) => { - result.set(p * i + k, q * j + l, v1 * v2); - return v2; + + this.withEachNonZero((i, j, v1) => { + const pi = p * i; + const qj = q * j; + other.withEachNonZero((k, l, v2) => { + result.set(pi + k, qj + l, v1 * v2); }); - return v1; }); return result; } + withEachNonZero(callback) { + const { state, table, values } = this.elements; + const nbStates = state.length; + const activeIndex = []; + for (let i = 0; i < nbStates; i++) { + if (state[i] === 1) activeIndex.push(i); + } + const columns = this.columns; + for (const i of activeIndex) { + const key = table[i]; + callback((key / columns) | 0, key % columns, values[i]); + } + } + forEachNonZero(callback) { this.elements.forEachPair((key, value) => { const i = (key / this.columns) | 0; From e4a7002b5cb5cbae205b8e14b388c904d45bc4c9 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Fri, 13 Jun 2025 16:09:36 -0500 Subject: [PATCH 02/33] feat: improve performance of mmul and kronecker product --- src/__tests__/test.js | 4 +++ src/index.js | 68 ++++++++++++++++++++++++++++++++----------- 2 files changed, 55 insertions(+), 17 deletions(-) diff --git a/src/__tests__/test.js b/src/__tests__/test.js index 56076b1..0bfca72 100644 --- a/src/__tests__/test.js +++ b/src/__tests__/test.js @@ -36,6 +36,10 @@ describe('Sparse Matrix', () => { expect(m3.cardinality).toBe(1); expect(m3.get(0, 1)).toBe(2); + expect(m3.to2DArray()).toStrictEqual([ + [0, 2], + [0, 0], + ]); }); it('kronecker', () => { diff --git a/src/index.js b/src/index.js index eff2077..114cd4b 100644 --- a/src/index.js +++ b/src/index.js @@ -158,14 +158,32 @@ export class SparseMatrix { const m = this.rows; const p = other.columns; + const { + columns: otherCols, + rows: otherRows, + values: otherValues, + } = other.getNonZeros(); + const { + columns: thisCols, + rows: thisRows, + values: thisValues, + } = this.getNonZeros(); + const result = new SparseMatrix(m, p); - this.withEachNonZero((i, j, v1) => { - other.withEachNonZero((k, l, v2) => { - if (j === k) { - result.set(i, l, result.get(i, l) + v1 * v2); + + const nbOtherActive = otherCols.length; + const nbThisActive = thisCols.length; + for (let t = 0; t < nbThisActive; t++) { + const i = thisRows[t]; + const j = thisCols[t]; + for (let o = 0; o < nbOtherActive; o++) { + if (j === otherRows[o]) { + const l = otherCols[o]; + result.set(i, l, result.get(i, l) + otherValues[o] * thisValues[t]); } - }); - }); + } + } + return result; } @@ -179,13 +197,31 @@ export class SparseMatrix { initialCapacity: this.cardinality * other.cardinality, }); - this.withEachNonZero((i, j, v1) => { - const pi = p * i; - const qj = q * j; - other.withEachNonZero((k, l, v2) => { - result.set(pi + k, qj + l, v1 * v2); - }); - }); + const { + columns: otherCols, + rows: otherRows, + values: otherValues, + } = other.getNonZeros(); + const { + columns: thisCols, + rows: thisRows, + values: thisValues, + } = this.getNonZeros(); + + const nbOtherActive = otherCols.length; + const nbThisActive = thisCols.length; + for (let t = 0; t < nbThisActive; t++) { + const pi = p * thisRows[t]; + const qj = q * thisCols[t]; + for (let o = 0; o < nbOtherActive; o++) { + result.set( + pi + otherRows[o], + qj + otherCols[o], + otherValues[o] * thisValues[t], + ); + } + } + return result; } @@ -229,12 +265,11 @@ export class SparseMatrix { const columns = new Array(cardinality); const values = new Array(cardinality); let idx = 0; - this.forEachNonZero((i, j, value) => { + this.withEachNonZero((i, j, value) => { rows[idx] = i; columns[idx] = j; values[idx] = value; idx++; - return value; }); return { rows, columns, values }; } @@ -254,9 +289,8 @@ export class SparseMatrix { let trans = new SparseMatrix(this.columns, this.rows, { initialCapacity: this.cardinality, }); - this.forEachNonZero((i, j, value) => { + this.withEachNonZero((i, j, value) => { trans.set(j, i, value); - return value; }); return trans; } From 9a44f16b1e674c6e86dffe88ad6b68689fd25da1 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Thu, 19 Jun 2025 23:27:15 -0500 Subject: [PATCH 03/33] feat: add new multiplication methods and improve getNonZeros function --- src/index.js | 139 ++++++++++++++++++++++++++++++++++++++++++++++++--- 1 file changed, 131 insertions(+), 8 deletions(-) diff --git a/src/index.js b/src/index.js index 114cd4b..90f185a 100644 --- a/src/index.js +++ b/src/index.js @@ -131,6 +131,7 @@ export class SparseMatrix { } else { this.elements.set(row * this.columns + column, value); } + return this; } @@ -147,7 +148,21 @@ export class SparseMatrix { return this; } - mmul(other) { + mulNew(other) { + if (typeof other !== 'number') { + throw new RangeError('the argument should be a number'); + } + + // if (this.cardinality / this.columns / this.rows > 0.1) + this.withEachNonZero((i, j, v) => { + // return v * other; + this.set(i, j, v * other); + }); + + return this; + } + + mmulNew(other) { if (this.columns !== other.rows) { // eslint-disable-next-line no-console console.warn( @@ -169,7 +184,7 @@ export class SparseMatrix { values: thisValues, } = this.getNonZeros(); - const result = new SparseMatrix(m, p); + const result = new SparseMatrix(m, p, { initialCapacity: m * p }); const nbOtherActive = otherCols.length; const nbThisActive = thisCols.length; @@ -187,6 +202,53 @@ export class SparseMatrix { return result; } + mmul(other) { + if (this.columns !== other.rows) { + // eslint-disable-next-line no-console + console.warn( + 'Number of columns of left matrix are not equal to number of rows of right matrix.', + ); + } + + const m = this.rows; + const p = other.columns; + + const { + columns: otherCols, + rows: otherRows, + values: otherValues, + } = other.getNonZeros('csr'); + const { + columns: thisCols, + rows: thisRows, + values: thisValues, + } = this.getNonZeros('csc'); + + const result = new SparseMatrix(m, p, { initialCapacity: m * p + 20 }); + + for (let t = 0; t < thisCols.length - 1; t++) { + const j = t; + const tValues = thisValues.subarray(thisCols[t], thisCols[t + 1]); + const tRows = thisRows.subarray(thisCols[t], thisCols[t + 1]); + let initOther = 0; + for (let o = initOther; o < otherRows.length - 1; o++) { + if (o === j) { + initOther++; + const oValues = otherValues.subarray(otherRows[o], otherRows[o + 1]); + const oCols = otherCols.subarray(otherRows[o], otherRows[o + 1]); + for (let f = 0; f < tValues.length; f++) { + for (let k = 0; k < oValues.length; k++) { + const i = tRows[f]; + const l = oCols[k]; + result.set(i, l, result.get(i, l) + tValues[f] * oValues[k]); + } + } + } else if (j < o) break; + } + } + return result; + } + kroneckerProduct(other) { const m = this.rows; const n = this.columns; @@ -194,7 +256,7 @@ export class SparseMatrix { const q = other.columns; const result = new SparseMatrix(m * p, n * q, { - initialCapacity: this.cardinality * other.cardinality, + initialCapacity: this.cardinality * other.cardinality + 20, }); const { @@ -259,11 +321,12 @@ export class SparseMatrix { return this; } - getNonZeros() { + //'csr' | 'csc' | undefined + getNonZeros(format) { const cardinality = this.cardinality; - const rows = new Array(cardinality); - const columns = new Array(cardinality); - const values = new Array(cardinality); + const rows = new Float64Array(cardinality); + const columns = new Float64Array(cardinality); + const values = new Float64Array(cardinality); let idx = 0; this.withEachNonZero((i, j, value) => { rows[idx] = i; @@ -271,7 +334,13 @@ export class SparseMatrix { values[idx] = value; idx++; }); - return { rows, columns, values }; + + const cooMatrix = { rows, columns, values }; + + if (!format) return cooMatrix; + + const csrMatrix = cooToCsr(cooMatrix, this.rows); + return format === 'csc' ? csrToCsc(csrMatrix, this.columns) : csrMatrix; } setThreshold(newThreshold) { @@ -452,3 +521,57 @@ function fillTemplateFunction(template, values) { } return template; } + +function csrToCsc(csrMatrix, numCols) { + const { + values: csrValues, + columns: csrColIndices, + rows: csrRowPtr, + } = csrMatrix; + // Initialize CSC arrays + const cscValues = new Float64Array(csrValues.length); + const cscRowIndices = new Float64Array(csrValues.length); + const cscColPtr = new Float64Array(numCols + 1); + + // Count non-zeros per column + for (let i = 0; i < csrColIndices.length; i++) { + cscColPtr[csrColIndices[i] + 1]++; + } + + // Compute column pointers (prefix sum) + for (let i = 1; i <= numCols; i++) { + cscColPtr[i] += cscColPtr[i - 1]; + } + + // Temporary copy for filling values + const next = cscColPtr.slice(); + + // Fill CSC arrays + for (let row = 0; row < csrRowPtr.length - 1; row++) { + for (let j = csrRowPtr[row], i = 0; j < csrRowPtr[row + 1]; j++, i++) { + const col = csrColIndices[j]; + const pos = next[col]; + cscValues[pos] = csrValues[j]; + cscRowIndices[pos] = row; + next[col]++; + } + // if (row === 1) break; + } + + return { rows: cscRowIndices, columns: cscColPtr, values: cscValues }; +} + +function cooToCsr(cooMatrix, nbRows = 9) { + const { values, columns, rows } = cooMatrix; + //could not be the same length + const csrRowPtr = new Float64Array(nbRows + 1); + const length = values.length; + let currentRow = rows[0]; + for (let index = 0; index < length; ) { + const prev = index; + while (currentRow === rows[index] && index < length) ++index; + csrRowPtr[currentRow + 1] = index; + currentRow += 1; + } + return { rows: csrRowPtr, columns, values }; +} From f0b16b20f8d924a27a6884598f585d1ef149f7e2 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Sat, 21 Jun 2025 00:16:15 -0500 Subject: [PATCH 04/33] feat: add benchmark for sparse matrix operations and performance comparison --- benchmark/benchmark.js | 204 ++++ benchmark/class/SparseMatrixOld.js | 1386 ++++++++++++++++++++++++++++ 2 files changed, 1590 insertions(+) create mode 100644 benchmark/benchmark.js create mode 100644 benchmark/class/SparseMatrixOld.js diff --git a/benchmark/benchmark.js b/benchmark/benchmark.js new file mode 100644 index 0000000..bb1eac5 --- /dev/null +++ b/benchmark/benchmark.js @@ -0,0 +1,204 @@ +import { SparseMatrix } from '../src/index.js'; +import { Matrix } from 'ml-matrix'; +import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; +import fs from 'fs'; + +function randomSparseMatrix(rows, cols, density = 0.01) { + const matrix = []; + for (let i = 0; i < rows; i++) { + const row = new Array(cols).fill(0); + for (let j = 0; j < cols; j++) { + if (Math.random() < density) { + row[j] = Math.random() * 10; + } + } + matrix.push(row); + } + return new SparseMatrix(matrix); +} + +function benchmark(fn, label, iterations = 5) { + const times = []; + for (let i = 0; i < iterations; i++) { + const t0 = performance.now(); + fn(); + const t1 = performance.now(); + times.push(t1 - t0); + } + const avg = times.reduce((a, b) => a + b, 0) / times.length; + console.log(`${label}: avg ${avg.toFixed(2)} ms over ${iterations} runs`); + return avg; +} + +function printWinner(label1, avg1, label2, avg2) { + let winner, loser, win, lose; + if (avg1 < avg2) { + winner = label1; + win = avg1; + loser = label2; + lose = avg2; + } else { + winner = label2; + win = avg2; + loser = label1; + lose = avg1; + } + + const percent = ((lose - win) / lose) * 100; + console.log( + ` -> ${winner} was ${(lose / win).toFixed(2)}x faster (${percent.toFixed( + 1, + )}% faster) than ${loser}\n`, + ); +} + +function runBenchmarks() { + const m = 128; + const n = 64; + const p = 128; + const densityA = 0.03; + const densityB = 0.03; + const A = randomSparseMatrix(m, n, densityA); + const B = randomSparseMatrix(n, p, densityB); + + let denseA = A.to2DArray(); + let denseB = B.to2DArray(); + + const AOld = new SparseMatrixOld(denseA); + const BOld = new SparseMatrixOld(denseB); + A.cardinality; + denseA = new Matrix(denseA); + denseB = new Matrix(denseB); + + // 1. add vs addNew + // const addAvg = benchmark(() => { + // const a = AOld.clone(); + // a.add(BOld); + // }, 'add'); + + // const addNewAvg = benchmark(() => { + // const a = A.clone(); + // a.add(B); + // }, 'addNew'); + + // printWinner('add', addAvg, 'addNew', addNewAvg); + + // 2. mmul vs mmulNew + + const mmulNewAvg = benchmark( + () => { + A.mmul(B); + }, + 'mmulNew', + 3, + ); + + const mmulAvg = benchmark( + () => { + AOld.mmul(BOld); + }, + 'mmul', + 3, + ); + + const denseAvg = benchmark( + () => { + denseA.mmul(denseB); + }, + 'denseMatrix', + 3, + ); + + printWinner('mmul', mmulAvg, 'mmulNew', mmulNewAvg); + + // 3. kroneckerProduct vs kroneckerProductNew + // const kronNewAvg = benchmark(() => { + // A.kroneckerProduct(B); + // }, 'kroneckerProductNew'); + // const kronAvg = benchmark(() => { + // AOld.kroneckerProduct(BOld); + // }, 'kroneckerProduct'); + + // printWinner('kroneckerProduct', kronAvg, 'kroneckerProductNew', kronNewAvg); + + // 4. matrix multiplication + // const mulAvg = benchmark(() => { + // A.mul(5); + // }, 'mul'); + + // const mulNewAvg = benchmark(() => { + // AOld.mul(5); + // }, 'mulNew'); + + // printWinner('mul', mulAvg, 'mulNew', mulNewAvg); +} + +function runSizeSweepBenchmark() { + const sizes = [8, 16, 32, 64, 128]; + const densities = [0.01, 0.015, 0.02, 0.025, 0.03, 0.35]; + const results = []; + + for (const densityA of densities) { + for (const densityB of densities) { + for (const m of sizes) { + for (const n of sizes) { + for (const p of sizes) { + // A: m x n, B: n x p + + const A = randomSparseMatrix(m, n, densityA); + const B = randomSparseMatrix(n, p, densityB); + + let denseA = A.to2DArray(); + let denseB = B.to2DArray(); + + const AOld = new SparseMatrixOld(denseA); + const BOld = new SparseMatrixOld(denseB); + + denseA = new Matrix(denseA); + denseB = new Matrix(denseB); + + const mmulNewAvg = benchmark( + () => { + A.mmul(B); + }, + 'mmulNew', + 3, + ); + + const mmulAvg = benchmark( + () => { + AOld.mmul(BOld); + }, + 'mmul', + 3, + ); + + const denseAvg = benchmark(() => { + denseA.mmul(denseB), 'denseMatrix', 3; + }); + + results.push({ + densityA, + densityB, + A_shape: [m, n], + B_shape: [n, p], + dense: denseAvg, + mmulNew: mmulNewAvg, + mmul: mmulAvg, + }); + } + } + } + } + } + + fs.writeFileSync( + './benchmark/size_sweep_results-dense.json', + JSON.stringify(results, null, 2), + ); + console.log('Size sweep benchmark results saved to size_sweep_results.json'); +} + +runBenchmarks(); +// Uncomment to run the size sweep benchmark +// runSizeSweepBenchmark(); diff --git a/benchmark/class/SparseMatrixOld.js b/benchmark/class/SparseMatrixOld.js new file mode 100644 index 0000000..f770b91 --- /dev/null +++ b/benchmark/class/SparseMatrixOld.js @@ -0,0 +1,1386 @@ +import HashTable from 'ml-hash-table'; + +/** @typedef {(row: number, column: number, value: number) => number | false} ForEachNonZeroCallback */ + +export class SparseMatrix { + constructor(rows, columns, options = {}) { + if (rows instanceof SparseMatrix) { + // clone + const other = rows; + this._init( + other.rows, + other.columns, + other.elements.clone(), + other.threshold, + ); + return; + } + + if (Array.isArray(rows)) { + const matrix = rows; + rows = matrix.length; + options = columns || {}; + columns = matrix[0].length; + this._init(rows, columns, new HashTable(options), options.threshold); + for (let i = 0; i < rows; i++) { + for (let j = 0; j < columns; j++) { + let value = matrix[i][j]; + if (this.threshold && Math.abs(value) < this.threshold) value = 0; + if (value !== 0) { + this.elements.set(i * columns + j, matrix[i][j]); + } + } + } + } else { + this._init(rows, columns, new HashTable(options), options.threshold); + } + } + + _init(rows, columns, elements, threshold) { + this.rows = rows; + this.columns = columns; + this.elements = elements; + this.threshold = threshold || 0; + } + + static eye(rows = 1, columns = rows) { + const min = Math.min(rows, columns); + const matrix = new SparseMatrix(rows, columns, { initialCapacity: min }); + for (let i = 0; i < min; i++) { + matrix.set(i, i, 1); + } + return matrix; + } + + clone() { + return new SparseMatrix(this); + } + + /** + * @returns {number[][]} + */ + to2DArray() { + const copy = new Array(this.rows); + for (let i = 0; i < this.rows; i++) { + copy[i] = new Array(this.columns); + for (let j = 0; j < this.columns; j++) { + copy[i][j] = this.get(i, j); + } + } + return copy; + } + + isSquare() { + return this.rows === this.columns; + } + + isSymmetric() { + if (!this.isSquare()) return false; + + let symmetric = true; + this.forEachNonZero((i, j, v) => { + if (this.get(j, i) !== v) { + symmetric = false; + return false; + } + return v; + }); + return symmetric; + } + + /** + * Search for the wither band in the main diagonals + * @returns {number} + */ + bandWidth() { + let min = this.columns; + let max = -1; + this.forEachNonZero((i, j, v) => { + let diff = i - j; + min = Math.min(min, diff); + max = Math.max(max, diff); + return v; + }); + return max - min; + } + + /** + * Test if a matrix is consider banded using a threshold + * @param {number} width + * @returns {boolean} + */ + isBanded(width) { + let bandWidth = this.bandWidth(); + return bandWidth <= width; + } + + /** + * @returns {number} + */ + get cardinality() { + return this.elements.size; + } + + /** + * @returns {number} + */ + get size() { + return this.rows * this.columns; + } + + /** + * @param {number} row + * @param {number} column + * @returns {number} + */ + get(row, column) { + return this.elements.get(row * this.columns + column); + } + + /** + * @param {number} row + * @param {number} column + * @param {number} value + * @returns {this} + */ + set(row, column, value) { + if (this.threshold && Math.abs(value) < this.threshold) value = 0; + if (value === 0) { + this.elements.remove(row * this.columns + column); + } else { + this.elements.set(row * this.columns + column, value); + } + return this; + } + + /** + * @param {SparseMatrix} other + * @returns {SparseMatrix} + */ + mmul(other) { + if (this.columns !== other.rows) { + // eslint-disable-next-line no-console + console.warn( + 'Number of columns of left matrix are not equal to number of rows of right matrix.', + ); + } + + const m = this.rows; + const p = other.columns; + + const result = new SparseMatrix(m, p); + this.forEachNonZero((i, j, v1) => { + other.forEachNonZero((k, l, v2) => { + if (j === k) { + result.set(i, l, result.get(i, l) + v1 * v2); + } + return v2; + }); + return v1; + }); + return result; + } + + /** + * @param {SparseMatrix} other + * @returns {SparseMatrix} + */ + kroneckerProduct(other) { + const m = this.rows; + const n = this.columns; + const p = other.rows; + const q = other.columns; + + const result = new SparseMatrix(m * p, n * q, { + initialCapacity: this.cardinality * other.cardinality, + }); + this.forEachNonZero((i, j, v1) => { + other.forEachNonZero((k, l, v2) => { + result.set(p * i + k, q * j + l, v1 * v2); + return v2; + }); + return v1; + }); + return result; + } + + /** + * Calls `callback` for each value in the matrix that is not zero. + * The callback can return: + * - `false` to stop the iteration immediately + * - the same value to keep it unchanged + * - 0 to remove the value + * - something else to change the value + * @param {ForEachNonZeroCallback} callback + * @returns {this} + */ + forEachNonZero(callback) { + this.elements.forEachPair((key, value) => { + const i = (key / this.columns) | 0; + const j = key % this.columns; + let r = callback(i, j, value); + if (r === false) return false; // stop iteration + if (this.threshold && Math.abs(r) < this.threshold) r = 0; + if (r !== value) { + if (r === 0) { + this.elements.remove(key, true); + } else { + this.elements.set(key, r); + } + } + return true; + }); + this.elements.maybeShrinkCapacity(); + return this; + } + + getNonZeros() { + const cardinality = this.cardinality; + /** @type {number[]} */ + const rows = new Array(cardinality); + /** @type {number[]} */ + const columns = new Array(cardinality); + /** @type {number[]} */ + const values = new Array(cardinality); + let idx = 0; + this.forEachNonZero((i, j, value) => { + rows[idx] = i; + columns[idx] = j; + values[idx] = value; + idx++; + return value; + }); + return { rows, columns, values }; + } + + /** + * @param {number} newThreshold + * @returns {this} + */ + setThreshold(newThreshold) { + if (newThreshold !== 0 && newThreshold !== this.threshold) { + this.threshold = newThreshold; + this.forEachNonZero((i, j, v) => v); + } + return this; + } + + /** + * @returns {SparseMatrix} - New transposed sparse matrix + */ + transpose() { + let trans = new SparseMatrix(this.columns, this.rows, { + initialCapacity: this.cardinality, + }); + this.forEachNonZero((i, j, value) => { + trans.set(j, i, value); + return value; + }); + return trans; + } + + /** + * @returns {boolean} + */ + isEmpty() { + return this.rows === 0 || this.columns === 0; + } + + // The following was generated by `makeMathMethods.mjs` + //#region Math methods. + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + add(value) { + if (typeof value === 'number') return this.addS(value); + return this.addM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + addS(value) { + this.forEachNonZero((i, j, v) => v + value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + addM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) + v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static add(matrix, value) { + return new SparseMatrix(matrix).addS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + sub(value) { + if (typeof value === 'number') return this.subS(value); + return this.subM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + subS(value) { + this.forEachNonZero((i, j, v) => v - value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + subM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) - v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static sub(matrix, value) { + return new SparseMatrix(matrix).subS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + subtract(value) { + if (typeof value === 'number') return this.subtractS(value); + return this.subtractM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + subtractS(value) { + this.forEachNonZero((i, j, v) => v - value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + subtractM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) - v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static subtract(matrix, value) { + return new SparseMatrix(matrix).subtractS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + mul(value) { + if (typeof value === 'number') return this.mulS(value); + return this.mulM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + mulS(value) { + this.forEachNonZero((i, j, v) => v * value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + mulM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) * v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static mul(matrix, value) { + return new SparseMatrix(matrix).mulS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + multiply(value) { + if (typeof value === 'number') return this.multiplyS(value); + return this.multiplyM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + multiplyS(value) { + this.forEachNonZero((i, j, v) => v * value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + multiplyM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) * v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static multiply(matrix, value) { + return new SparseMatrix(matrix).multiplyS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + div(value) { + if (typeof value === 'number') return this.divS(value); + return this.divM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + divS(value) { + this.forEachNonZero((i, j, v) => v / value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + divM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) / v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static div(matrix, value) { + return new SparseMatrix(matrix).divS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + divide(value) { + if (typeof value === 'number') return this.divideS(value); + return this.divideM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + divideS(value) { + this.forEachNonZero((i, j, v) => v / value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + divideM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) / v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static divide(matrix, value) { + return new SparseMatrix(matrix).divideS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + mod(value) { + if (typeof value === 'number') return this.modS(value); + return this.modM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + modS(value) { + this.forEachNonZero((i, j, v) => v % value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + modM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) % v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static mod(matrix, value) { + return new SparseMatrix(matrix).modS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + modulus(value) { + if (typeof value === 'number') return this.modulusS(value); + return this.modulusM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + modulusS(value) { + this.forEachNonZero((i, j, v) => v % value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + modulusM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) % v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static modulus(matrix, value) { + return new SparseMatrix(matrix).modulusS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + and(value) { + if (typeof value === 'number') return this.andS(value); + return this.andM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + andS(value) { + this.forEachNonZero((i, j, v) => v & value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + andM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) & v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static and(matrix, value) { + return new SparseMatrix(matrix).andS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + or(value) { + if (typeof value === 'number') return this.orS(value); + return this.orM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + orS(value) { + this.forEachNonZero((i, j, v) => v | value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + orM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) | v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static or(matrix, value) { + return new SparseMatrix(matrix).orS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + xor(value) { + if (typeof value === 'number') return this.xorS(value); + return this.xorM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + xorS(value) { + this.forEachNonZero((i, j, v) => v ^ value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + xorM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) ^ v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static xor(matrix, value) { + return new SparseMatrix(matrix).xorS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + leftShift(value) { + if (typeof value === 'number') return this.leftShiftS(value); + return this.leftShiftM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + leftShiftS(value) { + this.forEachNonZero((i, j, v) => v << value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + leftShiftM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) << v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static leftShift(matrix, value) { + return new SparseMatrix(matrix).leftShiftS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + signPropagatingRightShift(value) { + if (typeof value === 'number') { + return this.signPropagatingRightShiftS(value); + } + return this.signPropagatingRightShiftM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + signPropagatingRightShiftS(value) { + this.forEachNonZero((i, j, v) => v >> value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + signPropagatingRightShiftM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) >> v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static signPropagatingRightShift(matrix, value) { + return new SparseMatrix(matrix).signPropagatingRightShiftS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + rightShift(value) { + if (typeof value === 'number') return this.rightShiftS(value); + return this.rightShiftM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + rightShiftS(value) { + this.forEachNonZero((i, j, v) => v >>> value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + rightShiftM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) >>> v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static rightShift(matrix, value) { + return new SparseMatrix(matrix).rightShiftS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + zeroFillRightShift(value) { + if (typeof value === 'number') return this.zeroFillRightShiftS(value); + return this.zeroFillRightShiftM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + zeroFillRightShiftS(value) { + this.forEachNonZero((i, j, v) => v >>> value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + zeroFillRightShiftM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) >>> v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static zeroFillRightShift(matrix, value) { + return new SparseMatrix(matrix).zeroFillRightShiftS(value); + } + + /** + * @returns {this} + */ + not() { + this.forEachNonZero((i, j, v) => ~v); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static not(matrix) { + return new SparseMatrix(matrix).not(); + } + + /** + * @returns {this} + */ + abs() { + this.forEachNonZero((i, j, v) => Math.abs(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static abs(matrix) { + return new SparseMatrix(matrix).abs(); + } + + /** + * @returns {this} + */ + acos() { + this.forEachNonZero((i, j, v) => Math.acos(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static acos(matrix) { + return new SparseMatrix(matrix).acos(); + } + + /** + * @returns {this} + */ + acosh() { + this.forEachNonZero((i, j, v) => Math.acosh(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static acosh(matrix) { + return new SparseMatrix(matrix).acosh(); + } + + /** + * @returns {this} + */ + asin() { + this.forEachNonZero((i, j, v) => Math.asin(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static asin(matrix) { + return new SparseMatrix(matrix).asin(); + } + + /** + * @returns {this} + */ + asinh() { + this.forEachNonZero((i, j, v) => Math.asinh(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static asinh(matrix) { + return new SparseMatrix(matrix).asinh(); + } + + /** + * @returns {this} + */ + atan() { + this.forEachNonZero((i, j, v) => Math.atan(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static atan(matrix) { + return new SparseMatrix(matrix).atan(); + } + + /** + * @returns {this} + */ + atanh() { + this.forEachNonZero((i, j, v) => Math.atanh(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static atanh(matrix) { + return new SparseMatrix(matrix).atanh(); + } + + /** + * @returns {this} + */ + cbrt() { + this.forEachNonZero((i, j, v) => Math.cbrt(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static cbrt(matrix) { + return new SparseMatrix(matrix).cbrt(); + } + + /** + * @returns {this} + */ + ceil() { + this.forEachNonZero((i, j, v) => Math.ceil(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static ceil(matrix) { + return new SparseMatrix(matrix).ceil(); + } + + /** + * @returns {this} + */ + clz32() { + this.forEachNonZero((i, j, v) => Math.clz32(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static clz32(matrix) { + return new SparseMatrix(matrix).clz32(); + } + + /** + * @returns {this} + */ + cos() { + this.forEachNonZero((i, j, v) => Math.cos(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static cos(matrix) { + return new SparseMatrix(matrix).cos(); + } + + /** + * @returns {this} + */ + cosh() { + this.forEachNonZero((i, j, v) => Math.cosh(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static cosh(matrix) { + return new SparseMatrix(matrix).cosh(); + } + + /** + * @returns {this} + */ + exp() { + this.forEachNonZero((i, j, v) => Math.exp(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static exp(matrix) { + return new SparseMatrix(matrix).exp(); + } + + /** + * @returns {this} + */ + expm1() { + this.forEachNonZero((i, j, v) => Math.expm1(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static expm1(matrix) { + return new SparseMatrix(matrix).expm1(); + } + + /** + * @returns {this} + */ + floor() { + this.forEachNonZero((i, j, v) => Math.floor(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static floor(matrix) { + return new SparseMatrix(matrix).floor(); + } + + /** + * @returns {this} + */ + fround() { + this.forEachNonZero((i, j, v) => Math.fround(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static fround(matrix) { + return new SparseMatrix(matrix).fround(); + } + + /** + * @returns {this} + */ + log() { + this.forEachNonZero((i, j, v) => Math.log(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static log(matrix) { + return new SparseMatrix(matrix).log(); + } + + /** + * @returns {this} + */ + log1p() { + this.forEachNonZero((i, j, v) => Math.log1p(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static log1p(matrix) { + return new SparseMatrix(matrix).log1p(); + } + + /** + * @returns {this} + */ + log10() { + this.forEachNonZero((i, j, v) => Math.log10(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static log10(matrix) { + return new SparseMatrix(matrix).log10(); + } + + /** + * @returns {this} + */ + log2() { + this.forEachNonZero((i, j, v) => Math.log2(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static log2(matrix) { + return new SparseMatrix(matrix).log2(); + } + + /** + * @returns {this} + */ + round() { + this.forEachNonZero((i, j, v) => Math.round(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static round(matrix) { + return new SparseMatrix(matrix).round(); + } + + /** + * @returns {this} + */ + sign() { + this.forEachNonZero((i, j, v) => Math.sign(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static sign(matrix) { + return new SparseMatrix(matrix).sign(); + } + + /** + * @returns {this} + */ + sin() { + this.forEachNonZero((i, j, v) => Math.sin(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static sin(matrix) { + return new SparseMatrix(matrix).sin(); + } + + /** + * @returns {this} + */ + sinh() { + this.forEachNonZero((i, j, v) => Math.sinh(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static sinh(matrix) { + return new SparseMatrix(matrix).sinh(); + } + + /** + * @returns {this} + */ + sqrt() { + this.forEachNonZero((i, j, v) => Math.sqrt(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static sqrt(matrix) { + return new SparseMatrix(matrix).sqrt(); + } + + /** + * @returns {this} + */ + tan() { + this.forEachNonZero((i, j, v) => Math.tan(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static tan(matrix) { + return new SparseMatrix(matrix).tan(); + } + + /** + * @returns {this} + */ + tanh() { + this.forEachNonZero((i, j, v) => Math.tanh(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static tanh(matrix) { + return new SparseMatrix(matrix).tanh(); + } + + /** + * @returns {this} + */ + trunc() { + this.forEachNonZero((i, j, v) => Math.trunc(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static trunc(matrix) { + return new SparseMatrix(matrix).trunc(); + } + //#endregion +} + +SparseMatrix.prototype.klass = 'Matrix'; + +SparseMatrix.identity = SparseMatrix.eye; +SparseMatrix.prototype.tensorProduct = SparseMatrix.prototype.kroneckerProduct; From 7ac8a4c69e0b693c8b2ad6678cadd84a08a1f44c Mon Sep 17 00:00:00 2001 From: jobo322 Date: Sat, 21 Jun 2025 00:17:08 -0500 Subject: [PATCH 05/33] feat: add tests for getNonZeros method in SparseMatrix --- src/__tests__/getNonZeros.test.js | 37 +++++++++++++++++++++++++++++++ 1 file changed, 37 insertions(+) create mode 100644 src/__tests__/getNonZeros.test.js diff --git a/src/__tests__/getNonZeros.test.js b/src/__tests__/getNonZeros.test.js new file mode 100644 index 0000000..1b93740 --- /dev/null +++ b/src/__tests__/getNonZeros.test.js @@ -0,0 +1,37 @@ +import { SparseMatrix } from '../index'; + +describe('Sparse Matrix', () => { + it('getNonZeros', () => { + let m2 = new SparseMatrix( + [ + [0, 0, 0, 0, 0, 0], + [1, 0, 0, 2, 1, 1], + [0, 3, 0, 0, 5, 5], + [0, 0, 0, 0, 0, 0], + [1, 0, 0, 1, 9, 9], + ], + { initialCapacity: 12 }, + ); + + // Default (coordinate list) + expect(m2.getNonZeros()).toEqual({ + rows: Float64Array.from([4, 1, 4, 1, 1, 1, 4, 2, 4, 2, 2]), + columns: Float64Array.from([0, 0, 3, 3, 4, 5, 5, 1, 4, 4, 5]), + values: Float64Array.from([1, 1, 1, 2, 1, 1, 9, 3, 9, 5, 5]), + }); + + // CSR format + expect(m2.getNonZeros({ format: 'csr' })).toEqual({ + rows: Float64Array.from([0, 0, 4, 7, 7, 11]), + columns: Float64Array.from([0, 3, 4, 5, 1, 4, 5, 0, 3, 4, 5]), + values: Float64Array.from([1, 2, 1, 1, 3, 5, 5, 1, 1, 9, 9]), + }); + + //CSC format + expect(m2.getNonZeros({ format: 'csc' })).toEqual({ + rows: Float64Array.from([1, 4, 2, 1, 4, 1, 2, 4, 1, 2, 4]), + columns: Float64Array.from([0, 2, 3, 3, 5, 8, 11]), + values: Float64Array.from([1, 1, 3, 2, 1, 1, 5, 9, 1, 5, 9]), + }); + }); +}); From f2e2935e39414ec11394f520d679e05d607ff26f Mon Sep 17 00:00:00 2001 From: jobo322 Date: Sat, 21 Jun 2025 01:53:21 -0500 Subject: [PATCH 06/33] feat: update benchmark parameters and add comprehensive tests for SparseMatrix operations --- benchmark/benchmark.js | 17 ++- package.json | 3 +- src/__tests__/{test.js => matrix.test.js} | 52 +++++++ src/index.js | 176 ++++++++++++---------- 4 files changed, 160 insertions(+), 88 deletions(-) rename src/__tests__/{test.js => matrix.test.js} (67%) diff --git a/benchmark/benchmark.js b/benchmark/benchmark.js index bb1eac5..2557758 100644 --- a/benchmark/benchmark.js +++ b/benchmark/benchmark.js @@ -54,10 +54,12 @@ function printWinner(label1, avg1, label2, avg2) { function runBenchmarks() { const m = 128; - const n = 64; + const n = 128; const p = 128; - const densityA = 0.03; - const densityB = 0.03; + const densityA = 0.01; + const densityB = 0.01; + + const nbIterations = 3; const A = randomSparseMatrix(m, n, densityA); const B = randomSparseMatrix(n, p, densityB); @@ -66,10 +68,11 @@ function runBenchmarks() { const AOld = new SparseMatrixOld(denseA); const BOld = new SparseMatrixOld(denseB); - A.cardinality; + denseA = new Matrix(denseA); denseB = new Matrix(denseB); + denseA.mmul(denseB); // 1. add vs addNew // const addAvg = benchmark(() => { // const a = AOld.clone(); @@ -90,7 +93,7 @@ function runBenchmarks() { A.mmul(B); }, 'mmulNew', - 3, + nbIterations, ); const mmulAvg = benchmark( @@ -98,7 +101,7 @@ function runBenchmarks() { AOld.mmul(BOld); }, 'mmul', - 3, + nbIterations, ); const denseAvg = benchmark( @@ -106,7 +109,7 @@ function runBenchmarks() { denseA.mmul(denseB); }, 'denseMatrix', - 3, + nbIterations, ); printWinner('mmul', mmulAvg, 'mmulNew', mmulNewAvg); diff --git a/package.json b/package.json index a86a279..4501c71 100644 --- a/package.json +++ b/package.json @@ -12,7 +12,7 @@ "compile": "rollup -c", "eslint": "eslint src", "eslint-fix": "npm run eslint -- --fix", - "prepack": "npm run compile", + "prepack": "npm run compile", "prettier": "prettier --check src", "prettier-write": "prettier --write src", "test": "npm run test-only && npm run eslint && npm run prettier", @@ -37,6 +37,7 @@ "eslint": "^8.10.0", "eslint-config-cheminfo": "^7.2.2", "jest": "^27.5.1", + "ml-matrix": "^6.12.1", "prettier": "^2.5.1", "rollup": "^2.69.0" } diff --git a/src/__tests__/test.js b/src/__tests__/matrix.test.js similarity index 67% rename from src/__tests__/test.js rename to src/__tests__/matrix.test.js index 0bfca72..1885488 100644 --- a/src/__tests__/test.js +++ b/src/__tests__/matrix.test.js @@ -40,6 +40,29 @@ describe('Sparse Matrix', () => { [0, 2], [0, 0], ]); + + // Compare with dense multiplication + const denseM1 = m1.to2DArray(); + const denseM2 = m2.to2DArray(); + const expectedDense = denseMatrixMultiply(denseM1, denseM2); + expect(m3.to2DArray()).toStrictEqual(expectedDense); + }); + + it('mmul', () => { + const size = 32; + const density = 0.1; + const m1 = randomSparseMatrix(size, size, density); + const m2 = randomSparseMatrix(size, size, density); + let m3 = m1.mmul(m2); + + const denseM1 = m1.to2DArray(); + const denseM2 = m2.to2DArray(); + + const newSparse = new SparseMatrix(denseM1); + expect(newSparse.to2DArray()).toStrictEqual(denseM1); + const expectedDense = denseMatrixMultiply(denseM1, denseM2); + + expect(m3.to2DArray()).toStrictEqual(expectedDense); }); it('kronecker', () => { @@ -142,3 +165,32 @@ describe('Banded matrices', () => { expect(matrix4.isBanded(1)).toBe(true); }); }); + +function denseMatrixMultiply(A, B) { + const rowsA = A.length; + const colsA = A[0].length; + const colsB = B[0].length; + const result = Array.from({ length: rowsA }, () => Array(colsB).fill(0)); + for (let i = 0; i < rowsA; i++) { + for (let j = 0; j < colsB; j++) { + for (let k = 0; k < colsA; k++) { + result[i][j] += A[i][k] * B[k][j]; + } + } + } + return result; +} + +function randomSparseMatrix(rows, cols, density = 0.01) { + const matrix = []; + for (let i = 0; i < rows; i++) { + const row = new Float64Array(cols); + for (let j = 0; j < cols; j++) { + if (Math.random() < density) { + row[j] = Math.random() * 10; + } + } + matrix.push(row); + } + return new SparseMatrix(matrix); +} diff --git a/src/index.js b/src/index.js index 90f185a..e55139c 100644 --- a/src/index.js +++ b/src/index.js @@ -148,21 +148,25 @@ export class SparseMatrix { return this; } - mulNew(other) { + mul(other) { if (typeof other !== 'number') { throw new RangeError('the argument should be a number'); } - // if (this.cardinality / this.columns / this.rows > 0.1) + if (other === 0) { + return new SparseMatrix(this.rows, this.columns, { + initialCapacity: this.cardinality, + }); + } + this.withEachNonZero((i, j, v) => { - // return v * other; this.set(i, j, v * other); }); return this; } - mmulNew(other) { + mmul(other) { if (this.columns !== other.rows) { // eslint-disable-next-line no-console console.warn( @@ -173,80 +177,51 @@ export class SparseMatrix { const m = this.rows; const p = other.columns; + const result = matrixCreateEmpty(m, p); + const useCscCsr = useCscCsrFormat(this, other); + const { columns: otherCols, rows: otherRows, values: otherValues, - } = other.getNonZeros(); + } = other.getNonZeros(useCscCsr ? { format: 'csr' } : {}); const { columns: thisCols, rows: thisRows, values: thisValues, - } = this.getNonZeros(); - - const result = new SparseMatrix(m, p, { initialCapacity: m * p }); - - const nbOtherActive = otherCols.length; - const nbThisActive = thisCols.length; - for (let t = 0; t < nbThisActive; t++) { - const i = thisRows[t]; - const j = thisCols[t]; - for (let o = 0; o < nbOtherActive; o++) { - if (j === otherRows[o]) { - const l = otherCols[o]; - result.set(i, l, result.get(i, l) + otherValues[o] * thisValues[t]); + } = this.getNonZeros(useCscCsr ? { format: 'csc' } : {}); + + if (useCscCsr) { + const thisNbCols = this.columns; + for (let t = 0; t < thisNbCols; t++) { + const tValues = thisValues.subarray(thisCols[t], thisCols[t + 1]); + const tRows = thisRows.subarray(thisCols[t], thisCols[t + 1]); + const oValues = otherValues.subarray(otherRows[t], otherRows[t + 1]); + const oCols = otherCols.subarray(otherRows[t], otherRows[t + 1]); + for (let f = 0; f < tValues.length; f++) { + for (let k = 0; k < oValues.length; k++) { + const i = tRows[f]; + const l = oCols[k]; + result[i][l] += tValues[f] * oValues[k]; + } } } - } - - return result; - } - - mmul(other) { - if (this.columns !== other.rows) { - // eslint-disable-next-line no-console - console.warn( - 'Number of columns of left matrix are not equal to number of rows of right matrix.', - ); - } - - const m = this.rows; - const p = other.columns; - - const { - columns: otherCols, - rows: otherRows, - values: otherValues, - } = other.getNonZeros('csr'); - const { - columns: thisCols, - rows: thisRows, - values: thisValues, - } = this.getNonZeros('csc'); - - const result = new SparseMatrix(m, p, { initialCapacity: m * p + 20 }); - - for (let t = 0; t < thisCols.length - 1; t++) { - const j = t; - const tValues = thisValues.subarray(thisCols[t], thisCols[t + 1]); - const tRows = thisRows.subarray(thisCols[t], thisCols[t + 1]); - let initOther = 0; - for (let o = initOther; o < otherRows.length - 1; o++) { - if (o === j) { - initOther++; - const oValues = otherValues.subarray(otherRows[o], otherRows[o + 1]); - const oCols = otherCols.subarray(otherRows[o], otherRows[o + 1]); - for (let f = 0; f < tValues.length; f++) { - for (let k = 0; k < oValues.length; k++) { - const i = tRows[f]; - const l = oCols[k]; - result.set(i, l, result.get(i, l) + tValues[f] * oValues[k]); - } + } else { + const nbOtherActive = otherCols.length; + const nbThisActive = thisCols.length; + for (let t = 0; t < nbThisActive; t++) { + const i = thisRows[t]; + const j = thisCols[t]; + for (let o = 0; o < nbOtherActive; o++) { + if (j === otherRows[o]) { + const l = otherCols[o]; + result[i][l] += otherValues[o] * thisValues[t]; } - } else if (j < o) break; + } } } - return result; + + return new SparseMatrix(result); } kroneckerProduct(other) { @@ -270,8 +245,8 @@ export class SparseMatrix { values: thisValues, } = this.getNonZeros(); - const nbOtherActive = otherCols.length; const nbThisActive = thisCols.length; + const nbOtherActive = otherCols.length; for (let t = 0; t < nbThisActive; t++) { const pi = p * thisRows[t]; const qj = q * thisCols[t]; @@ -287,13 +262,18 @@ export class SparseMatrix { return result; } - withEachNonZero(callback) { + withEachNonZero(callback, sort = false) { const { state, table, values } = this.elements; const nbStates = state.length; - const activeIndex = []; - for (let i = 0; i < nbStates; i++) { - if (state[i] === 1) activeIndex.push(i); + const activeIndex = new Float64Array(this.cardinality); + for (let i = 0, j = 0; i < nbStates; i++) { + if (state[i] === 1) activeIndex[j++] = i; } + + if (sort) { + activeIndex.sort((a, b) => table[a] - table[b]); + } + const columns = this.columns; for (const i of activeIndex) { const key = table[i]; @@ -321,26 +301,34 @@ export class SparseMatrix { return this; } - //'csr' | 'csc' | undefined - getNonZeros(format) { + getNonZeros(options = {}) { const cardinality = this.cardinality; const rows = new Float64Array(cardinality); const columns = new Float64Array(cardinality); const values = new Float64Array(cardinality); + + const { format, sort = false } = options; + let idx = 0; this.withEachNonZero((i, j, value) => { rows[idx] = i; columns[idx] = j; values[idx] = value; idx++; - }); + }, sort || format); const cooMatrix = { rows, columns, values }; if (!format) return cooMatrix; + if (!['csr', 'csc'].includes(format.toLowerCase())) { + throw new Error(`format ${format} is not supported`); + } + const csrMatrix = cooToCsr(cooMatrix, this.rows); - return format === 'csc' ? csrToCsc(csrMatrix, this.columns) : csrMatrix; + return format.toLowerCase() === 'csc' + ? csrToCsc(csrMatrix, this.columns) + : csrMatrix; } setThreshold(newThreshold) { @@ -528,25 +516,21 @@ function csrToCsc(csrMatrix, numCols) { columns: csrColIndices, rows: csrRowPtr, } = csrMatrix; - // Initialize CSC arrays + const cscValues = new Float64Array(csrValues.length); const cscRowIndices = new Float64Array(csrValues.length); const cscColPtr = new Float64Array(numCols + 1); - // Count non-zeros per column for (let i = 0; i < csrColIndices.length; i++) { cscColPtr[csrColIndices[i] + 1]++; } - // Compute column pointers (prefix sum) for (let i = 1; i <= numCols; i++) { cscColPtr[i] += cscColPtr[i - 1]; } - // Temporary copy for filling values const next = cscColPtr.slice(); - // Fill CSC arrays for (let row = 0; row < csrRowPtr.length - 1; row++) { for (let j = csrRowPtr[row], i = 0; j < csrRowPtr[row + 1]; j++, i++) { const col = csrColIndices[j]; @@ -555,7 +539,6 @@ function csrToCsc(csrMatrix, numCols) { cscRowIndices[pos] = row; next[col]++; } - // if (row === 1) break; } return { rows: cscRowIndices, columns: cscColPtr, values: cscValues }; @@ -575,3 +558,36 @@ function cooToCsr(cooMatrix, nbRows = 9) { } return { rows: csrRowPtr, columns, values }; } + +function matrixCreateEmpty(nbRows, nbColumns) { + const newMatrix = []; + for (let row = 0; row < nbRows; row++) { + newMatrix.push(new Float64Array(nbColumns)); + } + return newMatrix; +} + +function useCscCsrFormat(A, B) { + // Validar dimensiones + if (!A.rows || !A.columns || !B.rows || !B.columns || A.columns !== B.rows) { + throw new Error('Dimensiones de matrices inválidas o no compatibles.'); + } + + const densityA = A.cardinality / (A.rows * A.columns); + const densityB = B.cardinality / (B.rows * B.columns); + + const M = A.rows; + const K = A.columns; + const N = B.columns; + + const isHighDensity = densityA >= 0.015 || densityB >= 0.015; + const isVeryHighDensity = densityB >= 0.03; + const isLargeMatrix = M >= 64 || K >= 64 || N >= 64; + const isLargeSquareMatrix = M === K && K === N && M >= 128; + + return ( + isVeryHighDensity || + (isHighDensity && isLargeMatrix) || + (isLargeSquareMatrix && densityA >= 0.03) + ); +} From 099f24f089a29ccd1c6393582c7a1ed9d27d9e5e Mon Sep 17 00:00:00 2001 From: jobo322 Date: Fri, 13 Jun 2025 15:04:54 -0500 Subject: [PATCH 07/33] feat: add new matrix addition method and improve performance of mmul and kroneckerProduct method --- src/__tests__/index.test.js | 18 ++++++++++++++++ src/index.js | 41 +++++++++++++++++++++++++------------ 2 files changed, 46 insertions(+), 13 deletions(-) diff --git a/src/__tests__/index.test.js b/src/__tests__/index.test.js index 0843eba..d53f522 100644 --- a/src/__tests__/index.test.js +++ b/src/__tests__/index.test.js @@ -3,6 +3,24 @@ import { describe, expect, it } from 'vitest'; import { SparseMatrix } from '../index.js'; describe('Sparse Matrix', () => { + it('add', () => { + let m1 = new SparseMatrix([ + [2, 0, 1], + [0, 0, 3], + [2, 0, 1], + ]); + let m2 = new SparseMatrix([ + [0, 1, 5], + [2, 0, 0], + [-2, 0, -1], + ]); + let m3 = m1.add(m2).to2DArray(); + expect(m3).toStrictEqual([ + [2, 1, 6], + [2, 0, 3], + [0, 0, 0], + ]); + }); it('mmul', () => { let m1 = new SparseMatrix([ [2, 0, 1], diff --git a/src/index.js b/src/index.js index f770b91..a1f16f3 100644 --- a/src/index.js +++ b/src/index.js @@ -18,11 +18,12 @@ export class SparseMatrix { if (Array.isArray(rows)) { const matrix = rows; - rows = matrix.length; - options = columns || {}; + const nbRows = matrix.length; + const nbColmuns = matrix[0].length; + options = columns || { initialCapacity: nbRows * nbColmuns }; columns = matrix[0].length; - this._init(rows, columns, new HashTable(options), options.threshold); - for (let i = 0; i < rows; i++) { + this._init(nbRows, columns, new HashTable(options), options.threshold); + for (let i = 0; i < nbRows; i++) { for (let j = 0; j < columns; j++) { let value = matrix[i][j]; if (this.threshold && Math.abs(value) < this.threshold) value = 0; @@ -31,6 +32,7 @@ export class SparseMatrix { } } } + this.elements.maybeShrinkCapacity(); } else { this._init(rows, columns, new HashTable(options), options.threshold); } @@ -169,14 +171,12 @@ export class SparseMatrix { const p = other.columns; const result = new SparseMatrix(m, p); - this.forEachNonZero((i, j, v1) => { - other.forEachNonZero((k, l, v2) => { + this.withEachNonZero((i, j, v1) => { + other.withEachNonZero((k, l, v2) => { if (j === k) { result.set(i, l, result.get(i, l) + v1 * v2); } - return v2; }); - return v1; }); return result; } @@ -194,16 +194,31 @@ export class SparseMatrix { const result = new SparseMatrix(m * p, n * q, { initialCapacity: this.cardinality * other.cardinality, }); - this.forEachNonZero((i, j, v1) => { - other.forEachNonZero((k, l, v2) => { - result.set(p * i + k, q * j + l, v1 * v2); - return v2; + + this.withEachNonZero((i, j, v1) => { + const pi = p * i; + const qj = q * j; + other.withEachNonZero((k, l, v2) => { + result.set(pi + k, qj + l, v1 * v2); }); - return v1; }); return result; } + withEachNonZero(callback) { + const { state, table, values } = this.elements; + const nbStates = state.length; + const activeIndex = []; + for (let i = 0; i < nbStates; i++) { + if (state[i] === 1) activeIndex.push(i); + } + const columns = this.columns; + for (const i of activeIndex) { + const key = table[i]; + callback((key / columns) | 0, key % columns, values[i]); + } + } + /** * Calls `callback` for each value in the matrix that is not zero. * The callback can return: From 104e3e135a83b322d014b8e622c3d93cdc94a253 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Fri, 13 Jun 2025 16:09:36 -0500 Subject: [PATCH 08/33] feat: improve performance of mmul and kronecker product --- src/__tests__/index.test.js | 4 +++ src/index.js | 68 +++++++++++++++++++++++++++---------- 2 files changed, 55 insertions(+), 17 deletions(-) diff --git a/src/__tests__/index.test.js b/src/__tests__/index.test.js index d53f522..d0c1caf 100644 --- a/src/__tests__/index.test.js +++ b/src/__tests__/index.test.js @@ -38,6 +38,10 @@ describe('Sparse Matrix', () => { expect(m3.cardinality).toBe(1); expect(m3.get(0, 1)).toBe(2); + expect(m3.to2DArray()).toStrictEqual([ + [0, 2], + [0, 0], + ]); }); it('kronecker', () => { diff --git a/src/index.js b/src/index.js index a1f16f3..c26ee3a 100644 --- a/src/index.js +++ b/src/index.js @@ -170,14 +170,32 @@ export class SparseMatrix { const m = this.rows; const p = other.columns; + const { + columns: otherCols, + rows: otherRows, + values: otherValues, + } = other.getNonZeros(); + const { + columns: thisCols, + rows: thisRows, + values: thisValues, + } = this.getNonZeros(); + const result = new SparseMatrix(m, p); - this.withEachNonZero((i, j, v1) => { - other.withEachNonZero((k, l, v2) => { - if (j === k) { - result.set(i, l, result.get(i, l) + v1 * v2); + + const nbOtherActive = otherCols.length; + const nbThisActive = thisCols.length; + for (let t = 0; t < nbThisActive; t++) { + const i = thisRows[t]; + const j = thisCols[t]; + for (let o = 0; o < nbOtherActive; o++) { + if (j === otherRows[o]) { + const l = otherCols[o]; + result.set(i, l, result.get(i, l) + otherValues[o] * thisValues[t]); } - }); - }); + } + } + return result; } @@ -195,13 +213,31 @@ export class SparseMatrix { initialCapacity: this.cardinality * other.cardinality, }); - this.withEachNonZero((i, j, v1) => { - const pi = p * i; - const qj = q * j; - other.withEachNonZero((k, l, v2) => { - result.set(pi + k, qj + l, v1 * v2); - }); - }); + const { + columns: otherCols, + rows: otherRows, + values: otherValues, + } = other.getNonZeros(); + const { + columns: thisCols, + rows: thisRows, + values: thisValues, + } = this.getNonZeros(); + + const nbOtherActive = otherCols.length; + const nbThisActive = thisCols.length; + for (let t = 0; t < nbThisActive; t++) { + const pi = p * thisRows[t]; + const qj = q * thisCols[t]; + for (let o = 0; o < nbOtherActive; o++) { + result.set( + pi + otherRows[o], + qj + otherCols[o], + otherValues[o] * thisValues[t], + ); + } + } + return result; } @@ -258,12 +294,11 @@ export class SparseMatrix { /** @type {number[]} */ const values = new Array(cardinality); let idx = 0; - this.forEachNonZero((i, j, value) => { + this.withEachNonZero((i, j, value) => { rows[idx] = i; columns[idx] = j; values[idx] = value; idx++; - return value; }); return { rows, columns, values }; } @@ -287,9 +322,8 @@ export class SparseMatrix { let trans = new SparseMatrix(this.columns, this.rows, { initialCapacity: this.cardinality, }); - this.forEachNonZero((i, j, value) => { + this.withEachNonZero((i, j, value) => { trans.set(j, i, value); - return value; }); return trans; } From 1bacd595f1acbfd71cb1067beacc8b8bbbdaa9a9 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Thu, 19 Jun 2025 23:27:15 -0500 Subject: [PATCH 09/33] feat: add new multiplication methods and improve getNonZeros function --- src/index.js | 146 +++++++++++++++++++++++++++++++++++++++++++++++---- 1 file changed, 135 insertions(+), 11 deletions(-) diff --git a/src/index.js b/src/index.js index c26ee3a..f4fe6ce 100644 --- a/src/index.js +++ b/src/index.js @@ -152,6 +152,7 @@ export class SparseMatrix { } else { this.elements.set(row * this.columns + column, value); } + return this; } @@ -159,7 +160,21 @@ export class SparseMatrix { * @param {SparseMatrix} other * @returns {SparseMatrix} */ - mmul(other) { + mulNew(other) { + if (typeof other !== 'number') { + throw new RangeError('the argument should be a number'); + } + + // if (this.cardinality / this.columns / this.rows > 0.1) + this.withEachNonZero((i, j, v) => { + // return v * other; + this.set(i, j, v * other); + }); + + return this; + } + + mmulNew(other) { if (this.columns !== other.rows) { // eslint-disable-next-line no-console console.warn( @@ -181,7 +196,7 @@ export class SparseMatrix { values: thisValues, } = this.getNonZeros(); - const result = new SparseMatrix(m, p); + const result = new SparseMatrix(m, p, { initialCapacity: m * p }); const nbOtherActive = otherCols.length; const nbThisActive = thisCols.length; @@ -199,6 +214,53 @@ export class SparseMatrix { return result; } + mmul(other) { + if (this.columns !== other.rows) { + // eslint-disable-next-line no-console + console.warn( + 'Number of columns of left matrix are not equal to number of rows of right matrix.', + ); + } + + const m = this.rows; + const p = other.columns; + + const { + columns: otherCols, + rows: otherRows, + values: otherValues, + } = other.getNonZeros('csr'); + const { + columns: thisCols, + rows: thisRows, + values: thisValues, + } = this.getNonZeros('csc'); + + const result = new SparseMatrix(m, p, { initialCapacity: m * p + 20 }); + + for (let t = 0; t < thisCols.length - 1; t++) { + const j = t; + const tValues = thisValues.subarray(thisCols[t], thisCols[t + 1]); + const tRows = thisRows.subarray(thisCols[t], thisCols[t + 1]); + let initOther = 0; + for (let o = initOther; o < otherRows.length - 1; o++) { + if (o === j) { + initOther++; + const oValues = otherValues.subarray(otherRows[o], otherRows[o + 1]); + const oCols = otherCols.subarray(otherRows[o], otherRows[o + 1]); + for (let f = 0; f < tValues.length; f++) { + for (let k = 0; k < oValues.length; k++) { + const i = tRows[f]; + const l = oCols[k]; + result.set(i, l, result.get(i, l) + tValues[f] * oValues[k]); + } + } + } else if (j < o) break; + } + } + return result; + } + /** * @param {SparseMatrix} other * @returns {SparseMatrix} @@ -210,7 +272,7 @@ export class SparseMatrix { const q = other.columns; const result = new SparseMatrix(m * p, n * q, { - initialCapacity: this.cardinality * other.cardinality, + initialCapacity: this.cardinality * other.cardinality + 20, }); const { @@ -285,14 +347,16 @@ export class SparseMatrix { return this; } - getNonZeros() { + //'csr' | 'csc' | undefined + getNonZeros(format) { const cardinality = this.cardinality; - /** @type {number[]} */ - const rows = new Array(cardinality); - /** @type {number[]} */ - const columns = new Array(cardinality); - /** @type {number[]} */ - const values = new Array(cardinality); + /** @type {Float64Array} */ + const rows = new Float64Array(cardinality); + /** @type {Float64Array} */ + const columns = new Float64Array(cardinality); + /** @type {Float64Array} */ + const values = new Float64Array(cardinality); + let idx = 0; this.withEachNonZero((i, j, value) => { rows[idx] = i; @@ -300,7 +364,13 @@ export class SparseMatrix { values[idx] = value; idx++; }); - return { rows, columns, values }; + + const cooMatrix = { rows, columns, values }; + + if (!format) return cooMatrix; + + const csrMatrix = cooToCsr(cooMatrix, this.rows); + return format === 'csc' ? csrToCsc(csrMatrix, this.columns) : csrMatrix; } /** @@ -1433,3 +1503,57 @@ SparseMatrix.prototype.klass = 'Matrix'; SparseMatrix.identity = SparseMatrix.eye; SparseMatrix.prototype.tensorProduct = SparseMatrix.prototype.kroneckerProduct; + +function csrToCsc(csrMatrix, numCols) { + const { + values: csrValues, + columns: csrColIndices, + rows: csrRowPtr, + } = csrMatrix; + // Initialize CSC arrays + const cscValues = new Float64Array(csrValues.length); + const cscRowIndices = new Float64Array(csrValues.length); + const cscColPtr = new Float64Array(numCols + 1); + + // Count non-zeros per column + for (let i = 0; i < csrColIndices.length; i++) { + cscColPtr[csrColIndices[i] + 1]++; + } + + // Compute column pointers (prefix sum) + for (let i = 1; i <= numCols; i++) { + cscColPtr[i] += cscColPtr[i - 1]; + } + + // Temporary copy for filling values + const next = cscColPtr.slice(); + + // Fill CSC arrays + for (let row = 0; row < csrRowPtr.length - 1; row++) { + for (let j = csrRowPtr[row], i = 0; j < csrRowPtr[row + 1]; j++, i++) { + const col = csrColIndices[j]; + const pos = next[col]; + cscValues[pos] = csrValues[j]; + cscRowIndices[pos] = row; + next[col]++; + } + // if (row === 1) break; + } + + return { rows: cscRowIndices, columns: cscColPtr, values: cscValues }; +} + +function cooToCsr(cooMatrix, nbRows = 9) { + const { values, columns, rows } = cooMatrix; + //could not be the same length + const csrRowPtr = new Float64Array(nbRows + 1); + const length = values.length; + let currentRow = rows[0]; + for (let index = 0; index < length; ) { + const prev = index; + while (currentRow === rows[index] && index < length) ++index; + csrRowPtr[currentRow + 1] = index; + currentRow += 1; + } + return { rows: csrRowPtr, columns, values }; +} From d647066e2039b98b25713967bf46b3d02930766f Mon Sep 17 00:00:00 2001 From: jobo322 Date: Sat, 21 Jun 2025 00:16:15 -0500 Subject: [PATCH 10/33] feat: add benchmark for sparse matrix operations and performance comparison --- benchmark/benchmark.js | 204 ++++ benchmark/class/SparseMatrixOld.js | 1386 ++++++++++++++++++++++++++++ 2 files changed, 1590 insertions(+) create mode 100644 benchmark/benchmark.js create mode 100644 benchmark/class/SparseMatrixOld.js diff --git a/benchmark/benchmark.js b/benchmark/benchmark.js new file mode 100644 index 0000000..bb1eac5 --- /dev/null +++ b/benchmark/benchmark.js @@ -0,0 +1,204 @@ +import { SparseMatrix } from '../src/index.js'; +import { Matrix } from 'ml-matrix'; +import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; +import fs from 'fs'; + +function randomSparseMatrix(rows, cols, density = 0.01) { + const matrix = []; + for (let i = 0; i < rows; i++) { + const row = new Array(cols).fill(0); + for (let j = 0; j < cols; j++) { + if (Math.random() < density) { + row[j] = Math.random() * 10; + } + } + matrix.push(row); + } + return new SparseMatrix(matrix); +} + +function benchmark(fn, label, iterations = 5) { + const times = []; + for (let i = 0; i < iterations; i++) { + const t0 = performance.now(); + fn(); + const t1 = performance.now(); + times.push(t1 - t0); + } + const avg = times.reduce((a, b) => a + b, 0) / times.length; + console.log(`${label}: avg ${avg.toFixed(2)} ms over ${iterations} runs`); + return avg; +} + +function printWinner(label1, avg1, label2, avg2) { + let winner, loser, win, lose; + if (avg1 < avg2) { + winner = label1; + win = avg1; + loser = label2; + lose = avg2; + } else { + winner = label2; + win = avg2; + loser = label1; + lose = avg1; + } + + const percent = ((lose - win) / lose) * 100; + console.log( + ` -> ${winner} was ${(lose / win).toFixed(2)}x faster (${percent.toFixed( + 1, + )}% faster) than ${loser}\n`, + ); +} + +function runBenchmarks() { + const m = 128; + const n = 64; + const p = 128; + const densityA = 0.03; + const densityB = 0.03; + const A = randomSparseMatrix(m, n, densityA); + const B = randomSparseMatrix(n, p, densityB); + + let denseA = A.to2DArray(); + let denseB = B.to2DArray(); + + const AOld = new SparseMatrixOld(denseA); + const BOld = new SparseMatrixOld(denseB); + A.cardinality; + denseA = new Matrix(denseA); + denseB = new Matrix(denseB); + + // 1. add vs addNew + // const addAvg = benchmark(() => { + // const a = AOld.clone(); + // a.add(BOld); + // }, 'add'); + + // const addNewAvg = benchmark(() => { + // const a = A.clone(); + // a.add(B); + // }, 'addNew'); + + // printWinner('add', addAvg, 'addNew', addNewAvg); + + // 2. mmul vs mmulNew + + const mmulNewAvg = benchmark( + () => { + A.mmul(B); + }, + 'mmulNew', + 3, + ); + + const mmulAvg = benchmark( + () => { + AOld.mmul(BOld); + }, + 'mmul', + 3, + ); + + const denseAvg = benchmark( + () => { + denseA.mmul(denseB); + }, + 'denseMatrix', + 3, + ); + + printWinner('mmul', mmulAvg, 'mmulNew', mmulNewAvg); + + // 3. kroneckerProduct vs kroneckerProductNew + // const kronNewAvg = benchmark(() => { + // A.kroneckerProduct(B); + // }, 'kroneckerProductNew'); + // const kronAvg = benchmark(() => { + // AOld.kroneckerProduct(BOld); + // }, 'kroneckerProduct'); + + // printWinner('kroneckerProduct', kronAvg, 'kroneckerProductNew', kronNewAvg); + + // 4. matrix multiplication + // const mulAvg = benchmark(() => { + // A.mul(5); + // }, 'mul'); + + // const mulNewAvg = benchmark(() => { + // AOld.mul(5); + // }, 'mulNew'); + + // printWinner('mul', mulAvg, 'mulNew', mulNewAvg); +} + +function runSizeSweepBenchmark() { + const sizes = [8, 16, 32, 64, 128]; + const densities = [0.01, 0.015, 0.02, 0.025, 0.03, 0.35]; + const results = []; + + for (const densityA of densities) { + for (const densityB of densities) { + for (const m of sizes) { + for (const n of sizes) { + for (const p of sizes) { + // A: m x n, B: n x p + + const A = randomSparseMatrix(m, n, densityA); + const B = randomSparseMatrix(n, p, densityB); + + let denseA = A.to2DArray(); + let denseB = B.to2DArray(); + + const AOld = new SparseMatrixOld(denseA); + const BOld = new SparseMatrixOld(denseB); + + denseA = new Matrix(denseA); + denseB = new Matrix(denseB); + + const mmulNewAvg = benchmark( + () => { + A.mmul(B); + }, + 'mmulNew', + 3, + ); + + const mmulAvg = benchmark( + () => { + AOld.mmul(BOld); + }, + 'mmul', + 3, + ); + + const denseAvg = benchmark(() => { + denseA.mmul(denseB), 'denseMatrix', 3; + }); + + results.push({ + densityA, + densityB, + A_shape: [m, n], + B_shape: [n, p], + dense: denseAvg, + mmulNew: mmulNewAvg, + mmul: mmulAvg, + }); + } + } + } + } + } + + fs.writeFileSync( + './benchmark/size_sweep_results-dense.json', + JSON.stringify(results, null, 2), + ); + console.log('Size sweep benchmark results saved to size_sweep_results.json'); +} + +runBenchmarks(); +// Uncomment to run the size sweep benchmark +// runSizeSweepBenchmark(); diff --git a/benchmark/class/SparseMatrixOld.js b/benchmark/class/SparseMatrixOld.js new file mode 100644 index 0000000..f770b91 --- /dev/null +++ b/benchmark/class/SparseMatrixOld.js @@ -0,0 +1,1386 @@ +import HashTable from 'ml-hash-table'; + +/** @typedef {(row: number, column: number, value: number) => number | false} ForEachNonZeroCallback */ + +export class SparseMatrix { + constructor(rows, columns, options = {}) { + if (rows instanceof SparseMatrix) { + // clone + const other = rows; + this._init( + other.rows, + other.columns, + other.elements.clone(), + other.threshold, + ); + return; + } + + if (Array.isArray(rows)) { + const matrix = rows; + rows = matrix.length; + options = columns || {}; + columns = matrix[0].length; + this._init(rows, columns, new HashTable(options), options.threshold); + for (let i = 0; i < rows; i++) { + for (let j = 0; j < columns; j++) { + let value = matrix[i][j]; + if (this.threshold && Math.abs(value) < this.threshold) value = 0; + if (value !== 0) { + this.elements.set(i * columns + j, matrix[i][j]); + } + } + } + } else { + this._init(rows, columns, new HashTable(options), options.threshold); + } + } + + _init(rows, columns, elements, threshold) { + this.rows = rows; + this.columns = columns; + this.elements = elements; + this.threshold = threshold || 0; + } + + static eye(rows = 1, columns = rows) { + const min = Math.min(rows, columns); + const matrix = new SparseMatrix(rows, columns, { initialCapacity: min }); + for (let i = 0; i < min; i++) { + matrix.set(i, i, 1); + } + return matrix; + } + + clone() { + return new SparseMatrix(this); + } + + /** + * @returns {number[][]} + */ + to2DArray() { + const copy = new Array(this.rows); + for (let i = 0; i < this.rows; i++) { + copy[i] = new Array(this.columns); + for (let j = 0; j < this.columns; j++) { + copy[i][j] = this.get(i, j); + } + } + return copy; + } + + isSquare() { + return this.rows === this.columns; + } + + isSymmetric() { + if (!this.isSquare()) return false; + + let symmetric = true; + this.forEachNonZero((i, j, v) => { + if (this.get(j, i) !== v) { + symmetric = false; + return false; + } + return v; + }); + return symmetric; + } + + /** + * Search for the wither band in the main diagonals + * @returns {number} + */ + bandWidth() { + let min = this.columns; + let max = -1; + this.forEachNonZero((i, j, v) => { + let diff = i - j; + min = Math.min(min, diff); + max = Math.max(max, diff); + return v; + }); + return max - min; + } + + /** + * Test if a matrix is consider banded using a threshold + * @param {number} width + * @returns {boolean} + */ + isBanded(width) { + let bandWidth = this.bandWidth(); + return bandWidth <= width; + } + + /** + * @returns {number} + */ + get cardinality() { + return this.elements.size; + } + + /** + * @returns {number} + */ + get size() { + return this.rows * this.columns; + } + + /** + * @param {number} row + * @param {number} column + * @returns {number} + */ + get(row, column) { + return this.elements.get(row * this.columns + column); + } + + /** + * @param {number} row + * @param {number} column + * @param {number} value + * @returns {this} + */ + set(row, column, value) { + if (this.threshold && Math.abs(value) < this.threshold) value = 0; + if (value === 0) { + this.elements.remove(row * this.columns + column); + } else { + this.elements.set(row * this.columns + column, value); + } + return this; + } + + /** + * @param {SparseMatrix} other + * @returns {SparseMatrix} + */ + mmul(other) { + if (this.columns !== other.rows) { + // eslint-disable-next-line no-console + console.warn( + 'Number of columns of left matrix are not equal to number of rows of right matrix.', + ); + } + + const m = this.rows; + const p = other.columns; + + const result = new SparseMatrix(m, p); + this.forEachNonZero((i, j, v1) => { + other.forEachNonZero((k, l, v2) => { + if (j === k) { + result.set(i, l, result.get(i, l) + v1 * v2); + } + return v2; + }); + return v1; + }); + return result; + } + + /** + * @param {SparseMatrix} other + * @returns {SparseMatrix} + */ + kroneckerProduct(other) { + const m = this.rows; + const n = this.columns; + const p = other.rows; + const q = other.columns; + + const result = new SparseMatrix(m * p, n * q, { + initialCapacity: this.cardinality * other.cardinality, + }); + this.forEachNonZero((i, j, v1) => { + other.forEachNonZero((k, l, v2) => { + result.set(p * i + k, q * j + l, v1 * v2); + return v2; + }); + return v1; + }); + return result; + } + + /** + * Calls `callback` for each value in the matrix that is not zero. + * The callback can return: + * - `false` to stop the iteration immediately + * - the same value to keep it unchanged + * - 0 to remove the value + * - something else to change the value + * @param {ForEachNonZeroCallback} callback + * @returns {this} + */ + forEachNonZero(callback) { + this.elements.forEachPair((key, value) => { + const i = (key / this.columns) | 0; + const j = key % this.columns; + let r = callback(i, j, value); + if (r === false) return false; // stop iteration + if (this.threshold && Math.abs(r) < this.threshold) r = 0; + if (r !== value) { + if (r === 0) { + this.elements.remove(key, true); + } else { + this.elements.set(key, r); + } + } + return true; + }); + this.elements.maybeShrinkCapacity(); + return this; + } + + getNonZeros() { + const cardinality = this.cardinality; + /** @type {number[]} */ + const rows = new Array(cardinality); + /** @type {number[]} */ + const columns = new Array(cardinality); + /** @type {number[]} */ + const values = new Array(cardinality); + let idx = 0; + this.forEachNonZero((i, j, value) => { + rows[idx] = i; + columns[idx] = j; + values[idx] = value; + idx++; + return value; + }); + return { rows, columns, values }; + } + + /** + * @param {number} newThreshold + * @returns {this} + */ + setThreshold(newThreshold) { + if (newThreshold !== 0 && newThreshold !== this.threshold) { + this.threshold = newThreshold; + this.forEachNonZero((i, j, v) => v); + } + return this; + } + + /** + * @returns {SparseMatrix} - New transposed sparse matrix + */ + transpose() { + let trans = new SparseMatrix(this.columns, this.rows, { + initialCapacity: this.cardinality, + }); + this.forEachNonZero((i, j, value) => { + trans.set(j, i, value); + return value; + }); + return trans; + } + + /** + * @returns {boolean} + */ + isEmpty() { + return this.rows === 0 || this.columns === 0; + } + + // The following was generated by `makeMathMethods.mjs` + //#region Math methods. + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + add(value) { + if (typeof value === 'number') return this.addS(value); + return this.addM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + addS(value) { + this.forEachNonZero((i, j, v) => v + value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + addM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) + v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static add(matrix, value) { + return new SparseMatrix(matrix).addS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + sub(value) { + if (typeof value === 'number') return this.subS(value); + return this.subM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + subS(value) { + this.forEachNonZero((i, j, v) => v - value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + subM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) - v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static sub(matrix, value) { + return new SparseMatrix(matrix).subS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + subtract(value) { + if (typeof value === 'number') return this.subtractS(value); + return this.subtractM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + subtractS(value) { + this.forEachNonZero((i, j, v) => v - value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + subtractM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) - v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static subtract(matrix, value) { + return new SparseMatrix(matrix).subtractS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + mul(value) { + if (typeof value === 'number') return this.mulS(value); + return this.mulM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + mulS(value) { + this.forEachNonZero((i, j, v) => v * value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + mulM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) * v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static mul(matrix, value) { + return new SparseMatrix(matrix).mulS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + multiply(value) { + if (typeof value === 'number') return this.multiplyS(value); + return this.multiplyM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + multiplyS(value) { + this.forEachNonZero((i, j, v) => v * value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + multiplyM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) * v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static multiply(matrix, value) { + return new SparseMatrix(matrix).multiplyS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + div(value) { + if (typeof value === 'number') return this.divS(value); + return this.divM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + divS(value) { + this.forEachNonZero((i, j, v) => v / value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + divM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) / v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static div(matrix, value) { + return new SparseMatrix(matrix).divS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + divide(value) { + if (typeof value === 'number') return this.divideS(value); + return this.divideM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + divideS(value) { + this.forEachNonZero((i, j, v) => v / value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + divideM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) / v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static divide(matrix, value) { + return new SparseMatrix(matrix).divideS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + mod(value) { + if (typeof value === 'number') return this.modS(value); + return this.modM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + modS(value) { + this.forEachNonZero((i, j, v) => v % value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + modM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) % v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static mod(matrix, value) { + return new SparseMatrix(matrix).modS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + modulus(value) { + if (typeof value === 'number') return this.modulusS(value); + return this.modulusM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + modulusS(value) { + this.forEachNonZero((i, j, v) => v % value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + modulusM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) % v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static modulus(matrix, value) { + return new SparseMatrix(matrix).modulusS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + and(value) { + if (typeof value === 'number') return this.andS(value); + return this.andM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + andS(value) { + this.forEachNonZero((i, j, v) => v & value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + andM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) & v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static and(matrix, value) { + return new SparseMatrix(matrix).andS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + or(value) { + if (typeof value === 'number') return this.orS(value); + return this.orM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + orS(value) { + this.forEachNonZero((i, j, v) => v | value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + orM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) | v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static or(matrix, value) { + return new SparseMatrix(matrix).orS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + xor(value) { + if (typeof value === 'number') return this.xorS(value); + return this.xorM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + xorS(value) { + this.forEachNonZero((i, j, v) => v ^ value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + xorM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) ^ v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static xor(matrix, value) { + return new SparseMatrix(matrix).xorS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + leftShift(value) { + if (typeof value === 'number') return this.leftShiftS(value); + return this.leftShiftM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + leftShiftS(value) { + this.forEachNonZero((i, j, v) => v << value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + leftShiftM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) << v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static leftShift(matrix, value) { + return new SparseMatrix(matrix).leftShiftS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + signPropagatingRightShift(value) { + if (typeof value === 'number') { + return this.signPropagatingRightShiftS(value); + } + return this.signPropagatingRightShiftM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + signPropagatingRightShiftS(value) { + this.forEachNonZero((i, j, v) => v >> value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + signPropagatingRightShiftM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) >> v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static signPropagatingRightShift(matrix, value) { + return new SparseMatrix(matrix).signPropagatingRightShiftS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + rightShift(value) { + if (typeof value === 'number') return this.rightShiftS(value); + return this.rightShiftM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + rightShiftS(value) { + this.forEachNonZero((i, j, v) => v >>> value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + rightShiftM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) >>> v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static rightShift(matrix, value) { + return new SparseMatrix(matrix).rightShiftS(value); + } + + /** + * @param {number|SparseMatrix} value + * @returns {this} + */ + zeroFillRightShift(value) { + if (typeof value === 'number') return this.zeroFillRightShiftS(value); + return this.zeroFillRightShiftM(value); + } + + /** + * @param {number} value + * @returns {this} + */ + zeroFillRightShiftS(value) { + this.forEachNonZero((i, j, v) => v >>> value); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {this} + */ + zeroFillRightShiftM(matrix) { + matrix.forEachNonZero((i, j, v) => { + this.set(i, j, this.get(i, j) >>> v); + return v; + }); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @param {number} value + * @returns {SparseMatrix} + */ + static zeroFillRightShift(matrix, value) { + return new SparseMatrix(matrix).zeroFillRightShiftS(value); + } + + /** + * @returns {this} + */ + not() { + this.forEachNonZero((i, j, v) => ~v); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static not(matrix) { + return new SparseMatrix(matrix).not(); + } + + /** + * @returns {this} + */ + abs() { + this.forEachNonZero((i, j, v) => Math.abs(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static abs(matrix) { + return new SparseMatrix(matrix).abs(); + } + + /** + * @returns {this} + */ + acos() { + this.forEachNonZero((i, j, v) => Math.acos(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static acos(matrix) { + return new SparseMatrix(matrix).acos(); + } + + /** + * @returns {this} + */ + acosh() { + this.forEachNonZero((i, j, v) => Math.acosh(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static acosh(matrix) { + return new SparseMatrix(matrix).acosh(); + } + + /** + * @returns {this} + */ + asin() { + this.forEachNonZero((i, j, v) => Math.asin(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static asin(matrix) { + return new SparseMatrix(matrix).asin(); + } + + /** + * @returns {this} + */ + asinh() { + this.forEachNonZero((i, j, v) => Math.asinh(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static asinh(matrix) { + return new SparseMatrix(matrix).asinh(); + } + + /** + * @returns {this} + */ + atan() { + this.forEachNonZero((i, j, v) => Math.atan(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static atan(matrix) { + return new SparseMatrix(matrix).atan(); + } + + /** + * @returns {this} + */ + atanh() { + this.forEachNonZero((i, j, v) => Math.atanh(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static atanh(matrix) { + return new SparseMatrix(matrix).atanh(); + } + + /** + * @returns {this} + */ + cbrt() { + this.forEachNonZero((i, j, v) => Math.cbrt(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static cbrt(matrix) { + return new SparseMatrix(matrix).cbrt(); + } + + /** + * @returns {this} + */ + ceil() { + this.forEachNonZero((i, j, v) => Math.ceil(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static ceil(matrix) { + return new SparseMatrix(matrix).ceil(); + } + + /** + * @returns {this} + */ + clz32() { + this.forEachNonZero((i, j, v) => Math.clz32(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static clz32(matrix) { + return new SparseMatrix(matrix).clz32(); + } + + /** + * @returns {this} + */ + cos() { + this.forEachNonZero((i, j, v) => Math.cos(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static cos(matrix) { + return new SparseMatrix(matrix).cos(); + } + + /** + * @returns {this} + */ + cosh() { + this.forEachNonZero((i, j, v) => Math.cosh(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static cosh(matrix) { + return new SparseMatrix(matrix).cosh(); + } + + /** + * @returns {this} + */ + exp() { + this.forEachNonZero((i, j, v) => Math.exp(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static exp(matrix) { + return new SparseMatrix(matrix).exp(); + } + + /** + * @returns {this} + */ + expm1() { + this.forEachNonZero((i, j, v) => Math.expm1(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static expm1(matrix) { + return new SparseMatrix(matrix).expm1(); + } + + /** + * @returns {this} + */ + floor() { + this.forEachNonZero((i, j, v) => Math.floor(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static floor(matrix) { + return new SparseMatrix(matrix).floor(); + } + + /** + * @returns {this} + */ + fround() { + this.forEachNonZero((i, j, v) => Math.fround(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static fround(matrix) { + return new SparseMatrix(matrix).fround(); + } + + /** + * @returns {this} + */ + log() { + this.forEachNonZero((i, j, v) => Math.log(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static log(matrix) { + return new SparseMatrix(matrix).log(); + } + + /** + * @returns {this} + */ + log1p() { + this.forEachNonZero((i, j, v) => Math.log1p(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static log1p(matrix) { + return new SparseMatrix(matrix).log1p(); + } + + /** + * @returns {this} + */ + log10() { + this.forEachNonZero((i, j, v) => Math.log10(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static log10(matrix) { + return new SparseMatrix(matrix).log10(); + } + + /** + * @returns {this} + */ + log2() { + this.forEachNonZero((i, j, v) => Math.log2(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static log2(matrix) { + return new SparseMatrix(matrix).log2(); + } + + /** + * @returns {this} + */ + round() { + this.forEachNonZero((i, j, v) => Math.round(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static round(matrix) { + return new SparseMatrix(matrix).round(); + } + + /** + * @returns {this} + */ + sign() { + this.forEachNonZero((i, j, v) => Math.sign(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static sign(matrix) { + return new SparseMatrix(matrix).sign(); + } + + /** + * @returns {this} + */ + sin() { + this.forEachNonZero((i, j, v) => Math.sin(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static sin(matrix) { + return new SparseMatrix(matrix).sin(); + } + + /** + * @returns {this} + */ + sinh() { + this.forEachNonZero((i, j, v) => Math.sinh(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static sinh(matrix) { + return new SparseMatrix(matrix).sinh(); + } + + /** + * @returns {this} + */ + sqrt() { + this.forEachNonZero((i, j, v) => Math.sqrt(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static sqrt(matrix) { + return new SparseMatrix(matrix).sqrt(); + } + + /** + * @returns {this} + */ + tan() { + this.forEachNonZero((i, j, v) => Math.tan(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static tan(matrix) { + return new SparseMatrix(matrix).tan(); + } + + /** + * @returns {this} + */ + tanh() { + this.forEachNonZero((i, j, v) => Math.tanh(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static tanh(matrix) { + return new SparseMatrix(matrix).tanh(); + } + + /** + * @returns {this} + */ + trunc() { + this.forEachNonZero((i, j, v) => Math.trunc(v)); + return this; + } + + /** + * @param {SparseMatrix} matrix + * @returns {SparseMatrix} + */ + static trunc(matrix) { + return new SparseMatrix(matrix).trunc(); + } + //#endregion +} + +SparseMatrix.prototype.klass = 'Matrix'; + +SparseMatrix.identity = SparseMatrix.eye; +SparseMatrix.prototype.tensorProduct = SparseMatrix.prototype.kroneckerProduct; From 12758f9f111800a5e177c21217c0917e8f058454 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Sat, 21 Jun 2025 00:17:08 -0500 Subject: [PATCH 11/33] feat: add tests for getNonZeros method in SparseMatrix --- src/__tests__/getNonZeros.test.js | 37 +++++++++++++++++++++++++++++++ 1 file changed, 37 insertions(+) create mode 100644 src/__tests__/getNonZeros.test.js diff --git a/src/__tests__/getNonZeros.test.js b/src/__tests__/getNonZeros.test.js new file mode 100644 index 0000000..1b93740 --- /dev/null +++ b/src/__tests__/getNonZeros.test.js @@ -0,0 +1,37 @@ +import { SparseMatrix } from '../index'; + +describe('Sparse Matrix', () => { + it('getNonZeros', () => { + let m2 = new SparseMatrix( + [ + [0, 0, 0, 0, 0, 0], + [1, 0, 0, 2, 1, 1], + [0, 3, 0, 0, 5, 5], + [0, 0, 0, 0, 0, 0], + [1, 0, 0, 1, 9, 9], + ], + { initialCapacity: 12 }, + ); + + // Default (coordinate list) + expect(m2.getNonZeros()).toEqual({ + rows: Float64Array.from([4, 1, 4, 1, 1, 1, 4, 2, 4, 2, 2]), + columns: Float64Array.from([0, 0, 3, 3, 4, 5, 5, 1, 4, 4, 5]), + values: Float64Array.from([1, 1, 1, 2, 1, 1, 9, 3, 9, 5, 5]), + }); + + // CSR format + expect(m2.getNonZeros({ format: 'csr' })).toEqual({ + rows: Float64Array.from([0, 0, 4, 7, 7, 11]), + columns: Float64Array.from([0, 3, 4, 5, 1, 4, 5, 0, 3, 4, 5]), + values: Float64Array.from([1, 2, 1, 1, 3, 5, 5, 1, 1, 9, 9]), + }); + + //CSC format + expect(m2.getNonZeros({ format: 'csc' })).toEqual({ + rows: Float64Array.from([1, 4, 2, 1, 4, 1, 2, 4, 1, 2, 4]), + columns: Float64Array.from([0, 2, 3, 3, 5, 8, 11]), + values: Float64Array.from([1, 1, 3, 2, 1, 1, 5, 9, 1, 5, 9]), + }); + }); +}); From d78911dd207943011ece4b0aaf2f053ed0a9d812 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Sat, 21 Jun 2025 01:53:21 -0500 Subject: [PATCH 12/33] feat: update benchmark parameters and add comprehensive tests for SparseMatrix operations --- benchmark/benchmark.js | 17 +-- src/__tests__/index.test.js | 52 +++++++++ src/__tests__/matrix.test.js | 198 +++++++++++++++++++++++++++++++++++ src/index.js | 190 ++++++++++++++++++--------------- 4 files changed, 367 insertions(+), 90 deletions(-) create mode 100644 src/__tests__/matrix.test.js diff --git a/benchmark/benchmark.js b/benchmark/benchmark.js index bb1eac5..2557758 100644 --- a/benchmark/benchmark.js +++ b/benchmark/benchmark.js @@ -54,10 +54,12 @@ function printWinner(label1, avg1, label2, avg2) { function runBenchmarks() { const m = 128; - const n = 64; + const n = 128; const p = 128; - const densityA = 0.03; - const densityB = 0.03; + const densityA = 0.01; + const densityB = 0.01; + + const nbIterations = 3; const A = randomSparseMatrix(m, n, densityA); const B = randomSparseMatrix(n, p, densityB); @@ -66,10 +68,11 @@ function runBenchmarks() { const AOld = new SparseMatrixOld(denseA); const BOld = new SparseMatrixOld(denseB); - A.cardinality; + denseA = new Matrix(denseA); denseB = new Matrix(denseB); + denseA.mmul(denseB); // 1. add vs addNew // const addAvg = benchmark(() => { // const a = AOld.clone(); @@ -90,7 +93,7 @@ function runBenchmarks() { A.mmul(B); }, 'mmulNew', - 3, + nbIterations, ); const mmulAvg = benchmark( @@ -98,7 +101,7 @@ function runBenchmarks() { AOld.mmul(BOld); }, 'mmul', - 3, + nbIterations, ); const denseAvg = benchmark( @@ -106,7 +109,7 @@ function runBenchmarks() { denseA.mmul(denseB); }, 'denseMatrix', - 3, + nbIterations, ); printWinner('mmul', mmulAvg, 'mmulNew', mmulNewAvg); diff --git a/src/__tests__/index.test.js b/src/__tests__/index.test.js index d0c1caf..1a6b38a 100644 --- a/src/__tests__/index.test.js +++ b/src/__tests__/index.test.js @@ -42,6 +42,29 @@ describe('Sparse Matrix', () => { [0, 2], [0, 0], ]); + + // Compare with dense multiplication + const denseM1 = m1.to2DArray(); + const denseM2 = m2.to2DArray(); + const expectedDense = denseMatrixMultiply(denseM1, denseM2); + expect(m3.to2DArray()).toStrictEqual(expectedDense); + }); + + it('mmul', () => { + const size = 32; + const density = 0.1; + const m1 = randomSparseMatrix(size, size, density); + const m2 = randomSparseMatrix(size, size, density); + let m3 = m1.mmul(m2); + + const denseM1 = m1.to2DArray(); + const denseM2 = m2.to2DArray(); + + const newSparse = new SparseMatrix(denseM1); + expect(newSparse.to2DArray()).toStrictEqual(denseM1); + const expectedDense = denseMatrixMultiply(denseM1, denseM2); + + expect(m3.to2DArray()).toStrictEqual(expectedDense); }); it('kronecker', () => { @@ -144,3 +167,32 @@ describe('Banded matrices', () => { expect(matrix4.isBanded(1)).toBe(true); }); }); + +function denseMatrixMultiply(A, B) { + const rowsA = A.length; + const colsA = A[0].length; + const colsB = B[0].length; + const result = Array.from({ length: rowsA }, () => Array(colsB).fill(0)); + for (let i = 0; i < rowsA; i++) { + for (let j = 0; j < colsB; j++) { + for (let k = 0; k < colsA; k++) { + result[i][j] += A[i][k] * B[k][j]; + } + } + } + return result; +} + +function randomSparseMatrix(rows, cols, density = 0.01) { + const matrix = []; + for (let i = 0; i < rows; i++) { + const row = new Float64Array(cols); + for (let j = 0; j < cols; j++) { + if (Math.random() < density) { + row[j] = Math.random() * 10; + } + } + matrix.push(row); + } + return new SparseMatrix(matrix); +} diff --git a/src/__tests__/matrix.test.js b/src/__tests__/matrix.test.js new file mode 100644 index 0000000..1a6b38a --- /dev/null +++ b/src/__tests__/matrix.test.js @@ -0,0 +1,198 @@ +import { describe, expect, it } from 'vitest'; + +import { SparseMatrix } from '../index.js'; + +describe('Sparse Matrix', () => { + it('add', () => { + let m1 = new SparseMatrix([ + [2, 0, 1], + [0, 0, 3], + [2, 0, 1], + ]); + let m2 = new SparseMatrix([ + [0, 1, 5], + [2, 0, 0], + [-2, 0, -1], + ]); + let m3 = m1.add(m2).to2DArray(); + expect(m3).toStrictEqual([ + [2, 1, 6], + [2, 0, 3], + [0, 0, 0], + ]); + }); + it('mmul', () => { + let m1 = new SparseMatrix([ + [2, 0, 1], + [0, 0, 3], + ]); + let m2 = new SparseMatrix([ + [0, 1], + [2, 0], + [0, 0], + ]); + let m3 = m1.mmul(m2); + + expect(m1.cardinality).toBe(3); + expect(m2.cardinality).toBe(2); + expect(m3.cardinality).toBe(1); + + expect(m3.get(0, 1)).toBe(2); + expect(m3.to2DArray()).toStrictEqual([ + [0, 2], + [0, 0], + ]); + + // Compare with dense multiplication + const denseM1 = m1.to2DArray(); + const denseM2 = m2.to2DArray(); + const expectedDense = denseMatrixMultiply(denseM1, denseM2); + expect(m3.to2DArray()).toStrictEqual(expectedDense); + }); + + it('mmul', () => { + const size = 32; + const density = 0.1; + const m1 = randomSparseMatrix(size, size, density); + const m2 = randomSparseMatrix(size, size, density); + let m3 = m1.mmul(m2); + + const denseM1 = m1.to2DArray(); + const denseM2 = m2.to2DArray(); + + const newSparse = new SparseMatrix(denseM1); + expect(newSparse.to2DArray()).toStrictEqual(denseM1); + const expectedDense = denseMatrixMultiply(denseM1, denseM2); + + expect(m3.to2DArray()).toStrictEqual(expectedDense); + }); + + it('kronecker', () => { + const matrix1 = new SparseMatrix([ + [1, 2], + [3, 4], + ]); + const matrix2 = new SparseMatrix([ + [0, 5], + [6, 7], + ]); + const product = matrix1.kroneckerProduct(matrix2); + expect(product.to2DArray()).toStrictEqual([ + [0, 5, 0, 10], + [6, 7, 12, 14], + [0, 15, 0, 20], + [18, 21, 24, 28], + ]); + }); + + it('isSymmetric', () => { + expect(new SparseMatrix(10, 10).isSymmetric()).toBe(true); + expect(new SparseMatrix(15, 10).isSymmetric()).toBe(false); + + let m = new SparseMatrix([ + [0, 1], + [1, 0], + ]); + expect(m.isSymmetric()).toBe(true); + + m = new SparseMatrix([ + [0, 1], + [0, 1], + ]); + expect(m.isSymmetric()).toBe(false); + }); + + it('transpose', () => { + const matrix = new SparseMatrix([ + [1, 2], + [3, 4], + ]); + expect(matrix.transpose().to2DArray()).toStrictEqual([ + [1, 3], + [2, 4], + ]); + }); +}); + +describe('Banded matrices', () => { + it('Check band size', () => { + const matrix1 = new SparseMatrix([ + [1, 0], + [0, 1], + ]); + const matrix2 = new SparseMatrix([ + [1, 0, 0], + [0, 1, 0], + ]); + const matrix3 = new SparseMatrix([ + [1, 0, 1], + [0, 1, 0], + ]); + const matrix4 = new SparseMatrix([ + [1, 0, 0], + [1, 1, 0], + ]); + const matrix5 = new SparseMatrix([ + [0, 0, 0], + [1, 0, 0], + [0, 1, 0], + ]); + expect(matrix1.bandWidth()).toBe(0); + expect(matrix2.bandWidth()).toBe(0); + expect(matrix3.bandWidth()).toBe(2); + expect(matrix4.bandWidth()).toBe(1); + expect(matrix5.bandWidth()).toBe(0); + }); + + it('isBanded', () => { + const matrix1 = new SparseMatrix([ + [1, 0], + [0, 1], + ]); + const matrix2 = new SparseMatrix([ + [1, 0, 0], + [0, 1, 0], + ]); + const matrix3 = new SparseMatrix([ + [1, 0, 1], + [0, 1, 0], + ]); + const matrix4 = new SparseMatrix([ + [1, 0, 0], + [1, 1, 0], + ]); + expect(matrix1.isBanded(1)).toBe(true); + expect(matrix2.isBanded(1)).toBe(true); + expect(matrix3.isBanded(1)).toBe(false); + expect(matrix4.isBanded(1)).toBe(true); + }); +}); + +function denseMatrixMultiply(A, B) { + const rowsA = A.length; + const colsA = A[0].length; + const colsB = B[0].length; + const result = Array.from({ length: rowsA }, () => Array(colsB).fill(0)); + for (let i = 0; i < rowsA; i++) { + for (let j = 0; j < colsB; j++) { + for (let k = 0; k < colsA; k++) { + result[i][j] += A[i][k] * B[k][j]; + } + } + } + return result; +} + +function randomSparseMatrix(rows, cols, density = 0.01) { + const matrix = []; + for (let i = 0; i < rows; i++) { + const row = new Float64Array(cols); + for (let j = 0; j < cols; j++) { + if (Math.random() < density) { + row[j] = Math.random() * 10; + } + } + matrix.push(row); + } + return new SparseMatrix(matrix); +} diff --git a/src/index.js b/src/index.js index f4fe6ce..0880ac0 100644 --- a/src/index.js +++ b/src/index.js @@ -156,62 +156,35 @@ export class SparseMatrix { return this; } - /** - * @param {SparseMatrix} other - * @returns {SparseMatrix} - */ - mulNew(other) { - if (typeof other !== 'number') { - throw new RangeError('the argument should be a number'); + add(other) { + if (typeof other === 'number') return this.addS(other); + if (this.rows !== other.rows || this.columns !== other.columns) { + throw new RangeError('Matrices dimensions must be equal'); } - // if (this.cardinality / this.columns / this.rows > 0.1) - this.withEachNonZero((i, j, v) => { - // return v * other; - this.set(i, j, v * other); + other.withEachNonZero((i, j, v) => { + this.set(i, j, v + this.get(i, j)); }); return this; } - mmulNew(other) { - if (this.columns !== other.rows) { - // eslint-disable-next-line no-console - console.warn( - 'Number of columns of left matrix are not equal to number of rows of right matrix.', - ); + mul(other) { + if (typeof other !== 'number') { + throw new RangeError('the argument should be a number'); } - const m = this.rows; - const p = other.columns; - - const { - columns: otherCols, - rows: otherRows, - values: otherValues, - } = other.getNonZeros(); - const { - columns: thisCols, - rows: thisRows, - values: thisValues, - } = this.getNonZeros(); - - const result = new SparseMatrix(m, p, { initialCapacity: m * p }); - - const nbOtherActive = otherCols.length; - const nbThisActive = thisCols.length; - for (let t = 0; t < nbThisActive; t++) { - const i = thisRows[t]; - const j = thisCols[t]; - for (let o = 0; o < nbOtherActive; o++) { - if (j === otherRows[o]) { - const l = otherCols[o]; - result.set(i, l, result.get(i, l) + otherValues[o] * thisValues[t]); - } - } + if (other === 0) { + return new SparseMatrix(this.rows, this.columns, { + initialCapacity: this.cardinality, + }); } - return result; + this.withEachNonZero((i, j, v) => { + this.set(i, j, v * other); + }); + + return this; } mmul(other) { @@ -225,40 +198,51 @@ export class SparseMatrix { const m = this.rows; const p = other.columns; + const result = matrixCreateEmpty(m, p); + const useCscCsr = useCscCsrFormat(this, other); + const { columns: otherCols, rows: otherRows, values: otherValues, - } = other.getNonZeros('csr'); + } = other.getNonZeros(useCscCsr ? { format: 'csr' } : {}); const { columns: thisCols, rows: thisRows, values: thisValues, - } = this.getNonZeros('csc'); - - const result = new SparseMatrix(m, p, { initialCapacity: m * p + 20 }); - - for (let t = 0; t < thisCols.length - 1; t++) { - const j = t; - const tValues = thisValues.subarray(thisCols[t], thisCols[t + 1]); - const tRows = thisRows.subarray(thisCols[t], thisCols[t + 1]); - let initOther = 0; - for (let o = initOther; o < otherRows.length - 1; o++) { - if (o === j) { - initOther++; - const oValues = otherValues.subarray(otherRows[o], otherRows[o + 1]); - const oCols = otherCols.subarray(otherRows[o], otherRows[o + 1]); - for (let f = 0; f < tValues.length; f++) { - for (let k = 0; k < oValues.length; k++) { - const i = tRows[f]; - const l = oCols[k]; - result.set(i, l, result.get(i, l) + tValues[f] * oValues[k]); - } + } = this.getNonZeros(useCscCsr ? { format: 'csc' } : {}); + + if (useCscCsr) { + const thisNbCols = this.columns; + for (let t = 0; t < thisNbCols; t++) { + const tValues = thisValues.subarray(thisCols[t], thisCols[t + 1]); + const tRows = thisRows.subarray(thisCols[t], thisCols[t + 1]); + const oValues = otherValues.subarray(otherRows[t], otherRows[t + 1]); + const oCols = otherCols.subarray(otherRows[t], otherRows[t + 1]); + for (let f = 0; f < tValues.length; f++) { + for (let k = 0; k < oValues.length; k++) { + const i = tRows[f]; + const l = oCols[k]; + result[i][l] += tValues[f] * oValues[k]; + } + } + } + } else { + const nbOtherActive = otherCols.length; + const nbThisActive = thisCols.length; + for (let t = 0; t < nbThisActive; t++) { + const i = thisRows[t]; + const j = thisCols[t]; + for (let o = 0; o < nbOtherActive; o++) { + if (j === otherRows[o]) { + const l = otherCols[o]; + result[i][l] += otherValues[o] * thisValues[t]; } - } else if (j < o) break; + } } } - return result; + + return new SparseMatrix(result); } /** @@ -286,8 +270,8 @@ export class SparseMatrix { values: thisValues, } = this.getNonZeros(); - const nbOtherActive = otherCols.length; const nbThisActive = thisCols.length; + const nbOtherActive = otherCols.length; for (let t = 0; t < nbThisActive; t++) { const pi = p * thisRows[t]; const qj = q * thisCols[t]; @@ -303,13 +287,18 @@ export class SparseMatrix { return result; } - withEachNonZero(callback) { + withEachNonZero(callback, sort = false) { const { state, table, values } = this.elements; const nbStates = state.length; - const activeIndex = []; - for (let i = 0; i < nbStates; i++) { - if (state[i] === 1) activeIndex.push(i); + const activeIndex = new Float64Array(this.cardinality); + for (let i = 0, j = 0; i < nbStates; i++) { + if (state[i] === 1) activeIndex[j++] = i; } + + if (sort) { + activeIndex.sort((a, b) => table[a] - table[b]); + } + const columns = this.columns; for (const i of activeIndex) { const key = table[i]; @@ -347,8 +336,7 @@ export class SparseMatrix { return this; } - //'csr' | 'csc' | undefined - getNonZeros(format) { + getNonZeros(options = {}) { const cardinality = this.cardinality; /** @type {Float64Array} */ const rows = new Float64Array(cardinality); @@ -357,20 +345,28 @@ export class SparseMatrix { /** @type {Float64Array} */ const values = new Float64Array(cardinality); + const { format, sort = false } = options; + let idx = 0; this.withEachNonZero((i, j, value) => { rows[idx] = i; columns[idx] = j; values[idx] = value; idx++; - }); + }, sort || format); const cooMatrix = { rows, columns, values }; if (!format) return cooMatrix; + if (!['csr', 'csc'].includes(format.toLowerCase())) { + throw new Error(`format ${format} is not supported`); + } + const csrMatrix = cooToCsr(cooMatrix, this.rows); - return format === 'csc' ? csrToCsc(csrMatrix, this.columns) : csrMatrix; + return format.toLowerCase() === 'csc' + ? csrToCsc(csrMatrix, this.columns) + : csrMatrix; } /** @@ -1510,25 +1506,21 @@ function csrToCsc(csrMatrix, numCols) { columns: csrColIndices, rows: csrRowPtr, } = csrMatrix; - // Initialize CSC arrays + const cscValues = new Float64Array(csrValues.length); const cscRowIndices = new Float64Array(csrValues.length); const cscColPtr = new Float64Array(numCols + 1); - // Count non-zeros per column for (let i = 0; i < csrColIndices.length; i++) { cscColPtr[csrColIndices[i] + 1]++; } - // Compute column pointers (prefix sum) for (let i = 1; i <= numCols; i++) { cscColPtr[i] += cscColPtr[i - 1]; } - // Temporary copy for filling values const next = cscColPtr.slice(); - // Fill CSC arrays for (let row = 0; row < csrRowPtr.length - 1; row++) { for (let j = csrRowPtr[row], i = 0; j < csrRowPtr[row + 1]; j++, i++) { const col = csrColIndices[j]; @@ -1537,7 +1529,6 @@ function csrToCsc(csrMatrix, numCols) { cscRowIndices[pos] = row; next[col]++; } - // if (row === 1) break; } return { rows: cscRowIndices, columns: cscColPtr, values: cscValues }; @@ -1557,3 +1548,36 @@ function cooToCsr(cooMatrix, nbRows = 9) { } return { rows: csrRowPtr, columns, values }; } + +function matrixCreateEmpty(nbRows, nbColumns) { + const newMatrix = []; + for (let row = 0; row < nbRows; row++) { + newMatrix.push(new Float64Array(nbColumns)); + } + return newMatrix; +} + +function useCscCsrFormat(A, B) { + // Validar dimensiones + if (!A.rows || !A.columns || !B.rows || !B.columns || A.columns !== B.rows) { + throw new Error('Dimensiones de matrices inválidas o no compatibles.'); + } + + const densityA = A.cardinality / (A.rows * A.columns); + const densityB = B.cardinality / (B.rows * B.columns); + + const M = A.rows; + const K = A.columns; + const N = B.columns; + + const isHighDensity = densityA >= 0.015 || densityB >= 0.015; + const isVeryHighDensity = densityB >= 0.03; + const isLargeMatrix = M >= 64 || K >= 64 || N >= 64; + const isLargeSquareMatrix = M === K && K === N && M >= 128; + + return ( + isVeryHighDensity || + (isHighDensity && isLargeMatrix) || + (isLargeSquareMatrix && densityA >= 0.03) + ); +} From 807ab3dac1583e50c9cb7d1c8c54963d4d0d089e Mon Sep 17 00:00:00 2001 From: jobo322 Date: Mon, 23 Jun 2025 06:41:01 -0500 Subject: [PATCH 13/33] chore: keep mmul simple --- src/index.js | 50 ++++++++++++++++---------------------------------- 1 file changed, 16 insertions(+), 34 deletions(-) diff --git a/src/index.js b/src/index.js index 0880ac0..bc15288 100644 --- a/src/index.js +++ b/src/index.js @@ -199,45 +199,29 @@ export class SparseMatrix { const p = other.columns; const result = matrixCreateEmpty(m, p); - const useCscCsr = useCscCsrFormat(this, other); const { columns: otherCols, rows: otherRows, values: otherValues, - } = other.getNonZeros(useCscCsr ? { format: 'csr' } : {}); + } = other.getNonZeros({ format: 'csr' }); const { columns: thisCols, rows: thisRows, values: thisValues, - } = this.getNonZeros(useCscCsr ? { format: 'csc' } : {}); - - if (useCscCsr) { - const thisNbCols = this.columns; - for (let t = 0; t < thisNbCols; t++) { - const tValues = thisValues.subarray(thisCols[t], thisCols[t + 1]); - const tRows = thisRows.subarray(thisCols[t], thisCols[t + 1]); - const oValues = otherValues.subarray(otherRows[t], otherRows[t + 1]); - const oCols = otherCols.subarray(otherRows[t], otherRows[t + 1]); - for (let f = 0; f < tValues.length; f++) { - for (let k = 0; k < oValues.length; k++) { - const i = tRows[f]; - const l = oCols[k]; - result[i][l] += tValues[f] * oValues[k]; - } - } - } - } else { - const nbOtherActive = otherCols.length; - const nbThisActive = thisCols.length; - for (let t = 0; t < nbThisActive; t++) { - const i = thisRows[t]; - const j = thisCols[t]; - for (let o = 0; o < nbOtherActive; o++) { - if (j === otherRows[o]) { - const l = otherCols[o]; - result[i][l] += otherValues[o] * thisValues[t]; - } + } = this.getNonZeros({ format: 'csc' }); + + const thisNbCols = this.columns; + for (let t = 0; t < thisNbCols; t++) { + const tStart = thisCols[t]; + const tEnd = thisCols[t + 1]; + const oStart = otherRows[t]; + const oEnd = otherRows[t + 1]; + for (let f = tStart; f < tEnd; f++) { + for (let k = oStart; k < oEnd; k++) { + const i = thisRows[f]; + const l = otherCols[k]; + result[i][l] += thisValues[f] * otherValues[k]; } } } @@ -355,15 +339,13 @@ export class SparseMatrix { idx++; }, sort || format); - const cooMatrix = { rows, columns, values }; - - if (!format) return cooMatrix; + if (!format) return { rows, columns, values }; if (!['csr', 'csc'].includes(format.toLowerCase())) { throw new Error(`format ${format} is not supported`); } - const csrMatrix = cooToCsr(cooMatrix, this.rows); + const csrMatrix = cooToCsr({ rows, columns, values }, this.rows); return format.toLowerCase() === 'csc' ? csrToCsc(csrMatrix, this.columns) : csrMatrix; From 09d0a9832027ac29776c72f8c6eee9c807d894ef Mon Sep 17 00:00:00 2001 From: jobo322 Date: Mon, 23 Jun 2025 09:40:52 -0500 Subject: [PATCH 14/33] chore: use benchmark package --- benchmark/benchmark.js | 67 ++++++++++++++--------- benchmark/squareMatrixBenchmark.js | 64 ++++++++++++++++++++++ benchmark/utils/benchmark.js | 16 ++++++ benchmark/utils/printWinner.js | 27 +++++++++ benchmark/utils/randomSparseMatrix.js | 24 ++++++++ benchmark/variateSizeDensityBenchmark.js | 70 ++++++++++++++++++++++++ src/index.js | 25 --------- 7 files changed, 243 insertions(+), 50 deletions(-) create mode 100644 benchmark/squareMatrixBenchmark.js create mode 100644 benchmark/utils/benchmark.js create mode 100644 benchmark/utils/printWinner.js create mode 100644 benchmark/utils/randomSparseMatrix.js create mode 100644 benchmark/variateSizeDensityBenchmark.js diff --git a/benchmark/benchmark.js b/benchmark/benchmark.js index 2557758..8d6d01b 100644 --- a/benchmark/benchmark.js +++ b/benchmark/benchmark.js @@ -1,23 +1,33 @@ import { SparseMatrix } from '../src/index.js'; import { Matrix } from 'ml-matrix'; import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; +// import { SparseMatrix as SparseMatrixOld } from '../src/Elements.js'; import fs from 'fs'; function randomSparseMatrix(rows, cols, density = 0.01) { - const matrix = []; - for (let i = 0; i < rows; i++) { - const row = new Array(cols).fill(0); - for (let j = 0; j < cols; j++) { - if (Math.random() < density) { - row[j] = Math.random() * 10; - } - } - matrix.push(row); + const total = rows * cols; + const cardinality = Math.round(total * density); + const positions = new Set(); + + // Generate unique random positions + while (positions.size < cardinality) { + positions.add(Math.floor(Math.random() * total)); } + + // Build the matrix with zeros + const matrix = Array.from({ length: rows }, () => new Float64Array(cols)); + + // Assign random values to the selected positions + for (const pos of positions) { + const i = Math.floor(pos / cols); + const j = pos % cols; + matrix[i][j] = Math.random() * 10; + } + return new SparseMatrix(matrix); } -function benchmark(fn, label, iterations = 5) { +function benchmark(fn, label, iterations = 5, logIt = false) { const times = []; for (let i = 0; i < iterations; i++) { const t0 = performance.now(); @@ -26,7 +36,9 @@ function benchmark(fn, label, iterations = 5) { times.push(t1 - t0); } const avg = times.reduce((a, b) => a + b, 0) / times.length; - console.log(`${label}: avg ${avg.toFixed(2)} ms over ${iterations} runs`); + if (logIt) { + console.log(`${label}: avg ${avg.toFixed(2)} ms over ${iterations} runs`); + } return avg; } @@ -53,9 +65,9 @@ function printWinner(label1, avg1, label2, avg2) { } function runBenchmarks() { - const m = 128; - const n = 128; - const p = 128; + const m = 256; + const n = 256; + const p = 256; const densityA = 0.01; const densityB = 0.01; @@ -137,8 +149,9 @@ function runBenchmarks() { } function runSizeSweepBenchmark() { - const sizes = [8, 16, 32, 64, 128]; - const densities = [0.01, 0.015, 0.02, 0.025, 0.03, 0.35]; + const nbIterations = 3; + const sizes = [32, 64, 128, 256]; + const densities = [0.01, 0.015, 0.02, 0.025, 0.03]; const results = []; for (const densityA of densities) { @@ -165,7 +178,7 @@ function runSizeSweepBenchmark() { A.mmul(B); }, 'mmulNew', - 3, + nbIterations, ); const mmulAvg = benchmark( @@ -173,12 +186,16 @@ function runSizeSweepBenchmark() { AOld.mmul(BOld); }, 'mmul', - 3, + nbIterations, ); - const denseAvg = benchmark(() => { - denseA.mmul(denseB), 'denseMatrix', 3; - }); + const denseAvg = benchmark( + () => { + denseA.mmul(denseB); + }, + 'denseMatrix', + nbIterations, + ); results.push({ densityA, @@ -186,8 +203,8 @@ function runSizeSweepBenchmark() { A_shape: [m, n], B_shape: [n, p], dense: denseAvg, - mmulNew: mmulNewAvg, - mmul: mmulAvg, + new: mmulNewAvg, + old: mmulAvg, }); } } @@ -202,6 +219,6 @@ function runSizeSweepBenchmark() { console.log('Size sweep benchmark results saved to size_sweep_results.json'); } -runBenchmarks(); +// runBenchmarks(); // Uncomment to run the size sweep benchmark -// runSizeSweepBenchmark(); +runSizeSweepBenchmark(); diff --git a/benchmark/squareMatrixBenchmark.js b/benchmark/squareMatrixBenchmark.js new file mode 100644 index 0000000..1e4e571 --- /dev/null +++ b/benchmark/squareMatrixBenchmark.js @@ -0,0 +1,64 @@ +import Benchmark from 'benchmark'; +import { Matrix } from 'ml-matrix'; +import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; +import fs from 'fs'; +import { randomSparseMatrix } from './utils/randomSparseMatrix.js'; + +function runSizeSweepBenchmarkCompact() { + const sizes = [32, 64, 128, 256]; + const densities = [0.01, 0.015, 0.02, 0.025, 0.03]; + const results = []; + + for (const densityA of densities) { + for (const m of sizes) { + const A = randomSparseMatrix(m, m, densityA); + const B = randomSparseMatrix(m, m, densityA); + let denseA = A.to2DArray(); + let denseB = B.to2DArray(); + const AOld = new SparseMatrixOld(denseA); + const BOld = new SparseMatrixOld(denseB); + denseA = new Matrix(denseA); + denseB = new Matrix(denseB); + // Warm up + A.mmul(B); + // Use Benchmark.js for each method + const suite = new Benchmark.Suite(); + let mmulNewAvg, mmulAvg, denseAvg; + suite + .add('mmulNew', function () { + A.mmul(B); + }) + .add('mmul', function () { + AOld.mmul(BOld); + }) + .add('denseMatrix', function () { + denseA.mmul(denseB); + }) + .on('cycle', function (event) { + // Optionally log: console.log(String(event.target)); + }) + .on('complete', function () { + this.forEach((bench) => { + if (bench.name === 'mmulNew') mmulNewAvg = 1000 / bench.hz; + if (bench.name === 'mmul') mmulAvg = 1000 / bench.hz; + if (bench.name === 'denseMatrix') denseAvg = 1000 / bench.hz; + }); + }) + .run({ async: false }); + results.push({ + density: densityA, + size: m, + dense: denseAvg, + new: mmulNewAvg, + old: mmulAvg, + }); + } + } + fs.writeFileSync( + './benchmark/size_sweep_results-dense.json', + JSON.stringify(results, null, 2), + ); + console.log('Size sweep benchmark results saved to size_sweep_results.json'); +} + +runSizeSweepBenchmarkCompact(); diff --git a/benchmark/utils/benchmark.js b/benchmark/utils/benchmark.js new file mode 100644 index 0000000..cc60e1c --- /dev/null +++ b/benchmark/utils/benchmark.js @@ -0,0 +1,16 @@ +export function benchmark(fn, label, iterations = 5, logIt = true) { + const times = []; + for (let i = 0; i < iterations; i++) { + const t0 = performance.now(); + fn(); + const t1 = performance.now(); + times.push(t1 - t0); + } + const avg = times.reduce((a, b) => a + b, 0) / times.length; + if (logIt) { + console.log( + `${label}: avg ${avg.toFixed(2)} ms over ${iterations} runs \n`, + ); + } + return avg; +} diff --git a/benchmark/utils/printWinner.js b/benchmark/utils/printWinner.js new file mode 100644 index 0000000..527e686 --- /dev/null +++ b/benchmark/utils/printWinner.js @@ -0,0 +1,27 @@ +export function printWinner(results) { + // results: Array of { label: string, avg: number } + if (!Array.isArray(results) || results.length < 2) { + console.log('Need at least two results to compare.'); + return; + } + + // Sort by avg ascending (fastest first) + const sorted = [...results].sort((a, b) => a.avg - b.avg); + const winner = sorted[0]; + + console.log( + ` -> ${winner.label} was the fastest (${winner.avg.toFixed(2)} ms avg)\n`, + ); + + for (let i = 1; i < sorted.length; i++) { + const slower = sorted[i]; + const speedup = (slower.avg / winner.avg).toFixed(2); + const percent = (((slower.avg - winner.avg) / slower.avg) * 100).toFixed(1); + console.log( + ` ${winner.label} was ${speedup}x faster (${percent}% faster) than ${ + slower.label + } (${slower.avg.toFixed(2)} ms avg)`, + ); + } + console.log(); +} diff --git a/benchmark/utils/randomSparseMatrix.js b/benchmark/utils/randomSparseMatrix.js new file mode 100644 index 0000000..61fe737 --- /dev/null +++ b/benchmark/utils/randomSparseMatrix.js @@ -0,0 +1,24 @@ +import { SparseMatrix } from '../../src/index.js'; + +export function randomSparseMatrix(rows, cols, density = 0.01) { + const total = rows * cols; + const cardinality = Math.round(total * density); + const positions = new Set(); + + // Generate unique random positions + while (positions.size < cardinality) { + positions.add(Math.floor(Math.random() * total)); + } + + // Build the matrix with zeros + const matrix = Array.from({ length: rows }, () => new Float64Array(cols)); + + // Assign random values to the selected positions + for (const pos of positions) { + const i = Math.floor(pos / cols); + const j = pos % cols; + matrix[i][j] = Math.random() * 10; + } + + return new SparseMatrix(matrix); +} diff --git a/benchmark/variateSizeDensityBenchmark.js b/benchmark/variateSizeDensityBenchmark.js new file mode 100644 index 0000000..39fe24c --- /dev/null +++ b/benchmark/variateSizeDensityBenchmark.js @@ -0,0 +1,70 @@ +import Benchmark from 'benchmark'; +import { Matrix } from 'ml-matrix'; +import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; +import fs from 'fs'; +import { randomSparseMatrix } from './utils/randomSparseMatrix.js'; + +function runSizeSweepBenchmark() { + const sizes = [32, 64, 128, 256]; + const densities = [0.01, 0.015, 0.02, 0.025, 0.03]; + const results = []; + + for (const densityA of densities) { + for (const densityB of densities) { + for (const m of sizes) { + for (const n of sizes) { + for (const p of sizes) { + const A = randomSparseMatrix(m, n, densityA); + const B = randomSparseMatrix(n, p, densityB); + let denseA = A.to2DArray(); + let denseB = B.to2DArray(); + const AOld = new SparseMatrixOld(denseA); + const BOld = new SparseMatrixOld(denseB); + denseA = new Matrix(denseA); + denseB = new Matrix(denseB); + // Use Benchmark.js for each method + const suite = new Benchmark.Suite(); + let mmulNewAvg, mmulAvg, denseAvg; + suite + .add('mmulNew', function () { + A.mmul(B); + }) + .add('mmul', function () { + AOld.mmul(BOld); + }) + .add('denseMatrix', function () { + denseA.mmul(denseB); + }) + .on('cycle', function (event) { + // Optionally log: console.log(String(event.target)); + }) + .on('complete', function () { + this.forEach((bench) => { + if (bench.name === 'mmulNew') mmulNewAvg = 1000 / bench.hz; + if (bench.name === 'mmul') mmulAvg = 1000 / bench.hz; + if (bench.name === 'denseMatrix') denseAvg = 1000 / bench.hz; + }); + }) + .run({ async: false }); + results.push({ + densityA, + densityB, + A_shape: [m, n], + B_shape: [n, p], + dense: denseAvg, + new: mmulNewAvg, + old: mmulAvg, + }); + } + } + } + } + } + fs.writeFileSync( + './benchmark/size_sweep_results-dense.json', + JSON.stringify(results, null, 2), + ); + console.log('Size sweep benchmark results saved to size_sweep_results.json'); +} + +runSizeSweepBenchmark(); diff --git a/src/index.js b/src/index.js index bc15288..23e6c4c 100644 --- a/src/index.js +++ b/src/index.js @@ -1538,28 +1538,3 @@ function matrixCreateEmpty(nbRows, nbColumns) { } return newMatrix; } - -function useCscCsrFormat(A, B) { - // Validar dimensiones - if (!A.rows || !A.columns || !B.rows || !B.columns || A.columns !== B.rows) { - throw new Error('Dimensiones de matrices inválidas o no compatibles.'); - } - - const densityA = A.cardinality / (A.rows * A.columns); - const densityB = B.cardinality / (B.rows * B.columns); - - const M = A.rows; - const K = A.columns; - const N = B.columns; - - const isHighDensity = densityA >= 0.015 || densityB >= 0.015; - const isVeryHighDensity = densityB >= 0.03; - const isLargeMatrix = M >= 64 || K >= 64 || N >= 64; - const isLargeSquareMatrix = M === K && K === N && M >= 128; - - return ( - isVeryHighDensity || - (isHighDensity && isLargeMatrix) || - (isLargeSquareMatrix && densityA >= 0.03) - ); -} From 72db1dd5b184d0c12ee346985c25bb9251377a64 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Mon, 23 Jun 2025 10:44:40 -0500 Subject: [PATCH 15/33] chore: simple benchmark run --- benchmark/benchmark.js | 281 +++++++++++++---------------------------- 1 file changed, 91 insertions(+), 190 deletions(-) diff --git a/benchmark/benchmark.js b/benchmark/benchmark.js index 8d6d01b..9836b52 100644 --- a/benchmark/benchmark.js +++ b/benchmark/benchmark.js @@ -1,68 +1,8 @@ -import { SparseMatrix } from '../src/index.js'; +import Benchmark from 'benchmark'; import { Matrix } from 'ml-matrix'; import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; -// import { SparseMatrix as SparseMatrixOld } from '../src/Elements.js'; -import fs from 'fs'; - -function randomSparseMatrix(rows, cols, density = 0.01) { - const total = rows * cols; - const cardinality = Math.round(total * density); - const positions = new Set(); - - // Generate unique random positions - while (positions.size < cardinality) { - positions.add(Math.floor(Math.random() * total)); - } - - // Build the matrix with zeros - const matrix = Array.from({ length: rows }, () => new Float64Array(cols)); - - // Assign random values to the selected positions - for (const pos of positions) { - const i = Math.floor(pos / cols); - const j = pos % cols; - matrix[i][j] = Math.random() * 10; - } - - return new SparseMatrix(matrix); -} - -function benchmark(fn, label, iterations = 5, logIt = false) { - const times = []; - for (let i = 0; i < iterations; i++) { - const t0 = performance.now(); - fn(); - const t1 = performance.now(); - times.push(t1 - t0); - } - const avg = times.reduce((a, b) => a + b, 0) / times.length; - if (logIt) { - console.log(`${label}: avg ${avg.toFixed(2)} ms over ${iterations} runs`); - } - return avg; -} - -function printWinner(label1, avg1, label2, avg2) { - let winner, loser, win, lose; - if (avg1 < avg2) { - winner = label1; - win = avg1; - loser = label2; - lose = avg2; - } else { - winner = label2; - win = avg2; - loser = label1; - lose = avg1; - } - - const percent = ((lose - win) / lose) * 100; - console.log( - ` -> ${winner} was ${(lose / win).toFixed(2)}x faster (${percent.toFixed( - 1, - )}% faster) than ${loser}\n`, - ); -} +import { randomSparseMatrix } from './utils/randomSparseMatrix.js'; +import { printWinner } from './utils/printWinner.js'; function runBenchmarks() { const m = 256; @@ -71,7 +11,6 @@ function runBenchmarks() { const densityA = 0.01; const densityB = 0.01; - const nbIterations = 3; const A = randomSparseMatrix(m, n, densityA); const B = randomSparseMatrix(n, p, densityB); @@ -84,141 +23,103 @@ function runBenchmarks() { denseA = new Matrix(denseA); denseB = new Matrix(denseB); - denseA.mmul(denseB); - // 1. add vs addNew - // const addAvg = benchmark(() => { - // const a = AOld.clone(); - // a.add(BOld); - // }, 'add'); - - // const addNewAvg = benchmark(() => { - // const a = A.clone(); - // a.add(B); - // }, 'addNew'); + const results = []; - // printWinner('add', addAvg, 'addNew', addNewAvg); + // 1. add vs addNew + const addSuite = new Benchmark.Suite('add'); + addSuite + .add('add', function () { + const a = AOld.clone(); + a.add(BOld); + }) + .add('addNew', function () { + const a = A.clone(); + a.add(B); + }) + .add('addDense', function () { + const a = denseA.clone(); + a.add(denseB); + }) + .on('cycle', function (event) { + console.log(String(event.target)); + }) + .on('complete', function () { + const fastest = this.filter('fastest').map('name'); + this.forEach((bench) => { + results.push({ label: bench.name, avg: 1000 / bench.hz }); + }); + printWinner(results); + }) + .run({ async: false }); // 2. mmul vs mmulNew - - const mmulNewAvg = benchmark( - () => { + const mmulSuite = new Benchmark.Suite('mmul'); + mmulSuite + .add('mmulNew', function () { A.mmul(B); - }, - 'mmulNew', - nbIterations, - ); - - const mmulAvg = benchmark( - () => { + }) + .add('mmul', function () { AOld.mmul(BOld); - }, - 'mmul', - nbIterations, - ); - - const denseAvg = benchmark( - () => { + }) + .add('denseMatrix', function () { denseA.mmul(denseB); - }, - 'denseMatrix', - nbIterations, - ); - - printWinner('mmul', mmulAvg, 'mmulNew', mmulNewAvg); + }) + .on('cycle', function (event) { + console.log(String(event.target)); + }) + .on('complete', function () { + const mmulResults = []; + this.forEach((bench) => { + mmulResults.push({ label: bench.name, avg: 1000 / bench.hz }); + }); + printWinner(mmulResults); + }) + .run({ async: false }); // 3. kroneckerProduct vs kroneckerProductNew - // const kronNewAvg = benchmark(() => { - // A.kroneckerProduct(B); - // }, 'kroneckerProductNew'); - // const kronAvg = benchmark(() => { - // AOld.kroneckerProduct(BOld); - // }, 'kroneckerProduct'); - - // printWinner('kroneckerProduct', kronAvg, 'kroneckerProductNew', kronNewAvg); + const kronSuite = new Benchmark.Suite('kroneckerProduct'); + kronSuite + .add('kroneckerProductNew', function () { + A.kroneckerProduct(B); + }) + .add('kroneckerProduct', function () { + AOld.kroneckerProduct(BOld); + }) + .add('kroneckerProductDense', function () { + denseA.kroneckerProduct(denseB); + }) + .on('cycle', function (event) { + console.log(String(event.target)); + }) + .on('complete', function () { + const kronResults = []; + this.forEach((bench) => { + kronResults.push({ label: bench.name, avg: 1000 / bench.hz }); + }); + printWinner(kronResults); + }) + .run({ async: false }); // 4. matrix multiplication - // const mulAvg = benchmark(() => { - // A.mul(5); - // }, 'mul'); - - // const mulNewAvg = benchmark(() => { - // AOld.mul(5); - // }, 'mulNew'); - - // printWinner('mul', mulAvg, 'mulNew', mulNewAvg); -} - -function runSizeSweepBenchmark() { - const nbIterations = 3; - const sizes = [32, 64, 128, 256]; - const densities = [0.01, 0.015, 0.02, 0.025, 0.03]; - const results = []; - - for (const densityA of densities) { - for (const densityB of densities) { - for (const m of sizes) { - for (const n of sizes) { - for (const p of sizes) { - // A: m x n, B: n x p - - const A = randomSparseMatrix(m, n, densityA); - const B = randomSparseMatrix(n, p, densityB); - - let denseA = A.to2DArray(); - let denseB = B.to2DArray(); - - const AOld = new SparseMatrixOld(denseA); - const BOld = new SparseMatrixOld(denseB); - - denseA = new Matrix(denseA); - denseB = new Matrix(denseB); - - const mmulNewAvg = benchmark( - () => { - A.mmul(B); - }, - 'mmulNew', - nbIterations, - ); - - const mmulAvg = benchmark( - () => { - AOld.mmul(BOld); - }, - 'mmul', - nbIterations, - ); - - const denseAvg = benchmark( - () => { - denseA.mmul(denseB); - }, - 'denseMatrix', - nbIterations, - ); - - results.push({ - densityA, - densityB, - A_shape: [m, n], - B_shape: [n, p], - dense: denseAvg, - new: mmulNewAvg, - old: mmulAvg, - }); - } - } - } - } - } - - fs.writeFileSync( - './benchmark/size_sweep_results-dense.json', - JSON.stringify(results, null, 2), - ); - console.log('Size sweep benchmark results saved to size_sweep_results.json'); + const mulSuite = new Benchmark.Suite('mul'); + mulSuite + .add('mul', function () { + A.mul(5); + }) + .add('mulNew', function () { + AOld.mul(5); + }) + .on('cycle', function (event) { + console.log(String(event.target)); + }) + .on('complete', function () { + const mulResults = []; + this.forEach((bench) => { + mulResults.push({ label: bench.name, avg: 1000 / bench.hz }); + }); + printWinner(mulResults); + }) + .run({ async: false }); } -// runBenchmarks(); -// Uncomment to run the size sweep benchmark -runSizeSweepBenchmark(); +runBenchmarks(); From 77eff3abd2547b3e9051e6d60ba7e5ccca16c7ad Mon Sep 17 00:00:00 2001 From: jobo322 Date: Mon, 23 Jun 2025 11:10:41 -0500 Subject: [PATCH 16/33] refactor: remove add method from SparseMatrix class --- src/index.js | 17 +---------------- 1 file changed, 1 insertion(+), 16 deletions(-) diff --git a/src/index.js b/src/index.js index 23e6c4c..02ab2c4 100644 --- a/src/index.js +++ b/src/index.js @@ -156,28 +156,13 @@ export class SparseMatrix { return this; } - add(other) { - if (typeof other === 'number') return this.addS(other); - if (this.rows !== other.rows || this.columns !== other.columns) { - throw new RangeError('Matrices dimensions must be equal'); - } - - other.withEachNonZero((i, j, v) => { - this.set(i, j, v + this.get(i, j)); - }); - - return this; - } - mul(other) { if (typeof other !== 'number') { throw new RangeError('the argument should be a number'); } if (other === 0) { - return new SparseMatrix(this.rows, this.columns, { - initialCapacity: this.cardinality, - }); + return new SparseMatrix(this.rows, this.columns); } this.withEachNonZero((i, j, v) => { From c05009002dab411deb5e4f9a24eb0fa29ee26897 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Mon, 23 Jun 2025 11:33:42 -0500 Subject: [PATCH 17/33] refactor: remove mul method from SparseMatrix class and update mulM method --- package.json | 1 + src/__tests__/getNonZeros.test.js | 3 +- src/__tests__/matrix.test.js | 198 ------------------------------ src/index.js | 31 ++--- 4 files changed, 14 insertions(+), 219 deletions(-) delete mode 100644 src/__tests__/matrix.test.js diff --git a/package.json b/package.json index df9f8e2..0c52aac 100644 --- a/package.json +++ b/package.json @@ -38,6 +38,7 @@ "devDependencies": { "@types/node": "^22.15.21", "@vitest/coverage-v8": "^3.1.4", + "benchmark": "^2.1.4", "eslint": "^9.27.0", "eslint-config-cheminfo-typescript": "^18.0.1", "ml-matrix": "^6.12.1", diff --git a/src/__tests__/getNonZeros.test.js b/src/__tests__/getNonZeros.test.js index 1b93740..6a0efd7 100644 --- a/src/__tests__/getNonZeros.test.js +++ b/src/__tests__/getNonZeros.test.js @@ -1,4 +1,5 @@ -import { SparseMatrix } from '../index'; +import { describe, it, expect } from 'vitest'; +import { SparseMatrix } from '../index.js'; describe('Sparse Matrix', () => { it('getNonZeros', () => { diff --git a/src/__tests__/matrix.test.js b/src/__tests__/matrix.test.js deleted file mode 100644 index 1a6b38a..0000000 --- a/src/__tests__/matrix.test.js +++ /dev/null @@ -1,198 +0,0 @@ -import { describe, expect, it } from 'vitest'; - -import { SparseMatrix } from '../index.js'; - -describe('Sparse Matrix', () => { - it('add', () => { - let m1 = new SparseMatrix([ - [2, 0, 1], - [0, 0, 3], - [2, 0, 1], - ]); - let m2 = new SparseMatrix([ - [0, 1, 5], - [2, 0, 0], - [-2, 0, -1], - ]); - let m3 = m1.add(m2).to2DArray(); - expect(m3).toStrictEqual([ - [2, 1, 6], - [2, 0, 3], - [0, 0, 0], - ]); - }); - it('mmul', () => { - let m1 = new SparseMatrix([ - [2, 0, 1], - [0, 0, 3], - ]); - let m2 = new SparseMatrix([ - [0, 1], - [2, 0], - [0, 0], - ]); - let m3 = m1.mmul(m2); - - expect(m1.cardinality).toBe(3); - expect(m2.cardinality).toBe(2); - expect(m3.cardinality).toBe(1); - - expect(m3.get(0, 1)).toBe(2); - expect(m3.to2DArray()).toStrictEqual([ - [0, 2], - [0, 0], - ]); - - // Compare with dense multiplication - const denseM1 = m1.to2DArray(); - const denseM2 = m2.to2DArray(); - const expectedDense = denseMatrixMultiply(denseM1, denseM2); - expect(m3.to2DArray()).toStrictEqual(expectedDense); - }); - - it('mmul', () => { - const size = 32; - const density = 0.1; - const m1 = randomSparseMatrix(size, size, density); - const m2 = randomSparseMatrix(size, size, density); - let m3 = m1.mmul(m2); - - const denseM1 = m1.to2DArray(); - const denseM2 = m2.to2DArray(); - - const newSparse = new SparseMatrix(denseM1); - expect(newSparse.to2DArray()).toStrictEqual(denseM1); - const expectedDense = denseMatrixMultiply(denseM1, denseM2); - - expect(m3.to2DArray()).toStrictEqual(expectedDense); - }); - - it('kronecker', () => { - const matrix1 = new SparseMatrix([ - [1, 2], - [3, 4], - ]); - const matrix2 = new SparseMatrix([ - [0, 5], - [6, 7], - ]); - const product = matrix1.kroneckerProduct(matrix2); - expect(product.to2DArray()).toStrictEqual([ - [0, 5, 0, 10], - [6, 7, 12, 14], - [0, 15, 0, 20], - [18, 21, 24, 28], - ]); - }); - - it('isSymmetric', () => { - expect(new SparseMatrix(10, 10).isSymmetric()).toBe(true); - expect(new SparseMatrix(15, 10).isSymmetric()).toBe(false); - - let m = new SparseMatrix([ - [0, 1], - [1, 0], - ]); - expect(m.isSymmetric()).toBe(true); - - m = new SparseMatrix([ - [0, 1], - [0, 1], - ]); - expect(m.isSymmetric()).toBe(false); - }); - - it('transpose', () => { - const matrix = new SparseMatrix([ - [1, 2], - [3, 4], - ]); - expect(matrix.transpose().to2DArray()).toStrictEqual([ - [1, 3], - [2, 4], - ]); - }); -}); - -describe('Banded matrices', () => { - it('Check band size', () => { - const matrix1 = new SparseMatrix([ - [1, 0], - [0, 1], - ]); - const matrix2 = new SparseMatrix([ - [1, 0, 0], - [0, 1, 0], - ]); - const matrix3 = new SparseMatrix([ - [1, 0, 1], - [0, 1, 0], - ]); - const matrix4 = new SparseMatrix([ - [1, 0, 0], - [1, 1, 0], - ]); - const matrix5 = new SparseMatrix([ - [0, 0, 0], - [1, 0, 0], - [0, 1, 0], - ]); - expect(matrix1.bandWidth()).toBe(0); - expect(matrix2.bandWidth()).toBe(0); - expect(matrix3.bandWidth()).toBe(2); - expect(matrix4.bandWidth()).toBe(1); - expect(matrix5.bandWidth()).toBe(0); - }); - - it('isBanded', () => { - const matrix1 = new SparseMatrix([ - [1, 0], - [0, 1], - ]); - const matrix2 = new SparseMatrix([ - [1, 0, 0], - [0, 1, 0], - ]); - const matrix3 = new SparseMatrix([ - [1, 0, 1], - [0, 1, 0], - ]); - const matrix4 = new SparseMatrix([ - [1, 0, 0], - [1, 1, 0], - ]); - expect(matrix1.isBanded(1)).toBe(true); - expect(matrix2.isBanded(1)).toBe(true); - expect(matrix3.isBanded(1)).toBe(false); - expect(matrix4.isBanded(1)).toBe(true); - }); -}); - -function denseMatrixMultiply(A, B) { - const rowsA = A.length; - const colsA = A[0].length; - const colsB = B[0].length; - const result = Array.from({ length: rowsA }, () => Array(colsB).fill(0)); - for (let i = 0; i < rowsA; i++) { - for (let j = 0; j < colsB; j++) { - for (let k = 0; k < colsA; k++) { - result[i][j] += A[i][k] * B[k][j]; - } - } - } - return result; -} - -function randomSparseMatrix(rows, cols, density = 0.01) { - const matrix = []; - for (let i = 0; i < rows; i++) { - const row = new Float64Array(cols); - for (let j = 0; j < cols; j++) { - if (Math.random() < density) { - row[j] = Math.random() * 10; - } - } - matrix.push(row); - } - return new SparseMatrix(matrix); -} diff --git a/src/index.js b/src/index.js index 02ab2c4..66029c2 100644 --- a/src/index.js +++ b/src/index.js @@ -156,22 +156,6 @@ export class SparseMatrix { return this; } - mul(other) { - if (typeof other !== 'number') { - throw new RangeError('the argument should be a number'); - } - - if (other === 0) { - return new SparseMatrix(this.rows, this.columns); - } - - this.withEachNonZero((i, j, v) => { - this.set(i, j, v * other); - }); - - return this; - } - mmul(other) { if (this.columns !== other.rows) { // eslint-disable-next-line no-console @@ -510,10 +494,18 @@ export class SparseMatrix { * @returns {this} */ mulM(matrix) { - matrix.forEachNonZero((i, j, v) => { - this.set(i, j, this.get(i, j) * v); - return v; + if (typeof matrix !== 'number') { + throw new RangeError('the argument should be a number'); + } + + if (matrix === 0) { + this.elements = new HashTable(); + } + + this.withEachNonZero((i, j, v) => { + this.set(i, j, v * matrix); }); + return this; } @@ -1503,7 +1495,6 @@ function csrToCsc(csrMatrix, numCols) { function cooToCsr(cooMatrix, nbRows = 9) { const { values, columns, rows } = cooMatrix; - //could not be the same length const csrRowPtr = new Float64Array(nbRows + 1); const length = values.length; let currentRow = rows[0]; From 44b982296009da811c943bb926a8059d5b4a449c Mon Sep 17 00:00:00 2001 From: jobo322 Date: Tue, 24 Jun 2025 10:44:59 -0500 Subject: [PATCH 18/33] chore: use mitota --- benchmark/benchmark.js | 125 ----------------------- benchmark/benchmarkLinePlot.js | 48 +++++++++ benchmark/benchmarkMitata.js | 45 ++++++++ benchmark/class/SparseMatrixOld.js | 2 +- benchmark/denseMatrix.js | 18 ---- benchmark/sparseMatrix.byhand.js | 41 -------- benchmark/sparseMatrix.js | 30 ------ benchmark/squareMatrixBenchmark.js | 64 ------------ benchmark/utils/benchmark.js | 16 --- benchmark/utils/printWinner.js | 27 ----- benchmark/utils/randomMatrix.js | 22 ++++ benchmark/variateSizeDensityBenchmark.js | 70 ------------- package.json | 2 +- src/index.js | 2 +- 14 files changed, 118 insertions(+), 394 deletions(-) delete mode 100644 benchmark/benchmark.js create mode 100644 benchmark/benchmarkLinePlot.js create mode 100644 benchmark/benchmarkMitata.js delete mode 100644 benchmark/denseMatrix.js delete mode 100644 benchmark/sparseMatrix.byhand.js delete mode 100644 benchmark/sparseMatrix.js delete mode 100644 benchmark/squareMatrixBenchmark.js delete mode 100644 benchmark/utils/benchmark.js delete mode 100644 benchmark/utils/printWinner.js create mode 100644 benchmark/utils/randomMatrix.js delete mode 100644 benchmark/variateSizeDensityBenchmark.js diff --git a/benchmark/benchmark.js b/benchmark/benchmark.js deleted file mode 100644 index 9836b52..0000000 --- a/benchmark/benchmark.js +++ /dev/null @@ -1,125 +0,0 @@ -import Benchmark from 'benchmark'; -import { Matrix } from 'ml-matrix'; -import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; -import { randomSparseMatrix } from './utils/randomSparseMatrix.js'; -import { printWinner } from './utils/printWinner.js'; - -function runBenchmarks() { - const m = 256; - const n = 256; - const p = 256; - const densityA = 0.01; - const densityB = 0.01; - - const A = randomSparseMatrix(m, n, densityA); - const B = randomSparseMatrix(n, p, densityB); - - let denseA = A.to2DArray(); - let denseB = B.to2DArray(); - - const AOld = new SparseMatrixOld(denseA); - const BOld = new SparseMatrixOld(denseB); - - denseA = new Matrix(denseA); - denseB = new Matrix(denseB); - - const results = []; - - // 1. add vs addNew - const addSuite = new Benchmark.Suite('add'); - addSuite - .add('add', function () { - const a = AOld.clone(); - a.add(BOld); - }) - .add('addNew', function () { - const a = A.clone(); - a.add(B); - }) - .add('addDense', function () { - const a = denseA.clone(); - a.add(denseB); - }) - .on('cycle', function (event) { - console.log(String(event.target)); - }) - .on('complete', function () { - const fastest = this.filter('fastest').map('name'); - this.forEach((bench) => { - results.push({ label: bench.name, avg: 1000 / bench.hz }); - }); - printWinner(results); - }) - .run({ async: false }); - - // 2. mmul vs mmulNew - const mmulSuite = new Benchmark.Suite('mmul'); - mmulSuite - .add('mmulNew', function () { - A.mmul(B); - }) - .add('mmul', function () { - AOld.mmul(BOld); - }) - .add('denseMatrix', function () { - denseA.mmul(denseB); - }) - .on('cycle', function (event) { - console.log(String(event.target)); - }) - .on('complete', function () { - const mmulResults = []; - this.forEach((bench) => { - mmulResults.push({ label: bench.name, avg: 1000 / bench.hz }); - }); - printWinner(mmulResults); - }) - .run({ async: false }); - - // 3. kroneckerProduct vs kroneckerProductNew - const kronSuite = new Benchmark.Suite('kroneckerProduct'); - kronSuite - .add('kroneckerProductNew', function () { - A.kroneckerProduct(B); - }) - .add('kroneckerProduct', function () { - AOld.kroneckerProduct(BOld); - }) - .add('kroneckerProductDense', function () { - denseA.kroneckerProduct(denseB); - }) - .on('cycle', function (event) { - console.log(String(event.target)); - }) - .on('complete', function () { - const kronResults = []; - this.forEach((bench) => { - kronResults.push({ label: bench.name, avg: 1000 / bench.hz }); - }); - printWinner(kronResults); - }) - .run({ async: false }); - - // 4. matrix multiplication - const mulSuite = new Benchmark.Suite('mul'); - mulSuite - .add('mul', function () { - A.mul(5); - }) - .add('mulNew', function () { - AOld.mul(5); - }) - .on('cycle', function (event) { - console.log(String(event.target)); - }) - .on('complete', function () { - const mulResults = []; - this.forEach((bench) => { - mulResults.push({ label: bench.name, avg: 1000 / bench.hz }); - }); - printWinner(mulResults); - }) - .run({ async: false }); -} - -runBenchmarks(); diff --git a/benchmark/benchmarkLinePlot.js b/benchmark/benchmarkLinePlot.js new file mode 100644 index 0000000..91f0224 --- /dev/null +++ b/benchmark/benchmarkLinePlot.js @@ -0,0 +1,48 @@ +import { run, bench, lineplot, do_not_optimize } from 'mitata'; +import { SparseMatrix } from '../src/index.js'; +import { randomSparseMatrix } from './utils/randomSparseMatrix.js'; +import { Matrix } from 'ml-matrix'; +import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; +import { randomMatrix } from './utils/randomMatrix.js'; +const density = 0.01; // Fixed density for this comparison; +const min = 32; +const max = 128; +// Prepare matrices once + +lineplot(() => { + bench('Sparse.mmul($size)', function* (ctx) { + const size = ctx.get('size'); + + // Prepare matrices once + const A = new SparseMatrix(randomMatrix(size, size, density)); + const B = new SparseMatrix(randomMatrix(size, size, density)); + // Benchmark the multiplication + yield () => do_not_optimize(A.mmul(B)); + }).range('size', min, max, 2); // 16, 32, 64, 128, 256 + + bench('SparseOld.mmul($size)', function* (ctx) { + const size = ctx.get('size'); + const A = randomMatrix(size, size, density); + const B = randomMatrix(size, size, density); + const AOld = new SparseMatrixOld(A); + const BOld = new SparseMatrixOld(B); + + // Benchmark the multiplication + yield () => do_not_optimize(AOld.mmul(BOld)); + }).range('size', min, max, 2); + + bench('Dense.mmul($size)', function* (ctx) { + const size = ctx.get('size'); + + // Prepare matrices once + const A = randomMatrix(size, size, density); + const B = randomMatrix(size, size, density); + const ADense = new Matrix(A); + const BDense = new Matrix(B); + + // Benchmark the multiplication + yield () => do_not_optimize(ADense.mmul(BDense)); + }).range('size', min, max, 2); +}); + +await run(); diff --git a/benchmark/benchmarkMitata.js b/benchmark/benchmarkMitata.js new file mode 100644 index 0000000..9a0ed7d --- /dev/null +++ b/benchmark/benchmarkMitata.js @@ -0,0 +1,45 @@ +import { run, bench, group, do_not_optimize } from 'mitata'; +import { Matrix } from 'ml-matrix'; +import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld2.js'; +import fs from 'fs'; +import { randomSparseMatrix } from './utils/randomSparseMatrix.js'; + +const sizes = [64, 128, 256]; +const densities = [0.01, 0.015, 0.02, 0.025, 0.03]; +const results = []; + +for (const density of densities) { + for (const size of sizes) { + const A = randomSparseMatrix(size, size, density); + const B = randomSparseMatrix(size, size, density); + let denseA = A.to2DArray(); + let denseB = B.to2DArray(); + const AOld = new SparseMatrixOld(denseA); + const BOld = new SparseMatrixOld(denseB); + denseA = new Matrix(denseA); + denseB = new Matrix(denseB); + + // Warm up + A.mmul(B); + + let mmulNewAvg, mmulAvg, denseAvg; + + group(`size:${size}-density:${density}`, () => { + bench('mmulNew', () => { + do_not_optimize(A.mmul(B)); + }); //.gc('inner'); + bench('mmul', () => { + do_not_optimize(AOld.mmul(BOld)); + }); //.gc('inner'); + bench('denseMatrix', () => { + do_not_optimize(denseA.mmul(denseB)); + }); //.gc('inner'); + }); + + // mitata will handle timing and reporting + // You can extract results from mitata's output or use the save() function + } +} + +await run({ silent: false }); +// await save({ file: './benchmark/mitata_results.json', format: 'json' }); diff --git a/benchmark/class/SparseMatrixOld.js b/benchmark/class/SparseMatrixOld.js index f770b91..e7ea2d3 100644 --- a/benchmark/class/SparseMatrixOld.js +++ b/benchmark/class/SparseMatrixOld.js @@ -1383,4 +1383,4 @@ export class SparseMatrix { SparseMatrix.prototype.klass = 'Matrix'; SparseMatrix.identity = SparseMatrix.eye; -SparseMatrix.prototype.tensorProduct = SparseMatrix.prototype.kroneckerProduct; +SparseMatrix.prototype.tensorProduct = SparseMatrix.prototype.kroneckerProduct; \ No newline at end of file diff --git a/benchmark/denseMatrix.js b/benchmark/denseMatrix.js deleted file mode 100644 index afab679..0000000 --- a/benchmark/denseMatrix.js +++ /dev/null @@ -1,18 +0,0 @@ -import { Matrix } from 'ml-matrix'; - -let size = 200; - -console.time('dense'); -const matrix = new Matrix(size, size).fill(1); -const matrix2 = new Matrix(size, size).fill(1); - -const product = matrix.mmul(matrix2); -// sum of all the elements -let sum = 0; -for (let i = 0; i < size; i++) { - for (let j = 0; j < size; j++) { - sum += product.get(i, j); - } -} -console.log(sum); -console.timeEnd('dense'); diff --git a/benchmark/sparseMatrix.byhand.js b/benchmark/sparseMatrix.byhand.js deleted file mode 100644 index 8261900..0000000 --- a/benchmark/sparseMatrix.byhand.js +++ /dev/null @@ -1,41 +0,0 @@ -import { SparseMatrix } from '../lib/index.js'; - -let size = 200; - -console.time('dense'); - -const matrix1 = new SparseMatrix(size, size); -fillSparseMatrix(matrix1); - -const matrix2 = new SparseMatrix(size, size); -fillSparseMatrix(matrix2); - -const product = new SparseMatrix(size, size); -for (let i = 0; i < size; i++) { - for (let j = 0; j < size; j++) { - for (let k = 0; k < size; k++) { - product.set( - i, - j, - product.get(i, j) + matrix1.get(i, k) * matrix2.get(k, j), - ); - } - } -} - -let sum = 0; -for (let i = 0; i < size; i++) { - for (let j = 0; j < size; j++) { - sum += product.get(i, j); - } -} -console.timeEnd('dense'); -console.log(sum); - -function fillSparseMatrix(matrix) { - for (let row = 0; row < matrix.rows; row++) { - for (let column = 0; column < matrix.columns; column++) { - matrix.set(row, column, 1); - } - } -} diff --git a/benchmark/sparseMatrix.js b/benchmark/sparseMatrix.js deleted file mode 100644 index 53b117e..0000000 --- a/benchmark/sparseMatrix.js +++ /dev/null @@ -1,30 +0,0 @@ -import { SparseMatrix } from '../lib/index.js'; - -let size = 200; - -console.time('dense'); - -const matrix = new SparseMatrix(size, size); -fillSparseMatrix(matrix); - -const matrix2 = new SparseMatrix(size, size); -fillSparseMatrix(matrix2); - -const product = matrix.mmul(matrix2); -// sum of all the elements -let sum = 0; -for (let i = 0; i < size; i++) { - for (let j = 0; j < size; j++) { - sum += product.get(i, j); - } -} -console.timeEnd('dense'); -console.log(sum); - -function fillSparseMatrix(matrix) { - for (let row = 0; row < matrix.rows; row++) { - for (let column = 0; column < matrix.columns; column++) { - matrix.set(row, column, 1); - } - } -} diff --git a/benchmark/squareMatrixBenchmark.js b/benchmark/squareMatrixBenchmark.js deleted file mode 100644 index 1e4e571..0000000 --- a/benchmark/squareMatrixBenchmark.js +++ /dev/null @@ -1,64 +0,0 @@ -import Benchmark from 'benchmark'; -import { Matrix } from 'ml-matrix'; -import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; -import fs from 'fs'; -import { randomSparseMatrix } from './utils/randomSparseMatrix.js'; - -function runSizeSweepBenchmarkCompact() { - const sizes = [32, 64, 128, 256]; - const densities = [0.01, 0.015, 0.02, 0.025, 0.03]; - const results = []; - - for (const densityA of densities) { - for (const m of sizes) { - const A = randomSparseMatrix(m, m, densityA); - const B = randomSparseMatrix(m, m, densityA); - let denseA = A.to2DArray(); - let denseB = B.to2DArray(); - const AOld = new SparseMatrixOld(denseA); - const BOld = new SparseMatrixOld(denseB); - denseA = new Matrix(denseA); - denseB = new Matrix(denseB); - // Warm up - A.mmul(B); - // Use Benchmark.js for each method - const suite = new Benchmark.Suite(); - let mmulNewAvg, mmulAvg, denseAvg; - suite - .add('mmulNew', function () { - A.mmul(B); - }) - .add('mmul', function () { - AOld.mmul(BOld); - }) - .add('denseMatrix', function () { - denseA.mmul(denseB); - }) - .on('cycle', function (event) { - // Optionally log: console.log(String(event.target)); - }) - .on('complete', function () { - this.forEach((bench) => { - if (bench.name === 'mmulNew') mmulNewAvg = 1000 / bench.hz; - if (bench.name === 'mmul') mmulAvg = 1000 / bench.hz; - if (bench.name === 'denseMatrix') denseAvg = 1000 / bench.hz; - }); - }) - .run({ async: false }); - results.push({ - density: densityA, - size: m, - dense: denseAvg, - new: mmulNewAvg, - old: mmulAvg, - }); - } - } - fs.writeFileSync( - './benchmark/size_sweep_results-dense.json', - JSON.stringify(results, null, 2), - ); - console.log('Size sweep benchmark results saved to size_sweep_results.json'); -} - -runSizeSweepBenchmarkCompact(); diff --git a/benchmark/utils/benchmark.js b/benchmark/utils/benchmark.js deleted file mode 100644 index cc60e1c..0000000 --- a/benchmark/utils/benchmark.js +++ /dev/null @@ -1,16 +0,0 @@ -export function benchmark(fn, label, iterations = 5, logIt = true) { - const times = []; - for (let i = 0; i < iterations; i++) { - const t0 = performance.now(); - fn(); - const t1 = performance.now(); - times.push(t1 - t0); - } - const avg = times.reduce((a, b) => a + b, 0) / times.length; - if (logIt) { - console.log( - `${label}: avg ${avg.toFixed(2)} ms over ${iterations} runs \n`, - ); - } - return avg; -} diff --git a/benchmark/utils/printWinner.js b/benchmark/utils/printWinner.js deleted file mode 100644 index 527e686..0000000 --- a/benchmark/utils/printWinner.js +++ /dev/null @@ -1,27 +0,0 @@ -export function printWinner(results) { - // results: Array of { label: string, avg: number } - if (!Array.isArray(results) || results.length < 2) { - console.log('Need at least two results to compare.'); - return; - } - - // Sort by avg ascending (fastest first) - const sorted = [...results].sort((a, b) => a.avg - b.avg); - const winner = sorted[0]; - - console.log( - ` -> ${winner.label} was the fastest (${winner.avg.toFixed(2)} ms avg)\n`, - ); - - for (let i = 1; i < sorted.length; i++) { - const slower = sorted[i]; - const speedup = (slower.avg / winner.avg).toFixed(2); - const percent = (((slower.avg - winner.avg) / slower.avg) * 100).toFixed(1); - console.log( - ` ${winner.label} was ${speedup}x faster (${percent}% faster) than ${ - slower.label - } (${slower.avg.toFixed(2)} ms avg)`, - ); - } - console.log(); -} diff --git a/benchmark/utils/randomMatrix.js b/benchmark/utils/randomMatrix.js new file mode 100644 index 0000000..3db3b66 --- /dev/null +++ b/benchmark/utils/randomMatrix.js @@ -0,0 +1,22 @@ +export function randomMatrix(rows, cols, density = 0.01) { + const total = rows * cols; + const cardinality = Math.ceil(total * density); + const positions = new Set(); + + // Generate unique random positions + while (positions.size < cardinality) { + positions.add(Math.floor(Math.random() * total)); + } + + // Build the matrix with zeros + const matrix = Array.from({ length: rows }, () => new Float64Array(cols)); + + // Assign random values to the selected positions + for (const pos of positions) { + const i = Math.floor(pos / cols); + const j = pos % cols; + matrix[i][j] = Math.random() * 10; + } + + return matrix; +} diff --git a/benchmark/variateSizeDensityBenchmark.js b/benchmark/variateSizeDensityBenchmark.js deleted file mode 100644 index 39fe24c..0000000 --- a/benchmark/variateSizeDensityBenchmark.js +++ /dev/null @@ -1,70 +0,0 @@ -import Benchmark from 'benchmark'; -import { Matrix } from 'ml-matrix'; -import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; -import fs from 'fs'; -import { randomSparseMatrix } from './utils/randomSparseMatrix.js'; - -function runSizeSweepBenchmark() { - const sizes = [32, 64, 128, 256]; - const densities = [0.01, 0.015, 0.02, 0.025, 0.03]; - const results = []; - - for (const densityA of densities) { - for (const densityB of densities) { - for (const m of sizes) { - for (const n of sizes) { - for (const p of sizes) { - const A = randomSparseMatrix(m, n, densityA); - const B = randomSparseMatrix(n, p, densityB); - let denseA = A.to2DArray(); - let denseB = B.to2DArray(); - const AOld = new SparseMatrixOld(denseA); - const BOld = new SparseMatrixOld(denseB); - denseA = new Matrix(denseA); - denseB = new Matrix(denseB); - // Use Benchmark.js for each method - const suite = new Benchmark.Suite(); - let mmulNewAvg, mmulAvg, denseAvg; - suite - .add('mmulNew', function () { - A.mmul(B); - }) - .add('mmul', function () { - AOld.mmul(BOld); - }) - .add('denseMatrix', function () { - denseA.mmul(denseB); - }) - .on('cycle', function (event) { - // Optionally log: console.log(String(event.target)); - }) - .on('complete', function () { - this.forEach((bench) => { - if (bench.name === 'mmulNew') mmulNewAvg = 1000 / bench.hz; - if (bench.name === 'mmul') mmulAvg = 1000 / bench.hz; - if (bench.name === 'denseMatrix') denseAvg = 1000 / bench.hz; - }); - }) - .run({ async: false }); - results.push({ - densityA, - densityB, - A_shape: [m, n], - B_shape: [n, p], - dense: denseAvg, - new: mmulNewAvg, - old: mmulAvg, - }); - } - } - } - } - } - fs.writeFileSync( - './benchmark/size_sweep_results-dense.json', - JSON.stringify(results, null, 2), - ); - console.log('Size sweep benchmark results saved to size_sweep_results.json'); -} - -runSizeSweepBenchmark(); diff --git a/package.json b/package.json index 0c52aac..8f4668c 100644 --- a/package.json +++ b/package.json @@ -38,9 +38,9 @@ "devDependencies": { "@types/node": "^22.15.21", "@vitest/coverage-v8": "^3.1.4", - "benchmark": "^2.1.4", "eslint": "^9.27.0", "eslint-config-cheminfo-typescript": "^18.0.1", + "mitata": "^1.0.34", "ml-matrix": "^6.12.1", "prettier": "^3.5.3", "rimraf": "^6.0.1", diff --git a/src/index.js b/src/index.js index 66029c2..7e6615d 100644 --- a/src/index.js +++ b/src/index.js @@ -209,7 +209,7 @@ export class SparseMatrix { const q = other.columns; const result = new SparseMatrix(m * p, n * q, { - initialCapacity: this.cardinality * other.cardinality + 20, + initialCapacity: this.cardinality * other.cardinality + 10, }); const { From 22d602d9ca7edadc7ce90ee42c5f849ef441d430 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Tue, 24 Jun 2025 11:25:49 -0500 Subject: [PATCH 19/33] Chore: use default value of initialCapacity from ml-hast-table --- benchmark/benchmarkLinePlot.js | 36 +++++++++++++++++----------------- benchmark/benchmarkMitata.js | 11 ----------- package.json | 1 + src/index.js | 3 +-- 4 files changed, 20 insertions(+), 31 deletions(-) diff --git a/benchmark/benchmarkLinePlot.js b/benchmark/benchmarkLinePlot.js index 91f0224..07c5f77 100644 --- a/benchmark/benchmarkLinePlot.js +++ b/benchmark/benchmarkLinePlot.js @@ -1,14 +1,14 @@ import { run, bench, lineplot, do_not_optimize } from 'mitata'; import { SparseMatrix } from '../src/index.js'; -import { randomSparseMatrix } from './utils/randomSparseMatrix.js'; -import { Matrix } from 'ml-matrix'; +import { xSequentialFillFromStep } from 'ml-spectra-processing'; import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; import { randomMatrix } from './utils/randomMatrix.js'; -const density = 0.01; // Fixed density for this comparison; -const min = 32; -const max = 128; -// Prepare matrices once +const density = 0.02; // Fixed density for this comparison; +// Prepare matrices once +const sizes = Array.from( + xSequentialFillFromStep({ from: 4, step: 4, size: 13 }), +); lineplot(() => { bench('Sparse.mmul($size)', function* (ctx) { const size = ctx.get('size'); @@ -18,7 +18,7 @@ lineplot(() => { const B = new SparseMatrix(randomMatrix(size, size, density)); // Benchmark the multiplication yield () => do_not_optimize(A.mmul(B)); - }).range('size', min, max, 2); // 16, 32, 64, 128, 256 + }).args('size', sizes); // 16, 32, 64, 128, 256 bench('SparseOld.mmul($size)', function* (ctx) { const size = ctx.get('size'); @@ -29,20 +29,20 @@ lineplot(() => { // Benchmark the multiplication yield () => do_not_optimize(AOld.mmul(BOld)); - }).range('size', min, max, 2); + }).args('size', sizes); - bench('Dense.mmul($size)', function* (ctx) { - const size = ctx.get('size'); + // bench('Dense.mmul($size)', function* (ctx) { + // const size = ctx.get('size'); - // Prepare matrices once - const A = randomMatrix(size, size, density); - const B = randomMatrix(size, size, density); - const ADense = new Matrix(A); - const BDense = new Matrix(B); + // // Prepare matrices once + // const A = randomMatrix(size, size, density); + // const B = randomMatrix(size, size, density); + // const ADense = new Matrix(A); + // const BDense = new Matrix(B); - // Benchmark the multiplication - yield () => do_not_optimize(ADense.mmul(BDense)); - }).range('size', min, max, 2); + // // Benchmark the multiplication + // yield () => do_not_optimize(ADense.mmul(BDense)); + // }).range('size', min, max, multiplier); }); await run(); diff --git a/benchmark/benchmarkMitata.js b/benchmark/benchmarkMitata.js index 9a0ed7d..8d694a5 100644 --- a/benchmark/benchmarkMitata.js +++ b/benchmark/benchmarkMitata.js @@ -1,12 +1,10 @@ import { run, bench, group, do_not_optimize } from 'mitata'; import { Matrix } from 'ml-matrix'; import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld2.js'; -import fs from 'fs'; import { randomSparseMatrix } from './utils/randomSparseMatrix.js'; const sizes = [64, 128, 256]; const densities = [0.01, 0.015, 0.02, 0.025, 0.03]; -const results = []; for (const density of densities) { for (const size of sizes) { @@ -19,11 +17,6 @@ for (const density of densities) { denseA = new Matrix(denseA); denseB = new Matrix(denseB); - // Warm up - A.mmul(B); - - let mmulNewAvg, mmulAvg, denseAvg; - group(`size:${size}-density:${density}`, () => { bench('mmulNew', () => { do_not_optimize(A.mmul(B)); @@ -35,11 +28,7 @@ for (const density of densities) { do_not_optimize(denseA.mmul(denseB)); }); //.gc('inner'); }); - - // mitata will handle timing and reporting - // You can extract results from mitata's output or use the save() function } } await run({ silent: false }); -// await save({ file: './benchmark/mitata_results.json', format: 'json' }); diff --git a/package.json b/package.json index 8f4668c..0801258 100644 --- a/package.json +++ b/package.json @@ -42,6 +42,7 @@ "eslint-config-cheminfo-typescript": "^18.0.1", "mitata": "^1.0.34", "ml-matrix": "^6.12.1", + "ml-spectra-processing": "^14.12.0", "prettier": "^3.5.3", "rimraf": "^6.0.1", "typescript": "^5.8.3", diff --git a/src/index.js b/src/index.js index 7e6615d..8de182d 100644 --- a/src/index.js +++ b/src/index.js @@ -19,8 +19,7 @@ export class SparseMatrix { if (Array.isArray(rows)) { const matrix = rows; const nbRows = matrix.length; - const nbColmuns = matrix[0].length; - options = columns || { initialCapacity: nbRows * nbColmuns }; + options = columns || {}; columns = matrix[0].length; this._init(nbRows, columns, new HashTable(options), options.threshold); for (let i = 0; i < nbRows; i++) { From dd0b69ca735ce24d77fb6bd4fd8f5b36c17899c2 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Tue, 24 Jun 2025 11:51:32 -0500 Subject: [PATCH 20/33] chore: fix eslint and prettier in benchmark --- benchmark/benchmarkLinePlot.js | 9 ++++++++- benchmark/benchmarkMitata.js | 6 ++++++ benchmark/class/SparseMatrixOld.js | 3 +-- 3 files changed, 15 insertions(+), 3 deletions(-) diff --git a/benchmark/benchmarkLinePlot.js b/benchmark/benchmarkLinePlot.js index 07c5f77..b413c2d 100644 --- a/benchmark/benchmarkLinePlot.js +++ b/benchmark/benchmarkLinePlot.js @@ -1,10 +1,17 @@ import { run, bench, lineplot, do_not_optimize } from 'mitata'; -import { SparseMatrix } from '../src/index.js'; import { xSequentialFillFromStep } from 'ml-spectra-processing'; + +import { SparseMatrix } from '../src/index.js'; + import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; import { randomMatrix } from './utils/randomMatrix.js'; + const density = 0.02; // Fixed density for this comparison; +/* eslint +func-names: 0 +camelcase: 0 +*/ // Prepare matrices once const sizes = Array.from( xSequentialFillFromStep({ from: 4, step: 4, size: 13 }), diff --git a/benchmark/benchmarkMitata.js b/benchmark/benchmarkMitata.js index 8d694a5..a99ade4 100644 --- a/benchmark/benchmarkMitata.js +++ b/benchmark/benchmarkMitata.js @@ -1,8 +1,14 @@ import { run, bench, group, do_not_optimize } from 'mitata'; import { Matrix } from 'ml-matrix'; + import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld2.js'; import { randomSparseMatrix } from './utils/randomSparseMatrix.js'; +/* eslint +func-names: 0 +camelcase: 0 +*/ + const sizes = [64, 128, 256]; const densities = [0.01, 0.015, 0.02, 0.025, 0.03]; diff --git a/benchmark/class/SparseMatrixOld.js b/benchmark/class/SparseMatrixOld.js index e7ea2d3..fca3b8e 100644 --- a/benchmark/class/SparseMatrixOld.js +++ b/benchmark/class/SparseMatrixOld.js @@ -159,7 +159,6 @@ export class SparseMatrix { */ mmul(other) { if (this.columns !== other.rows) { - // eslint-disable-next-line no-console console.warn( 'Number of columns of left matrix are not equal to number of rows of right matrix.', ); @@ -1383,4 +1382,4 @@ export class SparseMatrix { SparseMatrix.prototype.klass = 'Matrix'; SparseMatrix.identity = SparseMatrix.eye; -SparseMatrix.prototype.tensorProduct = SparseMatrix.prototype.kroneckerProduct; \ No newline at end of file +SparseMatrix.prototype.tensorProduct = SparseMatrix.prototype.kroneckerProduct; From c32bc1be8e8538c300cd04bd0e4f61054fec49e6 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Tue, 24 Jun 2025 11:51:52 -0500 Subject: [PATCH 21/33] chore: add jsdoc to functions --- src/index.js | 19 ++++++++++++++++++- 1 file changed, 18 insertions(+), 1 deletion(-) diff --git a/src/index.js b/src/index.js index 8de182d..346ba51 100644 --- a/src/index.js +++ b/src/index.js @@ -1458,6 +1458,12 @@ SparseMatrix.prototype.klass = 'Matrix'; SparseMatrix.identity = SparseMatrix.eye; SparseMatrix.prototype.tensorProduct = SparseMatrix.prototype.kroneckerProduct; +/** + * Converts a matrix from Compressed Sparse Row (CSR) format to Compressed Sparse Column (CSC) format. + * @param {Object} csrMatrix - The matrix in CSR format with properties: values, columns, rows. + * @param {number} numCols - The number of columns in the matrix. + * @returns {{rows: Float64Array, columns: Float64Array, values: Float64Array}} The matrix in CSC format. + */ function csrToCsc(csrMatrix, numCols) { const { values: csrValues, @@ -1492,13 +1498,18 @@ function csrToCsc(csrMatrix, numCols) { return { rows: cscRowIndices, columns: cscColPtr, values: cscValues }; } +/** + * Converts a matrix from Coordinate (COO) format to Compressed Sparse Row (CSR) format. + * @param {Object} cooMatrix - The matrix in COO format with properties: values, columns, rows. + * @param {number} [nbRows=9] - The number of rows in the matrix. + * @returns {{rows: Float64Array, columns: Float64Array, values: Float64Array}} The matrix in CSR format. + */ function cooToCsr(cooMatrix, nbRows = 9) { const { values, columns, rows } = cooMatrix; const csrRowPtr = new Float64Array(nbRows + 1); const length = values.length; let currentRow = rows[0]; for (let index = 0; index < length; ) { - const prev = index; while (currentRow === rows[index] && index < length) ++index; csrRowPtr[currentRow + 1] = index; currentRow += 1; @@ -1506,6 +1517,12 @@ function cooToCsr(cooMatrix, nbRows = 9) { return { rows: csrRowPtr, columns, values }; } +/** + * Creates an empty 2D matrix (array of Float64Array) with the given dimensions. + * @param {number} nbRows - Number of rows. + * @param {number} nbColumns - Number of columns. + * @returns {Float64Array[]} A 2D array representing the matrix. + */ function matrixCreateEmpty(nbRows, nbColumns) { const newMatrix = []; for (let row = 0; row < nbRows; row++) { From eac75af40a789b088c46ee902f063c515bb429df Mon Sep 17 00:00:00 2001 From: jobo322 Date: Tue, 24 Jun 2025 22:26:08 -0500 Subject: [PATCH 22/33] chore: exclude benchmark and lib from coverage report --- vitest.config.mjs | 14 ++++++++++++++ 1 file changed, 14 insertions(+) create mode 100644 vitest.config.mjs diff --git a/vitest.config.mjs b/vitest.config.mjs new file mode 100644 index 0000000..0bf34cd --- /dev/null +++ b/vitest.config.mjs @@ -0,0 +1,14 @@ +import { defineConfig, configDefaults } from 'vitest/config'; + +export default defineConfig({ + test: { + coverage: { + exclude: [ + ...configDefaults.coverage.exclude, + 'benchmark', + 'lib', + 'scripts', + ], + }, + }, +}); From 63b64ed942c70184e9ada3753c14eeb6e94bb367 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Tue, 24 Jun 2025 22:33:40 -0500 Subject: [PATCH 23/33] chore: ensure that getNonZeros sort data if format is required --- src/index.js | 21 ++++++++++++++++----- 1 file changed, 16 insertions(+), 5 deletions(-) diff --git a/src/index.js b/src/index.js index 346ba51..6d99c24 100644 --- a/src/index.js +++ b/src/index.js @@ -288,6 +288,17 @@ export class SparseMatrix { return this; } + /** + * Returns the non-zero elements of the matrix in coordinate (COO), CSR, or CSC format. + * + * @param {Object} [options={}] - Options for output format and sorting. + * @param {'csr'|'csc'} [options.format] - If specified, returns the result in CSR or CSC format. Otherwise, returns COO format. + * @param {boolean} [options.sort] - If true, sorts the non-zero elements by their indices. + * @returns {Object} If no format is specified, returns an object with Float64Array `rows`, `columns`, and `values` (COO format). + * If format is 'csr', returns { rows, columns, values } in CSR format. + * If format is 'csc', returns { rows, columns, values } in CSC format. + * @throws {Error} If an unsupported format is specified. + */ getNonZeros(options = {}) { const cardinality = this.cardinality; /** @type {Float64Array} */ @@ -297,7 +308,7 @@ export class SparseMatrix { /** @type {Float64Array} */ const values = new Float64Array(cardinality); - const { format, sort = false } = options; + const { format, sort = format !== undefined } = options; let idx = 0; this.withEachNonZero((i, j, value) => { @@ -305,11 +316,11 @@ export class SparseMatrix { columns[idx] = j; values[idx] = value; idx++; - }, sort || format); + }, sort); if (!format) return { rows, columns, values }; - if (!['csr', 'csc'].includes(format.toLowerCase())) { + if (!['csr', 'csc'].includes(format)) { throw new Error(`format ${format} is not supported`); } @@ -1501,10 +1512,10 @@ function csrToCsc(csrMatrix, numCols) { /** * Converts a matrix from Coordinate (COO) format to Compressed Sparse Row (CSR) format. * @param {Object} cooMatrix - The matrix in COO format with properties: values, columns, rows. - * @param {number} [nbRows=9] - The number of rows in the matrix. + * @param {number} nbRows - The number of rows in the matrix. * @returns {{rows: Float64Array, columns: Float64Array, values: Float64Array}} The matrix in CSR format. */ -function cooToCsr(cooMatrix, nbRows = 9) { +function cooToCsr(cooMatrix, nbRows) { const { values, columns, rows } = cooMatrix; const csrRowPtr = new Float64Array(nbRows + 1); const length = values.length; From 8cef7c05cc86850f8f7359b882ee239dfaf6541e Mon Sep 17 00:00:00 2001 From: jobo322 Date: Tue, 1 Jul 2025 23:07:39 -0500 Subject: [PATCH 24/33] chore: refactor getNonZeros method to simplify format handling and improve performance --- src/__tests__/getNonZeros.test.js | 10 +- src/index.js | 178 +++++++++++++++--------------- 2 files changed, 91 insertions(+), 97 deletions(-) diff --git a/src/__tests__/getNonZeros.test.js b/src/__tests__/getNonZeros.test.js index 6a0efd7..8715982 100644 --- a/src/__tests__/getNonZeros.test.js +++ b/src/__tests__/getNonZeros.test.js @@ -1,4 +1,5 @@ import { describe, it, expect } from 'vitest'; + import { SparseMatrix } from '../index.js'; describe('Sparse Matrix', () => { @@ -22,17 +23,10 @@ describe('Sparse Matrix', () => { }); // CSR format - expect(m2.getNonZeros({ format: 'csr' })).toEqual({ + expect(m2.getNonZeros({ format: true })).toEqual({ rows: Float64Array.from([0, 0, 4, 7, 7, 11]), columns: Float64Array.from([0, 3, 4, 5, 1, 4, 5, 0, 3, 4, 5]), values: Float64Array.from([1, 2, 1, 1, 3, 5, 5, 1, 1, 9, 9]), }); - - //CSC format - expect(m2.getNonZeros({ format: 'csc' })).toEqual({ - rows: Float64Array.from([1, 4, 2, 1, 4, 1, 2, 4, 1, 2, 4]), - columns: Float64Array.from([0, 2, 3, 3, 5, 8, 11]), - values: Float64Array.from([1, 1, 3, 2, 1, 1, 5, 9, 1, 5, 9]), - }); }); }); diff --git a/src/index.js b/src/index.js index 6d99c24..351e3e3 100644 --- a/src/index.js +++ b/src/index.js @@ -163,40 +163,97 @@ export class SparseMatrix { ); } + if (this.cardinality < 42 && other.cardinality < 42) { + return this._mmulSmall(other); + } else if (other.rows > 100 && other.cardinality < 110) { + return this._mmulLowDensity(other); + } + + return this._mmulMediumDensity(other); + } + + _mmulSmall(other) { const m = this.rows; const p = other.columns; + const { + columns: otherCols, + rows: otherRows, + values: otherValues, + } = other.getNonZeros(); - const result = matrixCreateEmpty(m, p); + const nbOtherActive = otherCols.length; + const result = new SparseMatrix(m, p); + this.withEachNonZero((i, j, v1) => { + for (let o = 0; o < nbOtherActive; o++) { + if (j === otherRows[o]) { + const l = otherCols[o]; + result.set(i, l, result.get(i, l) + otherValues[o] * v1); + } + } + }); + return result; + } + _mmulLowDensity(other) { + const m = this.rows; + const p = other.columns; const { columns: otherCols, rows: otherRows, values: otherValues, - } = other.getNonZeros({ format: 'csr' }); + } = other.getNonZeros(); const { columns: thisCols, rows: thisRows, values: thisValues, - } = this.getNonZeros({ format: 'csc' }); - - const thisNbCols = this.columns; - for (let t = 0; t < thisNbCols; t++) { - const tStart = thisCols[t]; - const tEnd = thisCols[t + 1]; - const oStart = otherRows[t]; - const oEnd = otherRows[t + 1]; - for (let f = tStart; f < tEnd; f++) { - for (let k = oStart; k < oEnd; k++) { - const i = thisRows[f]; - const l = otherCols[k]; - result[i][l] += thisValues[f] * otherValues[k]; + } = this.getNonZeros(); + + const result = new SparseMatrix(m, p); + const nbOtherActive = otherCols.length; + const nbThisActive = thisCols.length; + for (let t = 0; t < nbThisActive; t++) { + const i = thisRows[t]; + const j = thisCols[t]; + for (let o = 0; o < nbOtherActive; o++) { + if (j === otherRows[o]) { + const l = otherCols[o]; + result.set(i, l, result.get(i, l) + otherValues[o] * thisValues[t]); } } } - - return new SparseMatrix(result); + // console.log(result.cardinality); + return result; } + _mmulMediumDensity(other) { + const m = this.rows; + const p = other.columns; + const { + columns: otherCols, + rows: otherRows, + values: otherValues, + } = other.getNonZeros({ format: 'csr' }); + const { + columns: thisCols, + rows: thisRows, + values: thisValues, + } = this.getNonZeros(); + + const result = new SparseMatrix(m, p); + const nbThisActive = thisCols.length; + for (let t = 0; t < nbThisActive; t++) { + const i = thisRows[t]; + const j = thisCols[t]; + const oStart = otherRows[j]; + const oEnd = otherRows[j + 1]; + for (let k = oStart; k < oEnd; k++) { + const l = otherCols[k]; + result.set(i, l, result.get(i, l) + otherValues[k] * thisValues[t]); + } + } + + return result; + } /** * @param {SparseMatrix} other * @returns {SparseMatrix} @@ -292,11 +349,10 @@ export class SparseMatrix { * Returns the non-zero elements of the matrix in coordinate (COO), CSR, or CSC format. * * @param {Object} [options={}] - Options for output format and sorting. - * @param {'csr'|'csc'} [options.format] - If specified, returns the result in CSR or CSC format. Otherwise, returns COO format. + * @param {boolean} [options.format] - If specified, returns the result in CSR or CSC format. Otherwise, returns COO format. * @param {boolean} [options.sort] - If true, sorts the non-zero elements by their indices. * @returns {Object} If no format is specified, returns an object with Float64Array `rows`, `columns`, and `values` (COO format). - * If format is 'csr', returns { rows, columns, values } in CSR format. - * If format is 'csc', returns { rows, columns, values } in CSC format. + * If format is true, returns { rows, columns, values } in CSR format. * @throws {Error} If an unsupported format is specified. */ getNonZeros(options = {}) { @@ -318,16 +374,9 @@ export class SparseMatrix { idx++; }, sort); - if (!format) return { rows, columns, values }; - - if (!['csr', 'csc'].includes(format)) { - throw new Error(`format ${format} is not supported`); - } - - const csrMatrix = cooToCsr({ rows, columns, values }, this.rows); - return format.toLowerCase() === 'csc' - ? csrToCsc(csrMatrix, this.columns) - : csrMatrix; + return format + ? cooToCsr({ rows, columns, values }, this.rows) + : { rows, columns, values }; } /** @@ -1469,46 +1518,6 @@ SparseMatrix.prototype.klass = 'Matrix'; SparseMatrix.identity = SparseMatrix.eye; SparseMatrix.prototype.tensorProduct = SparseMatrix.prototype.kroneckerProduct; -/** - * Converts a matrix from Compressed Sparse Row (CSR) format to Compressed Sparse Column (CSC) format. - * @param {Object} csrMatrix - The matrix in CSR format with properties: values, columns, rows. - * @param {number} numCols - The number of columns in the matrix. - * @returns {{rows: Float64Array, columns: Float64Array, values: Float64Array}} The matrix in CSC format. - */ -function csrToCsc(csrMatrix, numCols) { - const { - values: csrValues, - columns: csrColIndices, - rows: csrRowPtr, - } = csrMatrix; - - const cscValues = new Float64Array(csrValues.length); - const cscRowIndices = new Float64Array(csrValues.length); - const cscColPtr = new Float64Array(numCols + 1); - - for (let i = 0; i < csrColIndices.length; i++) { - cscColPtr[csrColIndices[i] + 1]++; - } - - for (let i = 1; i <= numCols; i++) { - cscColPtr[i] += cscColPtr[i - 1]; - } - - const next = cscColPtr.slice(); - - for (let row = 0; row < csrRowPtr.length - 1; row++) { - for (let j = csrRowPtr[row], i = 0; j < csrRowPtr[row + 1]; j++, i++) { - const col = csrColIndices[j]; - const pos = next[col]; - cscValues[pos] = csrValues[j]; - cscRowIndices[pos] = row; - next[col]++; - } - } - - return { rows: cscRowIndices, columns: cscColPtr, values: cscValues }; -} - /** * Converts a matrix from Coordinate (COO) format to Compressed Sparse Row (CSR) format. * @param {Object} cooMatrix - The matrix in COO format with properties: values, columns, rows. @@ -1518,26 +1527,17 @@ function csrToCsc(csrMatrix, numCols) { function cooToCsr(cooMatrix, nbRows) { const { values, columns, rows } = cooMatrix; const csrRowPtr = new Float64Array(nbRows + 1); - const length = values.length; - let currentRow = rows[0]; - for (let index = 0; index < length; ) { - while (currentRow === rows[index] && index < length) ++index; - csrRowPtr[currentRow + 1] = index; - currentRow += 1; + + // Count non-zeros per row + const numberOfNonZeros = rows.length; + for (let i = 0; i < numberOfNonZeros; i++) { + csrRowPtr[rows[i] + 1]++; } - return { rows: csrRowPtr, columns, values }; -} -/** - * Creates an empty 2D matrix (array of Float64Array) with the given dimensions. - * @param {number} nbRows - Number of rows. - * @param {number} nbColumns - Number of columns. - * @returns {Float64Array[]} A 2D array representing the matrix. - */ -function matrixCreateEmpty(nbRows, nbColumns) { - const newMatrix = []; - for (let row = 0; row < nbRows; row++) { - newMatrix.push(new Float64Array(nbColumns)); + // Compute cumulative sum + for (let i = 1; i <= nbRows; i++) { + csrRowPtr[i] += csrRowPtr[i - 1]; } - return newMatrix; + + return { rows: csrRowPtr, columns, values }; } From 80686cd3df634a2fd8965192ccb6839f9fb14740 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Wed, 2 Jul 2025 05:55:21 -0500 Subject: [PATCH 25/33] chore: remove benchmark script --- benchmark/benchmarkMitata.js | 40 ------------------------------------ 1 file changed, 40 deletions(-) delete mode 100644 benchmark/benchmarkMitata.js diff --git a/benchmark/benchmarkMitata.js b/benchmark/benchmarkMitata.js deleted file mode 100644 index a99ade4..0000000 --- a/benchmark/benchmarkMitata.js +++ /dev/null @@ -1,40 +0,0 @@ -import { run, bench, group, do_not_optimize } from 'mitata'; -import { Matrix } from 'ml-matrix'; - -import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld2.js'; -import { randomSparseMatrix } from './utils/randomSparseMatrix.js'; - -/* eslint -func-names: 0 -camelcase: 0 -*/ - -const sizes = [64, 128, 256]; -const densities = [0.01, 0.015, 0.02, 0.025, 0.03]; - -for (const density of densities) { - for (const size of sizes) { - const A = randomSparseMatrix(size, size, density); - const B = randomSparseMatrix(size, size, density); - let denseA = A.to2DArray(); - let denseB = B.to2DArray(); - const AOld = new SparseMatrixOld(denseA); - const BOld = new SparseMatrixOld(denseB); - denseA = new Matrix(denseA); - denseB = new Matrix(denseB); - - group(`size:${size}-density:${density}`, () => { - bench('mmulNew', () => { - do_not_optimize(A.mmul(B)); - }); //.gc('inner'); - bench('mmul', () => { - do_not_optimize(AOld.mmul(BOld)); - }); //.gc('inner'); - bench('denseMatrix', () => { - do_not_optimize(denseA.mmul(denseB)); - }); //.gc('inner'); - }); - } -} - -await run({ silent: false }); From 255e35161f66985404140d3d8cd3e5eabaa76f5b Mon Sep 17 00:00:00 2001 From: jobo322 Date: Wed, 2 Jul 2025 08:00:21 -0500 Subject: [PATCH 26/33] chore: add benchmarks for cardinality and size comparisons --- benchmark/benchmarkLinePlot.js | 10 +- benchmark/combinationOfCardinalityAndSize.js | 142 +++++++++++++++ benchmark/densityDecreaseWhenSizeIncrease.js | 37 ++++ benchmark/plot.html | 174 +++++++++++++++++++ 4 files changed, 360 insertions(+), 3 deletions(-) create mode 100644 benchmark/combinationOfCardinalityAndSize.js create mode 100644 benchmark/densityDecreaseWhenSizeIncrease.js create mode 100644 benchmark/plot.html diff --git a/benchmark/benchmarkLinePlot.js b/benchmark/benchmarkLinePlot.js index b413c2d..452dca7 100644 --- a/benchmark/benchmarkLinePlot.js +++ b/benchmark/benchmarkLinePlot.js @@ -25,7 +25,9 @@ lineplot(() => { const B = new SparseMatrix(randomMatrix(size, size, density)); // Benchmark the multiplication yield () => do_not_optimize(A.mmul(B)); - }).args('size', sizes); // 16, 32, 64, 128, 256 + }) + .gc('inner') + .args('size', sizes); // 16, 32, 64, 128, 256 bench('SparseOld.mmul($size)', function* (ctx) { const size = ctx.get('size'); @@ -36,7 +38,9 @@ lineplot(() => { // Benchmark the multiplication yield () => do_not_optimize(AOld.mmul(BOld)); - }).args('size', sizes); + }) + .gc('inner') + .args('size', sizes); // bench('Dense.mmul($size)', function* (ctx) { // const size = ctx.get('size'); @@ -49,7 +53,7 @@ lineplot(() => { // // Benchmark the multiplication // yield () => do_not_optimize(ADense.mmul(BDense)); - // }).range('size', min, max, multiplier); + // }).gc('inner').range('size', min, max, multiplier); }); await run(); diff --git a/benchmark/combinationOfCardinalityAndSize.js b/benchmark/combinationOfCardinalityAndSize.js new file mode 100644 index 0000000..f15a727 --- /dev/null +++ b/benchmark/combinationOfCardinalityAndSize.js @@ -0,0 +1,142 @@ +import { run, bench, group, do_not_optimize } from 'mitata'; +import { xSequentialFillFromStep } from 'ml-spectra-processing'; + +import { SparseMatrix } from '../src/index.js'; + +import { randomMatrix } from './utils/randomMatrix.js'; +import { writeFile } from 'node:fs/promises'; +import path from 'node:path'; + +/* eslint +func-names: 0 +camelcase: 0 +*/ +// Prepare matrices once +const cardinalities = Array.from( + xSequentialFillFromStep({ from: 10, step: 5, size: 25 }), +); + +// const dimensions = Array.from( +// xSequentialFillFromStep({ from: 700, step: 100, size: 13 }), +// ); + +const dimensions = [512]; +group('comparation internal method', () => { + bench('hibrid($cardinality,$dimension)', function* (ctx) { + const cardinality = ctx.get('cardinality'); + const size = ctx.get('dimension'); + // Prepare matrices once + let A = new SparseMatrix(randomMatrix(size, size, cardinality)); + let B = new SparseMatrix(randomMatrix(size, size, cardinality)); + + // Benchmark the multiplication + yield () => do_not_optimize(A.mmul(B)); + do_not_optimize(A); + do_not_optimize(B); + }) + .gc('inner') + .args('cardinality', cardinalities) //.range('size', 32, 1024, 2); //.args('size', sizes); + .args('dimension', dimensions); + + bench('small($cardinality,$dimension)', function* (ctx) { + const cardinality = ctx.get('cardinality'); + const size = ctx.get('dimension'); + // Prepare matrices once + let A = new SparseMatrix(randomMatrix(size, size, cardinality)); + let B = new SparseMatrix(randomMatrix(size, size, cardinality)); + + // Benchmark the multiplication + yield () => do_not_optimize(A._mmulSmall(B)); + // Explicit cleanup + do_not_optimize(A); + do_not_optimize(B); + }) + .gc('inner') + .args('cardinality', cardinalities) //.range('size', 32, 1024, 2); //.args('size', sizes); + .args('dimension', dimensions); + + bench('low($cardinality,$dimension)', function* (ctx) { + const cardinality = ctx.get('cardinality'); + const size = ctx.get('dimension'); + // Prepare matrices once + let A = new SparseMatrix(randomMatrix(size, size, cardinality)); + let B = new SparseMatrix(randomMatrix(size, size, cardinality)); + + // Benchmark the multiplication + yield () => do_not_optimize(A._mmulLowDensity(B)); + // Explicit cleanup + do_not_optimize(A); + do_not_optimize(B); + }) + .gc('inner') + .args('cardinality', cardinalities) //.range('size', 32, 1024, 2); //.args('size', sizes); + .args('dimension', dimensions); + + bench('medium($cardinality,$dimension)', function* (ctx) { + const cardinality = ctx.get('cardinality'); + const size = ctx.get('dimension'); + // Prepare matrices once + let A = new SparseMatrix(randomMatrix(size, size, cardinality)); + let B = new SparseMatrix(randomMatrix(size, size, cardinality)); + + // Benchmark the multiplication + yield () => { + do_not_optimize(A._mmulMediumDensity(B)); + }; + + // Explicit cleanup + do_not_optimize(A); + do_not_optimize(B); + }) + .gc('inner') + .args('cardinality', cardinalities) //.range('size', 32, 1024, 2); //.args('size', sizes); + .args('dimension', dimensions); +}); + +// Run benchmarks and capture results +const results = await run({ + // Force GC between every benchmark + gc: true, + // // More samples for statistical significance + // min_samples: 20, + // max_samples: 200, + // // Longer warmup to stabilize CPU state + // warmup_samples: 10, + // warmup_threshold: 100, // ms + // Longer minimum time for stable measurements + min_cpu_time: 2000, // 2 seconds minimum + // Batch settings to reduce variance + batch_samples: 5, + batch_threshold: 10, // ms + // Enable colors + colors: true, +}); + +// Process and store results +const processedResults = []; + +for (const benchmark of results.benchmarks) { + for (const run of benchmark.runs) { + if (run.stats) { + processedResults.push({ + name: benchmark.alias, + cardinality: run.args.cardinality, + dimension: run.args.dimension, + avg: run.stats.avg, + min: run.stats.min, + max: run.stats.max, + p50: run.stats.p50, + p75: run.stats.p75, + p99: run.stats.p99, + samples: run.stats.samples.length, + ticks: run.stats.ticks, + }); + } + } +} + +// Save results to JSON file +await writeFile( + path.join(import.meta.dirname, `benchmark-results.json`), + JSON.stringify(processedResults, null, 2), +); diff --git a/benchmark/densityDecreaseWhenSizeIncrease.js b/benchmark/densityDecreaseWhenSizeIncrease.js new file mode 100644 index 0000000..7d828a2 --- /dev/null +++ b/benchmark/densityDecreaseWhenSizeIncrease.js @@ -0,0 +1,37 @@ +import { run, bench, group, do_not_optimize } from 'mitata'; +// import { Matrix } from 'ml-matrix'; + +import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; +import { randomMatrix } from './utils/randomMatrix.js'; +import { SparseMatrix } from '../src/index.js'; + +/* eslint +func-names: 0 +camelcase: 0 +*/ + +const sizes = [8, 16, 32, 256, 512, 1024]; +const densities = [0.125, 0.0625, 0.03125, 0.0039, 0.00197, 0.001]; + +for (let i = 0; i < sizes.length; i++) { + const size = sizes[i]; + const density = densities[i]; + const denseA = randomMatrix(size, size, density); + const denseB = randomMatrix(size, size, density); + const AOld = new SparseMatrixOld(denseA); + const BOld = new SparseMatrixOld(denseB); + + const A = new SparseMatrix(denseA); + const B = new SparseMatrix(denseB); + + group(`size:${size}-density:${density}`, () => { + bench('mmulNew', () => { + do_not_optimize(A.mmul(B)); + }).gc('inner'); + bench('mmul', () => { + do_not_optimize(AOld.mmul(BOld)); + }).gc('inner'); + }); +} + +await run({ silent: false }); diff --git a/benchmark/plot.html b/benchmark/plot.html new file mode 100644 index 0000000..3309757 --- /dev/null +++ b/benchmark/plot.html @@ -0,0 +1,174 @@ + + + + + Benchmark Results Plot + + + + +

Benchmark Results: avg vs cardinality

+ +
+
+ + + \ No newline at end of file From 87966e4eab84bf9835c78676c288e9bab667725f Mon Sep 17 00:00:00 2001 From: jobo322 Date: Wed, 2 Jul 2025 08:38:37 -0500 Subject: [PATCH 27/33] chore: fix randommatrix util --- benchmark/combinationOfCardinalityAndSize.js | 18 +++++++++--------- benchmark/utils/randomMatrix.js | 3 +-- src/index.js | 2 +- 3 files changed, 11 insertions(+), 12 deletions(-) diff --git a/benchmark/combinationOfCardinalityAndSize.js b/benchmark/combinationOfCardinalityAndSize.js index f15a727..0105161 100644 --- a/benchmark/combinationOfCardinalityAndSize.js +++ b/benchmark/combinationOfCardinalityAndSize.js @@ -1,4 +1,4 @@ -import { run, bench, group, do_not_optimize } from 'mitata'; +import { run, bench, do_not_optimize, lineplot } from 'mitata'; import { xSequentialFillFromStep } from 'ml-spectra-processing'; import { SparseMatrix } from '../src/index.js'; @@ -13,7 +13,7 @@ camelcase: 0 */ // Prepare matrices once const cardinalities = Array.from( - xSequentialFillFromStep({ from: 10, step: 5, size: 25 }), + xSequentialFillFromStep({ from: 10, step: 5, size: 2 }), ); // const dimensions = Array.from( @@ -21,7 +21,7 @@ const cardinalities = Array.from( // ); const dimensions = [512]; -group('comparation internal method', () => { +lineplot(() => { bench('hibrid($cardinality,$dimension)', function* (ctx) { const cardinality = ctx.get('cardinality'); const size = ctx.get('dimension'); @@ -104,12 +104,12 @@ const results = await run({ // warmup_samples: 10, // warmup_threshold: 100, // ms // Longer minimum time for stable measurements - min_cpu_time: 2000, // 2 seconds minimum - // Batch settings to reduce variance - batch_samples: 5, - batch_threshold: 10, // ms - // Enable colors - colors: true, + // min_cpu_time: 2000, // 2 seconds minimum + // // Batch settings to reduce variance + // batch_samples: 5, + // batch_threshold: 10, // ms + // // Enable colors + // colors: true, }); // Process and store results diff --git a/benchmark/utils/randomMatrix.js b/benchmark/utils/randomMatrix.js index 3db3b66..41f4923 100644 --- a/benchmark/utils/randomMatrix.js +++ b/benchmark/utils/randomMatrix.js @@ -1,6 +1,5 @@ -export function randomMatrix(rows, cols, density = 0.01) { +export function randomMatrix(rows, cols, cardinality) { const total = rows * cols; - const cardinality = Math.ceil(total * density); const positions = new Set(); // Generate unique random positions diff --git a/src/index.js b/src/index.js index 351e3e3..1404021 100644 --- a/src/index.js +++ b/src/index.js @@ -165,7 +165,7 @@ export class SparseMatrix { if (this.cardinality < 42 && other.cardinality < 42) { return this._mmulSmall(other); - } else if (other.rows > 100 && other.cardinality < 110) { + } else if (other.rows > 100 && other.cardinality < 100) { return this._mmulLowDensity(other); } From 163c079536dbe593a6ee94020f968aff42f71362 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Micha=C3=ABl=20Zasso?= Date: Thu, 3 Jul 2025 17:09:19 +0200 Subject: [PATCH 28/33] chore: fix lint and remove redundant vitest config --- .github/workflows/typedoc.yml | 33 +- benchmark/benchmarkLinePlot.js | 6 +- benchmark/combinationOfCardinalityAndSize.js | 11 +- benchmark/densityDecreaseWhenSizeIncrease.js | 10 +- benchmark/plot.html | 350 ++++++++++--------- eslint.config.mjs => eslint.config.js | 2 + vitest.config.mjs | 14 - 7 files changed, 209 insertions(+), 217 deletions(-) rename eslint.config.mjs => eslint.config.js (86%) delete mode 100644 vitest.config.mjs diff --git a/.github/workflows/typedoc.yml b/.github/workflows/typedoc.yml index ff91151..f7816a7 100644 --- a/.github/workflows/typedoc.yml +++ b/.github/workflows/typedoc.yml @@ -1,32 +1,15 @@ -name: Deploy TypeDoc on GitHub pages +name: TypeDoc on: workflow_dispatch: release: types: [published] -env: - NODE_VERSION: 22.x - ENTRY_FILE: 'src/index.js' - jobs: - deploy: - runs-on: ubuntu-latest - steps: - - uses: actions/checkout@v4 - - uses: actions/setup-node@v4 - with: - node-version: ${{ env.NODE_VERSION }} - - name: Install dependencies - run: npm install - - name: Build documentation - uses: zakodium/typedoc-action@v2 - with: - entry: ${{ env.ENTRY_FILE }} - - name: Deploy to GitHub pages - uses: JamesIves/github-pages-deploy-action@releases/v4 - with: - token: ${{ secrets.BOT_TOKEN }} - branch: gh-pages - folder: docs - clean: true + typedoc: + # Documentation: https://github.com/zakodium/workflows#typedoc + uses: zakodium/workflows/.github/workflows/typedoc.yml@typedoc-v1 + with: + entry: 'src/index.js' + secrets: + github-token: ${{ secrets.BOT_TOKEN }} diff --git a/benchmark/benchmarkLinePlot.js b/benchmark/benchmarkLinePlot.js index 452dca7..175b042 100644 --- a/benchmark/benchmarkLinePlot.js +++ b/benchmark/benchmarkLinePlot.js @@ -1,4 +1,4 @@ -import { run, bench, lineplot, do_not_optimize } from 'mitata'; +import { bench, do_not_optimize, lineplot, run } from 'mitata'; import { xSequentialFillFromStep } from 'ml-spectra-processing'; import { SparseMatrix } from '../src/index.js'; @@ -8,10 +8,6 @@ import { randomMatrix } from './utils/randomMatrix.js'; const density = 0.02; // Fixed density for this comparison; -/* eslint -func-names: 0 -camelcase: 0 -*/ // Prepare matrices once const sizes = Array.from( xSequentialFillFromStep({ from: 4, step: 4, size: 13 }), diff --git a/benchmark/combinationOfCardinalityAndSize.js b/benchmark/combinationOfCardinalityAndSize.js index 0105161..2f6545d 100644 --- a/benchmark/combinationOfCardinalityAndSize.js +++ b/benchmark/combinationOfCardinalityAndSize.js @@ -1,16 +1,13 @@ -import { run, bench, do_not_optimize, lineplot } from 'mitata'; +import { writeFile } from 'node:fs/promises'; +import path from 'node:path'; + +import { bench, do_not_optimize, lineplot, run } from 'mitata'; import { xSequentialFillFromStep } from 'ml-spectra-processing'; import { SparseMatrix } from '../src/index.js'; import { randomMatrix } from './utils/randomMatrix.js'; -import { writeFile } from 'node:fs/promises'; -import path from 'node:path'; -/* eslint -func-names: 0 -camelcase: 0 -*/ // Prepare matrices once const cardinalities = Array.from( xSequentialFillFromStep({ from: 10, step: 5, size: 2 }), diff --git a/benchmark/densityDecreaseWhenSizeIncrease.js b/benchmark/densityDecreaseWhenSizeIncrease.js index 7d828a2..8ed6101 100644 --- a/benchmark/densityDecreaseWhenSizeIncrease.js +++ b/benchmark/densityDecreaseWhenSizeIncrease.js @@ -1,14 +1,10 @@ -import { run, bench, group, do_not_optimize } from 'mitata'; +import { bench, do_not_optimize, group, run } from 'mitata'; + // import { Matrix } from 'ml-matrix'; +import { SparseMatrix } from '../src/index.js'; import { SparseMatrix as SparseMatrixOld } from './class/SparseMatrixOld.js'; import { randomMatrix } from './utils/randomMatrix.js'; -import { SparseMatrix } from '../src/index.js'; - -/* eslint -func-names: 0 -camelcase: 0 -*/ const sizes = [8, 16, 32, 256, 512, 1024]; const densities = [0.125, 0.0625, 0.03125, 0.0039, 0.00197, 0.001]; diff --git a/benchmark/plot.html b/benchmark/plot.html index 3309757..4710a4d 100644 --- a/benchmark/plot.html +++ b/benchmark/plot.html @@ -1,174 +1,206 @@ - + - - - Benchmark Results Plot - - - - -

Benchmark Results: avg vs cardinality

- -
-
- + + + +

Benchmark Results: avg vs cardinality

+ +
+
+ - - \ No newline at end of file + async function main() { + let data = await loadData(); + if (data) { + renderCharts(data); + } else { + // Show file input for manual loading + document.getElementById('fileInput').style.display = ''; + document.getElementById('error').textContent = + "Could not load 'benchmark-results.json'. Please select the file manually."; + document + .getElementById('fileInput') + .addEventListener('change', function (e) { + const file = e.target.files[0]; + if (!file) return; + const reader = new FileReader(); + reader.onload = function (evt) { + try { + const json = JSON.parse(evt.target.result); + document.getElementById('error').textContent = ''; + renderCharts(json); + } catch (err) { + document.getElementById('error').textContent = + 'Invalid JSON file.'; + } + }; + reader.readAsText(file); + }); + } + } + + main(); + + + diff --git a/eslint.config.mjs b/eslint.config.js similarity index 86% rename from eslint.config.mjs rename to eslint.config.js index 79d8fb8..3ca4395 100644 --- a/eslint.config.mjs +++ b/eslint.config.js @@ -4,6 +4,8 @@ import ts from 'eslint-config-cheminfo-typescript/base'; export default defineConfig(globalIgnores(['benchmark/class/*', 'lib']), ts, { files: ['benchmark/**'], rules: { + camelcase: 'off', + 'func-names': 'off', 'no-console': 'off', 'no-unused-vars': 'off', }, diff --git a/vitest.config.mjs b/vitest.config.mjs deleted file mode 100644 index 0bf34cd..0000000 --- a/vitest.config.mjs +++ /dev/null @@ -1,14 +0,0 @@ -import { defineConfig, configDefaults } from 'vitest/config'; - -export default defineConfig({ - test: { - coverage: { - exclude: [ - ...configDefaults.coverage.exclude, - 'benchmark', - 'lib', - 'scripts', - ], - }, - }, -}); From 85054017617ae8c3cb5b07da7400cce333d897a2 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Thu, 3 Jul 2025 15:35:44 -0500 Subject: [PATCH 29/33] chore(getNonZeros): rename format options to csr and move cooToCsr as an utility function --- src/__tests__/getNonZeros.test.js | 2 +- src/index.js | 92 +++++++++++++++++-------------- src/utils/cooToCsr.ts | 23 ++++++++ 3 files changed, 74 insertions(+), 43 deletions(-) create mode 100644 src/utils/cooToCsr.ts diff --git a/src/__tests__/getNonZeros.test.js b/src/__tests__/getNonZeros.test.js index 6fc3182..0520c30 100644 --- a/src/__tests__/getNonZeros.test.js +++ b/src/__tests__/getNonZeros.test.js @@ -23,7 +23,7 @@ describe('Sparse Matrix', () => { }); // CSR format - expect(m2.getNonZeros({ format: true })).toEqual({ + expect(m2.getNonZeros({ csr: true })).toEqual({ rows: Float64Array.from([0, 0, 4, 7, 7, 11]), columns: Float64Array.from([0, 3, 4, 5, 1, 4, 5, 0, 3, 4, 5]), values: Float64Array.from([1, 2, 1, 1, 3, 5, 5, 1, 1, 9, 9]), diff --git a/src/index.js b/src/index.js index 5148384..ee24f5a 100644 --- a/src/index.js +++ b/src/index.js @@ -1,4 +1,5 @@ import HashTable from 'ml-hash-table'; +import { cooToCsr } from './utils/cooToCsr.js'; /** @typedef {(row: number, column: number, value: number) => void} WithEachNonZeroCallback */ /** @typedef {(row: number, column: number, value: number) => number | false} ForEachNonZeroCallback */ @@ -157,13 +158,14 @@ export class SparseMatrix { } /** + * Matrix multiplication, does not modify the current instance. * @param {SparseMatrix} other - * @returns {SparseMatrix} + * @returns {SparseMatrix} returns a new matrix instance. + * @throws {Error} If the number of columns of this matrix does not match the number of rows of the other matrix. */ mmul(other) { if (this.columns !== other.rows) { - // eslint-disable-next-line no-console - console.warn( + throw new RangeError( 'Number of columns of left matrix are not equal to number of rows of right matrix.', ); } @@ -177,6 +179,13 @@ export class SparseMatrix { return this._mmulMediumDensity(other); } + /** + * Matrix multiplication optimized for very small matrices (both cardinalities < 42). + * + * @private + * @param {SparseMatrix} other - The right-hand side matrix to multiply with. + * @returns {SparseMatrix} - The resulting matrix after multiplication. + */ _mmulSmall(other) { const m = this.rows; const p = other.columns; @@ -199,6 +208,13 @@ export class SparseMatrix { return result; } + /** + * Matrix multiplication optimized for low-density right-hand side matrices (other.rows > 100 and other.cardinality < 100). + * + * @private + * @param {SparseMatrix} other - The right-hand side matrix to multiply with. + * @returns {SparseMatrix} - The resulting matrix after multiplication. + */ _mmulLowDensity(other) { const m = this.rows; const p = other.columns; @@ -230,20 +246,27 @@ export class SparseMatrix { return result; } + /** + * Matrix multiplication for medium-density matrices using CSR format for the right-hand side. + * + * @private + * @param {SparseMatrix} other - The right-hand side matrix to multiply with. + * @returns {SparseMatrix} - The resulting matrix after multiplication. + */ _mmulMediumDensity(other) { - const m = this.rows; - const p = other.columns; - const { - columns: otherCols, - rows: otherRows, - values: otherValues, - } = other.getNonZeros({ format: 'csr' }); const { columns: thisCols, rows: thisRows, values: thisValues, } = this.getNonZeros(); + const { + columns: otherCols, + rows: otherRows, + values: otherValues, + } = other.getNonZeros({ csr: true }); + const m = this.rows; + const p = other.columns; const result = new SparseMatrix(m, p); const nbThisActive = thisCols.length; for (let t = 0; t < nbThisActive; t++) { @@ -358,14 +381,23 @@ export class SparseMatrix { } /** - * Returns the non-zero elements of the matrix in coordinate (COO), CSR, or CSC format. + * Returns the non-zero elements of the matrix in coordinates (COO) or CSR format. + * + * **COO (Coordinate) format:** + * Stores the non-zero elements as three arrays: `rows`, `columns`, and `values`, where each index corresponds to a non-zero entry at (row, column) with the given value. + * + * **CSR (Compressed Sparse Row) format:** + * Stores the matrix using three arrays: + * - `rows`: Row pointer array of length `numRows + 1`, where each entry indicates the start of a row in the `columns` and `values` arrays. + * - `columns`: Column indices of non-zero elements. + * - `values`: Non-zero values. + * This format is efficient for row slicing and matrix-vector multiplication. * * @param {Object} [options={}] - Options for output format and sorting. - * @param {boolean} [options.format] - If specified, returns the result in CSR or CSC format. Otherwise, returns COO format. + * @param {boolean} [options.csr] - If true, returns the result in CSR format. Otherwise, returns COO format. * @param {boolean} [options.sort] - If true, sorts the non-zero elements by their indices. - * @returns {Object} If no format is specified, returns an object with Float64Array `rows`, `columns`, and `values` (COO format). - * If format is true, returns { rows, columns, values } in CSR format. - * @throws {Error} If an unsupported format is specified. + * @returns {Object} If `csr` is not specified, returns an object with Float64Array `rows`, `columns`, and `values` (COO format). + * If `csr` is true, returns `{ rows, columns, values }` in CSR format. */ getNonZeros(options = {}) { const cardinality = this.cardinality; @@ -373,7 +405,7 @@ export class SparseMatrix { const columns = new Float64Array(cardinality); const values = new Float64Array(cardinality); - const { format, sort = format !== undefined } = options; + const { csr, sort } = options; let idx = 0; this.withEachNonZero((i, j, value) => { @@ -381,9 +413,9 @@ export class SparseMatrix { columns[idx] = j; values[idx] = value; idx++; - }, sort); + }, sort || csr); - return format + return csr ? cooToCsr({ rows, columns, values }, this.rows) : { rows, columns, values }; } @@ -1518,27 +1550,3 @@ SparseMatrix.prototype.klass = 'Matrix'; SparseMatrix.identity = SparseMatrix.eye; SparseMatrix.prototype.tensorProduct = SparseMatrix.prototype.kroneckerProduct; - -/** - * Converts a matrix from Coordinate (COO) format to Compressed Sparse Row (CSR) format. - * @param {Object} cooMatrix - The matrix in COO format with properties: values, columns, rows. - * @param {number} nbRows - The number of rows in the matrix. - * @returns {{rows: Float64Array, columns: Float64Array, values: Float64Array}} The matrix in CSR format. - */ -function cooToCsr(cooMatrix, nbRows) { - const { values, columns, rows } = cooMatrix; - const csrRowPtr = new Float64Array(nbRows + 1); - - // Count non-zeros per row - const numberOfNonZeros = rows.length; - for (let i = 0; i < numberOfNonZeros; i++) { - csrRowPtr[rows[i] + 1]++; - } - - // Compute cumulative sum - for (let i = 1; i <= nbRows; i++) { - csrRowPtr[i] += csrRowPtr[i - 1]; - } - - return { rows: csrRowPtr, columns, values }; -} diff --git a/src/utils/cooToCsr.ts b/src/utils/cooToCsr.ts new file mode 100644 index 0000000..0f4a2dd --- /dev/null +++ b/src/utils/cooToCsr.ts @@ -0,0 +1,23 @@ +/** + * Converts a matrix from Coordinate (COO) format to Compressed Sparse Row (CSR) format. + * @param {Object} cooMatrix - The matrix in COO format with properties: values, columns, rows. + * @param {number} nbRows - The number of rows in the matrix. + * @returns {{rows: Float64Array, columns: Float64Array, values: Float64Array}} The matrix in CSR format. + */ +export function cooToCsr(cooMatrix, nbRows) { + const { values, columns, rows } = cooMatrix; + const csrRowPtr = new Float64Array(nbRows + 1); + + // Count non-zeros per row + const numberOfNonZeros = rows.length; + for (let i = 0; i < numberOfNonZeros; i++) { + csrRowPtr[rows[i] + 1]++; + } + + // Compute cumulative sum + for (let i = 1; i <= nbRows; i++) { + csrRowPtr[i] += csrRowPtr[i - 1]; + } + + return { rows: csrRowPtr, columns, values }; +} From 10aee91f42fa016330c427eb805f9ec0406a4dc7 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Thu, 3 Jul 2025 16:13:56 -0500 Subject: [PATCH 30/33] test: refactor matrix multiplication tests and improve random matrix generation --- src/__tests__/index.test.js | 96 ++++++++++++++++---------- src/index.js | 3 +- src/utils/{cooToCsr.ts => cooToCsr.js} | 0 3 files changed, 59 insertions(+), 40 deletions(-) rename src/utils/{cooToCsr.ts => cooToCsr.js} (100%) diff --git a/src/__tests__/index.test.js b/src/__tests__/index.test.js index 1a6b38a..81ad9b3 100644 --- a/src/__tests__/index.test.js +++ b/src/__tests__/index.test.js @@ -1,3 +1,4 @@ +import { Matrix } from 'ml-matrix'; import { describe, expect, it } from 'vitest'; import { SparseMatrix } from '../index.js'; @@ -38,33 +39,46 @@ describe('Sparse Matrix', () => { expect(m3.cardinality).toBe(1); expect(m3.get(0, 1)).toBe(2); - expect(m3.to2DArray()).toStrictEqual([ + expectMatrixClose(m3.to2DArray(), [ [0, 2], [0, 0], ]); // Compare with dense multiplication - const denseM1 = m1.to2DArray(); - const denseM2 = m2.to2DArray(); - const expectedDense = denseMatrixMultiply(denseM1, denseM2); - expect(m3.to2DArray()).toStrictEqual(expectedDense); + const denseM1 = new Matrix(m1.to2DArray()); + const denseM2 = new Matrix(m2.to2DArray()); + const expectedDense = denseM1.mmul(denseM2); + expectMatrixClose(m3.to2DArray(), expectedDense.to2DArray()); }); it('mmul', () => { const size = 32; const density = 0.1; - const m1 = randomSparseMatrix(size, size, density); - const m2 = randomSparseMatrix(size, size, density); - let m3 = m1.mmul(m2); - - const denseM1 = m1.to2DArray(); - const denseM2 = m2.to2DArray(); - - const newSparse = new SparseMatrix(denseM1); - expect(newSparse.to2DArray()).toStrictEqual(denseM1); - const expectedDense = denseMatrixMultiply(denseM1, denseM2); + const A = randomMatrix(size, size, density * size ** 2); + const B = randomMatrix(size, size, density * size ** 2); + const m1 = new SparseMatrix(A); + const m2 = new SparseMatrix(B); + const m3 = m1.mmul(m2); + + const denseM1 = new Matrix(A); + const denseM2 = new Matrix(B); + const expectedDense = denseM1.mmul(denseM2); + expectMatrixClose(m3.to2DArray(), expectedDense.to2DArray()); + }); - expect(m3.to2DArray()).toStrictEqual(expectedDense); + it('mmul with low density', () => { + const size = 128; + const cardinality = 64; + const A = randomMatrix(size, size, cardinality); + const B = randomMatrix(size, size, cardinality); + const m1 = new SparseMatrix(A); + const m2 = new SparseMatrix(B); + const m3 = m1.mmul(m2); + + const denseM1 = new Matrix(A); + const denseM2 = new Matrix(B); + const expectedDense = denseM1.mmul(denseM2); + expectMatrixClose(m3.to2DArray(), expectedDense.to2DArray()); }); it('kronecker', () => { @@ -168,31 +182,37 @@ describe('Banded matrices', () => { }); }); -function denseMatrixMultiply(A, B) { - const rowsA = A.length; - const colsA = A[0].length; - const colsB = B[0].length; - const result = Array.from({ length: rowsA }, () => Array(colsB).fill(0)); - for (let i = 0; i < rowsA; i++) { - for (let j = 0; j < colsB; j++) { - for (let k = 0; k < colsA; k++) { - result[i][j] += A[i][k] * B[k][j]; - } +/** + * Helper to compare two 2D arrays element-wise using toBeCloseTo. + */ +function expectMatrixClose(received, expected, precision = 6) { + expect(received.length).toBe(expected.length); + for (let i = 0; i < received.length; i++) { + expect(received[i].length).toBe(expected[i].length); + for (let j = 0; j < received[i].length; j++) { + expect(received[i][j]).toBeCloseTo(expected[i][j], precision); } } - return result; } -function randomSparseMatrix(rows, cols, density = 0.01) { - const matrix = []; - for (let i = 0; i < rows; i++) { - const row = new Float64Array(cols); - for (let j = 0; j < cols; j++) { - if (Math.random() < density) { - row[j] = Math.random() * 10; - } - } - matrix.push(row); +function randomMatrix(rows, cols, cardinality) { + const total = rows * cols; + const positions = new Set(); + + // Generate unique random positions + while (positions.size < cardinality) { + positions.add(Math.floor(Math.random() * total)); } - return new SparseMatrix(matrix); + + // Build the matrix with zeros + const matrix = Array.from({ length: rows }, () => new Float64Array(cols)); + + // Assign random values to the selected positions + for (const pos of positions) { + const i = Math.floor(pos / cols); + const j = pos % cols; + matrix[i][j] = Math.random() * 10; + } + + return matrix; } diff --git a/src/index.js b/src/index.js index ee24f5a..5c23318 100644 --- a/src/index.js +++ b/src/index.js @@ -1,4 +1,5 @@ import HashTable from 'ml-hash-table'; + import { cooToCsr } from './utils/cooToCsr.js'; /** @typedef {(row: number, column: number, value: number) => void} WithEachNonZeroCallback */ @@ -169,13 +170,11 @@ export class SparseMatrix { 'Number of columns of left matrix are not equal to number of rows of right matrix.', ); } - if (this.cardinality < 42 && other.cardinality < 42) { return this._mmulSmall(other); } else if (other.rows > 100 && other.cardinality < 100) { return this._mmulLowDensity(other); } - return this._mmulMediumDensity(other); } diff --git a/src/utils/cooToCsr.ts b/src/utils/cooToCsr.js similarity index 100% rename from src/utils/cooToCsr.ts rename to src/utils/cooToCsr.js From e287ca19d5decce5c7d1d20c4fb009e81f427d36 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Thu, 3 Jul 2025 16:30:41 -0500 Subject: [PATCH 31/33] feat: move internal specialized matrix multiplication as pure functions --- src/index.js | 112 +++------------------------------ src/utils/mmulLowDensity.js | 41 ++++++++++++ src/utils/mmulMediumDensity.js | 41 ++++++++++++ src/utils/mmulSmall.js | 33 ++++++++++ 4 files changed, 122 insertions(+), 105 deletions(-) create mode 100644 src/utils/mmulLowDensity.js create mode 100644 src/utils/mmulMediumDensity.js create mode 100644 src/utils/mmulSmall.js diff --git a/src/index.js b/src/index.js index 5c23318..6c23c82 100644 --- a/src/index.js +++ b/src/index.js @@ -1,6 +1,9 @@ import HashTable from 'ml-hash-table'; import { cooToCsr } from './utils/cooToCsr.js'; +import { mmulLowDensity } from './utils/mmulLowDensity.js'; +import { mmulMediumDensity } from './utils/mmulMediumDensity.js'; +import { mmulSmall } from './utils/mmulSmall.js'; /** @typedef {(row: number, column: number, value: number) => void} WithEachNonZeroCallback */ /** @typedef {(row: number, column: number, value: number) => number | false} ForEachNonZeroCallback */ @@ -170,117 +173,16 @@ export class SparseMatrix { 'Number of columns of left matrix are not equal to number of rows of right matrix.', ); } + if (this.cardinality < 42 && other.cardinality < 42) { - return this._mmulSmall(other); + return mmulSmall(this, other); } else if (other.rows > 100 && other.cardinality < 100) { - return this._mmulLowDensity(other); + return mmulLowDensity(this, other); } - return this._mmulMediumDensity(other); - } - - /** - * Matrix multiplication optimized for very small matrices (both cardinalities < 42). - * - * @private - * @param {SparseMatrix} other - The right-hand side matrix to multiply with. - * @returns {SparseMatrix} - The resulting matrix after multiplication. - */ - _mmulSmall(other) { - const m = this.rows; - const p = other.columns; - const { - columns: otherCols, - rows: otherRows, - values: otherValues, - } = other.getNonZeros(); - - const nbOtherActive = otherCols.length; - const result = new SparseMatrix(m, p); - this.withEachNonZero((i, j, v1) => { - for (let o = 0; o < nbOtherActive; o++) { - if (j === otherRows[o]) { - const l = otherCols[o]; - result.set(i, l, result.get(i, l) + otherValues[o] * v1); - } - } - }); - return result; - } - - /** - * Matrix multiplication optimized for low-density right-hand side matrices (other.rows > 100 and other.cardinality < 100). - * - * @private - * @param {SparseMatrix} other - The right-hand side matrix to multiply with. - * @returns {SparseMatrix} - The resulting matrix after multiplication. - */ - _mmulLowDensity(other) { - const m = this.rows; - const p = other.columns; - const { - columns: otherCols, - rows: otherRows, - values: otherValues, - } = other.getNonZeros(); - const { - columns: thisCols, - rows: thisRows, - values: thisValues, - } = this.getNonZeros(); - const result = new SparseMatrix(m, p); - const nbOtherActive = otherCols.length; - const nbThisActive = thisCols.length; - for (let t = 0; t < nbThisActive; t++) { - const i = thisRows[t]; - const j = thisCols[t]; - for (let o = 0; o < nbOtherActive; o++) { - if (j === otherRows[o]) { - const l = otherCols[o]; - result.set(i, l, result.get(i, l) + otherValues[o] * thisValues[t]); - } - } - } - // console.log(result.cardinality); - return result; + return mmulMediumDensity(this, other); } - /** - * Matrix multiplication for medium-density matrices using CSR format for the right-hand side. - * - * @private - * @param {SparseMatrix} other - The right-hand side matrix to multiply with. - * @returns {SparseMatrix} - The resulting matrix after multiplication. - */ - _mmulMediumDensity(other) { - const { - columns: thisCols, - rows: thisRows, - values: thisValues, - } = this.getNonZeros(); - const { - columns: otherCols, - rows: otherRows, - values: otherValues, - } = other.getNonZeros({ csr: true }); - - const m = this.rows; - const p = other.columns; - const result = new SparseMatrix(m, p); - const nbThisActive = thisCols.length; - for (let t = 0; t < nbThisActive; t++) { - const i = thisRows[t]; - const j = thisCols[t]; - const oStart = otherRows[j]; - const oEnd = otherRows[j + 1]; - for (let k = oStart; k < oEnd; k++) { - const l = otherCols[k]; - result.set(i, l, result.get(i, l) + otherValues[k] * thisValues[t]); - } - } - - return result; - } /** * @param {SparseMatrix} other * @returns {SparseMatrix} diff --git a/src/utils/mmulLowDensity.js b/src/utils/mmulLowDensity.js new file mode 100644 index 0000000..0310357 --- /dev/null +++ b/src/utils/mmulLowDensity.js @@ -0,0 +1,41 @@ +import { SparseMatrix } from '../index.js'; + +/** + * Multiplies two sparse matrices, optimized for cases where the right-hand side matrix + * has low density (number of rows > 100 and cardinality < 100). + * + * @private + * @param {SparseMatrix} left - The left-hand side matrix. + * @param {SparseMatrix} right - The right-hand side matrix to multiply with. + * @returns {SparseMatrix} The resulting matrix after multiplication. + */ +export function mmulLowDensity(left, right) { + const { + columns: otherCols, + rows: otherRows, + values: otherValues, + } = right.getNonZeros(); + const { + columns: thisCols, + rows: thisRows, + values: thisValues, + } = left.getNonZeros(); + + const m = left.rows; + const p = right.columns; + const output = new SparseMatrix(m, p); + + const nbOtherActive = otherCols.length; + const nbThisActive = thisCols.length; + for (let t = 0; t < nbThisActive; t++) { + const i = thisRows[t]; + const j = thisCols[t]; + for (let o = 0; o < nbOtherActive; o++) { + if (j === otherRows[o]) { + const l = otherCols[o]; + output.set(i, l, output.get(i, l) + otherValues[o] * thisValues[t]); + } + } + } + return output; +} diff --git a/src/utils/mmulMediumDensity.js b/src/utils/mmulMediumDensity.js new file mode 100644 index 0000000..2bee6f8 --- /dev/null +++ b/src/utils/mmulMediumDensity.js @@ -0,0 +1,41 @@ +import { SparseMatrix } from '../index.js'; + +/** + * Multiplies two sparse matrices where the right-hand side matrix is of medium density. + * Uses CSR format for the right-hand side for efficient multiplication. + * + * @private + * @param {SparseMatrix} left - The left-hand side matrix. + * @param {SparseMatrix} right - The right-hand side matrix to multiply with (in CSR format). + * @returns {SparseMatrix} The resulting matrix after multiplication. + */ + +export function mmulMediumDensity(left, right) { + const m = left.rows; + const p = right.columns; + const result = new SparseMatrix(m, p); + const { + columns: thisCols, + rows: thisRows, + values: thisValues, + } = left.getNonZeros(); + const { + columns: otherCols, + rows: otherRows, + values: otherValues, + } = right.getNonZeros({ csr: true }); + + const nbThisActive = thisCols.length; + for (let t = 0; t < nbThisActive; t++) { + const i = thisRows[t]; + const j = thisCols[t]; + const oStart = otherRows[j]; + const oEnd = otherRows[j + 1]; + for (let k = oStart; k < oEnd; k++) { + const l = otherCols[k]; + result.set(i, l, result.get(i, l) + otherValues[k] * thisValues[t]); + } + } + + return result; +} diff --git a/src/utils/mmulSmall.js b/src/utils/mmulSmall.js new file mode 100644 index 0000000..ba23481 --- /dev/null +++ b/src/utils/mmulSmall.js @@ -0,0 +1,33 @@ +import { SparseMatrix } from '../index.js'; + +/** + * Multiplies two very small sparse matrices (both with cardinalities < 42). + * + * @private + * @param {SparseMatrix} left - The left-hand side matrix. + * @param {SparseMatrix} right - The right-hand side matrix to multiply with. + * @returns {SparseMatrix} The resulting matrix after multiplication. + */ +export function mmulSmall(left, right) { + const { + columns: otherCols, + rows: otherRows, + values: otherValues, + } = right.getNonZeros(); + + const nbOtherActive = otherCols.length; + + const m = left.rows; + const p = right.columns; + const output = new SparseMatrix(m, p); + + left.withEachNonZero((i, j, v1) => { + for (let o = 0; o < nbOtherActive; o++) { + if (j === otherRows[o]) { + const l = otherCols[o]; + output.set(i, l, output.get(i, l) + otherValues[o] * v1); + } + } + }); + return output; +} From f101f9df5c10dfa71754f0245428b5c3abd6a0e6 Mon Sep 17 00:00:00 2001 From: jobo322 Date: Thu, 3 Jul 2025 16:35:05 -0500 Subject: [PATCH 32/33] chore: update benchmark --- benchmark/combinationOfCardinalityAndSize.js | 73 +++++++++++++------- 1 file changed, 49 insertions(+), 24 deletions(-) diff --git a/benchmark/combinationOfCardinalityAndSize.js b/benchmark/combinationOfCardinalityAndSize.js index 2f6545d..d9217fa 100644 --- a/benchmark/combinationOfCardinalityAndSize.js +++ b/benchmark/combinationOfCardinalityAndSize.js @@ -1,23 +1,31 @@ import { writeFile } from 'node:fs/promises'; import path from 'node:path'; -import { bench, do_not_optimize, lineplot, run } from 'mitata'; +import { run, bench, do_not_optimize, lineplot } from 'mitata'; +import { Matrix } from 'ml-matrix'; import { xSequentialFillFromStep } from 'ml-spectra-processing'; import { SparseMatrix } from '../src/index.js'; import { randomMatrix } from './utils/randomMatrix.js'; - +import { mmulMediumDensity } from '../src/utils/mmulMediumDensity.js'; +import { mmulSmall } from '../src/utils/mmulSmall.js'; +import { mmulLowDensity } from '../src/utils/mmulLowDensity.js'; + +/* eslint +func-names: 0 +camelcase: 0 +*/ // Prepare matrices once const cardinalities = Array.from( - xSequentialFillFromStep({ from: 10, step: 5, size: 2 }), + xSequentialFillFromStep({ from: 10, step: 50, size: 9 }), ); // const dimensions = Array.from( -// xSequentialFillFromStep({ from: 700, step: 100, size: 13 }), +// xSequentialFillFromStep({ from: 30, step: 100, size: 20 }), // ); -const dimensions = [512]; +const dimensions = [32]; lineplot(() => { bench('hibrid($cardinality,$dimension)', function* (ctx) { const cardinality = ctx.get('cardinality'); @@ -25,7 +33,7 @@ lineplot(() => { // Prepare matrices once let A = new SparseMatrix(randomMatrix(size, size, cardinality)); let B = new SparseMatrix(randomMatrix(size, size, cardinality)); - + A.mmul(B); // Benchmark the multiplication yield () => do_not_optimize(A.mmul(B)); do_not_optimize(A); @@ -41,9 +49,9 @@ lineplot(() => { // Prepare matrices once let A = new SparseMatrix(randomMatrix(size, size, cardinality)); let B = new SparseMatrix(randomMatrix(size, size, cardinality)); - + mmulSmall(A, B); // Benchmark the multiplication - yield () => do_not_optimize(A._mmulSmall(B)); + yield () => do_not_optimize(mmulSmall(A, B)); // Explicit cleanup do_not_optimize(A); do_not_optimize(B); @@ -58,9 +66,26 @@ lineplot(() => { // Prepare matrices once let A = new SparseMatrix(randomMatrix(size, size, cardinality)); let B = new SparseMatrix(randomMatrix(size, size, cardinality)); + mmulLowDensity(A, B); + // Benchmark the multiplication + yield () => do_not_optimize(mmulLowDensity(A, B)); + // Explicit cleanup + do_not_optimize(A); + do_not_optimize(B); + }) + .gc('inner') + .args('cardinality', cardinalities) //.range('size', 32, 1024, 2); //.args('size', sizes); + .args('dimension', dimensions); + bench('dense($cardinality,$dimension)', function* (ctx) { + const cardinality = ctx.get('cardinality'); + const size = ctx.get('dimension'); + // Prepare matrices once + let A = new Matrix(randomMatrix(size, size, cardinality)); + let B = new Matrix(randomMatrix(size, size, cardinality)); + A.mmul(B); // Benchmark the multiplication - yield () => do_not_optimize(A._mmulLowDensity(B)); + yield () => do_not_optimize(A.mmul(B)); // Explicit cleanup do_not_optimize(A); do_not_optimize(B); @@ -75,10 +100,10 @@ lineplot(() => { // Prepare matrices once let A = new SparseMatrix(randomMatrix(size, size, cardinality)); let B = new SparseMatrix(randomMatrix(size, size, cardinality)); - + mmulMediumDensity(A, B); // Benchmark the multiplication yield () => { - do_not_optimize(A._mmulMediumDensity(B)); + do_not_optimize(mmulMediumDensity(A, B)); }; // Explicit cleanup @@ -94,19 +119,19 @@ lineplot(() => { const results = await run({ // Force GC between every benchmark gc: true, - // // More samples for statistical significance - // min_samples: 20, - // max_samples: 200, - // // Longer warmup to stabilize CPU state - // warmup_samples: 10, - // warmup_threshold: 100, // ms + // More samples for statistical significance + min_samples: 20, + max_samples: 200, + // Longer warmup to stabilize CPU state + warmup_samples: 10, + warmup_threshold: 100, // ms // Longer minimum time for stable measurements - // min_cpu_time: 2000, // 2 seconds minimum - // // Batch settings to reduce variance - // batch_samples: 5, - // batch_threshold: 10, // ms - // // Enable colors - // colors: true, + min_cpu_time: 2000, // 2 seconds minimum + // Batch settings to reduce variance + batch_samples: 5, + batch_threshold: 10, // ms + // Enable colors + colors: true, }); // Process and store results @@ -134,6 +159,6 @@ for (const benchmark of results.benchmarks) { // Save results to JSON file await writeFile( - path.join(import.meta.dirname, `benchmark-results.json`), + path.join(import.meta.dirname, `benchmark-results-${dimensions[0]}.json`), JSON.stringify(processedResults, null, 2), ); From 5a7b46ed489eb4f773acb120682339c2a21604c4 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Micha=C3=ABl=20Zasso?= Date: Fri, 4 Jul 2025 10:58:34 +0200 Subject: [PATCH 33/33] fix eslint import order --- benchmark/combinationOfCardinalityAndSize.js | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/benchmark/combinationOfCardinalityAndSize.js b/benchmark/combinationOfCardinalityAndSize.js index d9217fa..e2b3907 100644 --- a/benchmark/combinationOfCardinalityAndSize.js +++ b/benchmark/combinationOfCardinalityAndSize.js @@ -6,11 +6,11 @@ import { Matrix } from 'ml-matrix'; import { xSequentialFillFromStep } from 'ml-spectra-processing'; import { SparseMatrix } from '../src/index.js'; - -import { randomMatrix } from './utils/randomMatrix.js'; +import { mmulLowDensity } from '../src/utils/mmulLowDensity.js'; import { mmulMediumDensity } from '../src/utils/mmulMediumDensity.js'; import { mmulSmall } from '../src/utils/mmulSmall.js'; -import { mmulLowDensity } from '../src/utils/mmulLowDensity.js'; + +import { randomMatrix } from './utils/randomMatrix.js'; /* eslint func-names: 0