Questions about Goals and Future Direction for Matrix Multiplication

[epsilonError]

Hello Everyone,

I worked on some tests and type fixes for multiplyMatrices.js and wanted to share some edge cases and ask about future implementation goals and priorities.

I previously mentioned I was going to open an issue to outline what I found while working on tests and types, but there's enough to talk about that I thought it would be better to try opening a discussion instead.

After some consensus is reached, I'll breakdown the decisions into individual issues and more focused pull requests.

Reference Implementation

I'm testing possible implementation details and running tests with an adaptation of color.js's multiplyMatrices.js in a private typescript project. The typescript code is available in a github gist (hover to see when last updated).

I already have testing and code coverage setup for that project so it's easier for me to play around with things there. And I plan on porting the relevant tests to hTest and color.js in future.

The private project also includes tests for simulating prota-/deutera-/trita-nopia on sRGB colors using color.js and compares reference results from DaltonLens-Python. So I've got quite a lot of tests already. (After making the simulation work with non-sRGB colors I plan on open sourcing it, but that is a long way off and not part of this discussion.)

Observations

Priority

What I gleaned from a previous version of the API docs is that multiplyMatrices.js isn't focused on ensuring the correctness of every arbitrary matrix product. The priority seemed to be 1) correctness for common use cases, and 2) performance to reduce the overhead of the frequent transformation of coordinates from one Color Space to another.

Arguments & Return Types

In color.js (as of 2024-08-10) there are 25 non-test uses of multiplyMatrices(), roughly the types work out to be:

Assigned/Result Type	A	B	Returns
[number number number]		11	7
[number number number][]
number[]		10	16
number[][]	25	2	2
any		2

My Guess of Type	Times	Out of
any seems to be 3-Array:	2	2
number[] seems to be 1×3:
number[] seems to be 3×1:	10	26
number[] seems to be 3-Array:	10	26
number[][] seems to be 3×3:	25	29

I did a rough scan of the code to get these numbers, they aren't definitive.

Current multiplyMatrices Implementation Notes

1. Lines 65–69 seem unreachable. At that point, AM is a row-ordered matrix so the row argument of the first map will always be an Array. (Code Coverage never reached this code in my tests.)
2. Line 73 uses (col[i] || 0) to handle possible undefined access of the col variable. This is guaranteed to be an invalid product if the || 0 is used. Another possible choice is using NaN instead of 0.{See comment below.}
3. Since the for loop on line 72 loops until the length of the row in the A matrix, the returned value's shape is always constrained by A. Iterating based on col.length, and using (row[i] || 0), returns the same correct answers and is also a reasonable choice.
4. The general types in the function overloads work well, but it could be worthwhile to add a TransformationMatrix (3×3) type and use it along with a non-null Coords type to specify precise return types.
5. Along with 3 and 4 above, the shape of the return type will always be constrained by either A's rows or B's columns. Which means either A or B is guaranteed to constrain the final returned shape.
6. A vector × vector generally returns a scalar, but multiplyMatrices returns a [number] result.
7. Overly specified vectors arguments are permissible and fulfill the parameter type requirements (i.e., [[]], [[1]], [[1],[2],[3]], etc.)
8. Overly specified vector arguments don't fulfill the return type guarantee. They are all number[][] type but are mathematically equivalent to vectors. Notationally 1 is the same as [1] is the same as [[1]], but they are not the same computationally. Computational equivalence requires properly packing and unpacking these nested arrays.
9. Empty B vectors throw an error (can't call .length on undefined). We could throw our own errors, attempt to return an obviously incorrect value (but type/shape correct) result, or leave as is for performance reasons
10. If B is [] then calling map on it will assign [] to BM, rather than the expected [[]].
11. m is actually assigned the value of n, when A is a number[].

The Future

My Next Steps Recommendation

• Add more tests to account for notes-1,2,3,5,9,10,11
• Add a TransformationMatrix type and use it and/or non-null Coords type to highlight "happy path" usages of multiplyMatrices (ref notes-4,5)

Needs Discussion

• How do we want to pack and/or unpack mathematically equivalent nested arrays? (Some Options: Unpack to simplest equivalent. Pack to match provided arguments return type.) [related to notes-6,7,8]
• How much error checking and handling do we want to do? (Some Checking Options: None, An empty edge-case, All possible empty edge-cases, full dimension checking. Some Handling Options: Throw error for some/each check, return obviously incorrect values (0/NaN) in the proper return shape) [related to notes-9,2,3]
• Since Coords includes null and numbers as acceptable values, should multiplyMatrices accept them as well? (I'm not sure how multiplyMatrices could handle this correctly for all cases. So it's probably better for each Color Space to resolve null values.) [related to notes-4,2,9]

For Now

Thanks for reading and let know what you think. I'll probably focus on adding more test in the near future. But I'm open to hearing other perspectives and recommendations.

[epsilonError]

2. Line 73 uses (col[i] || 0) to handle possible undefined access of the col variable. This is guaranteed to be an invalid product if the || 0 is used. Another possible choice is using NaN instead of 0.{See comment below.}

NaN is not a reasonable choice if || is used for comparison. 0 is a falsy value and so all 0's in a matrix would be turned into NaN. So using NaN wouldn't be a good choice.

[lloydk]

Thanks for the detailed analysis. Looks like some of the types need to be cleaned up.

I was doing some profiling of a benchmark that used Color.js on the weekend and found that multiplyMatrices was the function that contributed the most amount of time to the benchmark. Using a specialized function to transform a [number, number, number] by a 3x3 matrix resulted in a significant performance improvement. Here is the PR with the changes and benchmark timings: #585

Since Coords includes null and numbers as acceptable values, should multiplyMatrices accept them as well? (I'm not sure how multiplyMatrices could handle this correctly for all cases. So it's probably better for each Color Space to resolve null values.)

I agree that null should be resolved to 0 before calling and of the matrix multiplication functions.

[epsilonError]

I was doing some profiling of a benchmark that used Color.js on the weekend and found that multiplyMatrices was the function that contributed the most amount of time to the benchmark. Using a specialized function to transform a [number, number, number] by a 3x3 matrix resulted in a significant performance improvement. Here is the PR with the changes and benchmark timings: #585

That's a good use of the specialized function, and with type-checking a great tradeoff for the performance.

It makes sense to me that multiplyMatrices is a bottleneck, it's generic enough to handle a lot of things and is consistently called whenever a color space is traversed.

As just a nit-pick, I prefer the semantics of multiplyMatrices' parameters over transform's. Even though they are functionally not much different. transform reads to me as transform input by matrix while optionally using some existing vector. While multiplyMatrices hints at a chain of matrices all lined up together and finally a column vector at the end reifies the abstract linear transformations into a specific position. But I'm a sucker for a reduce/scan pattern so I like imagining 𝐀𝐁𝐂𝐃𝐱ᵀ over 𝑓(𝐱, 𝐀𝐁𝐂𝐃) in my head.

I'm happy to keep working on multiplyMatrices if a generic way to do matrix multiplication with scalars, vectors, and/or matrices is needed. But I suspect a specialized function for 3x3 matrices and vectors will always outperform it.

[epsilonError]

I found some benchmarks on matrix math. It only benchmarks element-wise addition of matrices but could be useful in structuring future specialized functions.

[epsilonError]

I took a look at the CSS 4 Editor's Draft's Color Conversion Code section
for some background intentions for matrix multiplication.

My understanding of the linked article on column-major ordering, the wikipedia article on matrix multiplication,
and my memory of undergraduate linear algebra, matrices work something like this:

Matrices are 2 dimensional rectangular collections and we're only interested in collections of numbers.

So matrices can be described as tables where each cell in the table can contain a collected number.

Matrix	col 1	col 2	col 3
row 1	𝑒₁₁	𝑒₁₂	𝑒₁₃
row 2	𝑒₂₁	𝑒₂₂	𝑒₂₃
row 3	𝑒₃₁	𝑒₃₂	𝑒₃₃

Mathematically

Mathematically, matrices can be accessed or manipulated by element, row, or column. All 3 are equally important and each can be preferred depending how the matrix will be used.

The trouble comes with how to use vectors and matrices together. A matrix can be described as made up of both rows and columns at the same time. But a vector is 1 dimensional and can only be described as either a row or column.

The transpose operator easily changes a row to a column and column to a row. So for consistency it needs to be stated if a vector is by default a representation of a row or a column. This is where row-major order and column-major order come in. The -major order states the default format of vectors.

Row vectors and column vectors each require a specific multiplication order to be equivalent:

• Row vectors are pre-multiplied, i.e. they are multiplied first (𝐯⋅𝐌).
• Column vectors are post-multiplied, i.e. they are multiplied last (𝐌⋅𝐯).

This importance of multiplication ordering especially applies to successive multiplications, but multiplyMatrices only takes two arguments to multiply together so we can ignore that for now.

Computer Implementation

Computer implementations can extend the notational representation of row-major order and column-major order and actually use the -major ordering as a low level structure of the data in memory. There are performance reason and optimizations that can be made for how the data is accessed and stored in-memory, but we can ignore any of that since Color.js is beholden to JavaScript. (We could probably use TypedArrays to make our own low-level format, but I would suggest just nixing that idea right away.)

In Color.js, matrices are written as nested arrays. And with lines typed horizontally the inner array is always written as a row. So the js matrices are written as if they are row-major ordered.

But plain arrays in js don't contain enough information to say if they represent row, or column vectors. To absolutely specify the type requires representing the vector as a nested array. So, to force the array [1, 2, 3] into a column or row vector would be:

Array	Row Vector	Column Vector
[1, 2, 3]	[[1, 2, 3]]	[[1], [2], [3]]

multiplyMatrices.js

The mathematical and implementational details gives some options for how multiplyMatrices could work:

1. Treat an array as appropriate vector based on argument position (First arg = row, second arg = column)
2. Only support/force use of row-major order or column-major order (always expect/treat 1D arrays as one of the major orders, may require an additional transpose function)
3. Pervasively allow use of both row-major order and column-major order (would mean adapting vector based on argument position above, and creating both pre- and post- multiplication specialized functions for any performance-focused functions)
4. Focus performant multiplication implementations on either row-major order or column-major order, not both

And if the multiplication functionality should additionally support,
5. vector times vector returning a scalar
6. performance focused vector multiplications (vector times vector, and/or vector times scalar)

Color.js multiplication order

Lastly, another factor needs be considered when deciding between pre-, post- (or both) multiplication. JavaScript numbers are not equivalent to mathematical Real Numbers. The available precision in them should be good enough for Color.js's needs. But since they aren't Reals the commutative property of multiplication can breakdown, meaning the order in which numbers are multiplied can effect the final result. And this effect can compound with the number of multiplications performed.

So picking/prioritizing pre/post-multiplication is necessary for highly consistent results. You can't avoid the possible compounding precision loss, but it is possible to get decently consistent results by judiciously choosing between pre/post-multiplication.

[facelessuser]

5. vector times vector returning a scalar

This would seem to be the far more likely need in this library, essentially a vector dot. If you want to multiply it some other way, matrixMultiply would likely not be the choice, and some generic multiplication function could be used instead.

1. Treat an array as appropriate vector based on argument position (First arg = row, second arg = column)

The above is essentially what Numpy does (a scientific library in Python).

The library is row major by default, but matmul(mat, vec) treats vec like a column vector, and matmul(vec, mat) treats the operation as a row major operation.

>>> import numpy as np
>>> np.matmul([1, 2, 3], [[1, 2, 3], [4, 5, 6], [7, 8, 9]])
array([30, 36, 42])
>>> np.matmul([[1, 2, 3], [4, 5, 6], [7, 8, 9]], [1, 2, 3])
array([14, 32, 50])
>>> np.matmul([1, 2, 3], [4, 5, 6])
32

I think they don't allow scalar operations in this function as a scalar is obviously not a matrix, but we could do whatever we want.

This is just my thoughts. I'm familiar with this approach, so I'm biased, but I'll adapt to whatever.

[epsilonError]

Just submitted a pull request that includes the flexible function interface like Numpy's matmul uses. I did also include multiplying vectors and returning a scalar in it as well.

At this point, I'm thinking multiplyMatrices will end up being a slow generic matrix multiplication operator. Probably mostly used in a prototype/rapid iteration way.

I suspect the speed-up @lloydk got from using specialized Vector3 and Matrix3x3 multiplication means the standard use will be specialized functions. And the existing dot3 for Vector3 dot products, and multiply_v3_m3x3 already supersedes almost all uses of multiplyMatrices.

Once a specialized multiply_m3x3_v3 and multiply_m3x3 are implemented, there isn't really anything for multiplyMatrices to do in Color.js. But I do think exposing a general matrix multiplication function would be useful for people who don't know about the specialized functions. Especially if the docs on multiplyMatrices point out the existence of the specialized functions.