Linear Algebra extension: svd, inv, *, /, \ for QuantityArray

[ggebbie]

Thanks for the cool and promising package!

On Julia Discourse there seemed to be some interest in extending this package for "multidimensional analysis" (following George W. Hart) which I mostly interpret to be linear algebra utilities for the QuantityArray type. I don't know whether this fits with the vision of this package or whether it is best completed in a separate package. Any comments/advice on this PR are appreciated.

Here I make a start using the existing DynamicQuantitiesLinearAlgebraExt to add functionality for: inv, (\), svd, SVD, Diagonal, eigen, and det for QuantityArrays.
All functions here use the "algebraic" interpretation of dimensional analysis.

Because QuantityArrays represent a dimensionally uniform matrix, the number of applicable linear algebra methods is fairly restricted. In the parlance of Multidimensional Analysis, an array would be any combination of values and dimensions/units in a tabular form. A matrix or multipliable matrix would be a restricted dimensional structure such that multiplication and other linear algebra functions are possible. In this light, the QuantityArray name for uniform matrices is somewhat confusing.

I think that other multipliable matrices (i.e., non-uniform) would best be represented by a different struct with the dimensional range and domain as two of the fields. This goes beyond an extension to this package and I plan to experiment on this idea in a separate package.

[github-actions]

Benchmark Results

	main	a7acc5a...	main/a7acc5a341e7ca...
Quantity/creation/Quantity(x)	2.79 ± 0 ns	3.42 ± 0.011 ns	0.818
Quantity/creation/Quantity(x, length=y)	3.41 ± 0.01 ns	3.11 ± 0.01 ns	1.1
Quantity/with_numbers/*real	3.11 ± 0.01 ns	3.11 ± 0.01 ns	1
Quantity/with_numbers/^int	8.37 ± 2.5 ns	8.05 ± 2.2 ns	1.04
Quantity/with_numbers/^int * real	8.05 ± 2.2 ns	8.67 ± 2.5 ns	0.929
Quantity/with_quantity/+y	4.04 ± 0.01 ns	4.04 ± 0.001 ns	1
Quantity/with_quantity//y	3.11 ± 0.001 ns	3.42 ± 0.011 ns	0.909
Quantity/with_self/dimension	2.79 ± 0.001 ns	3.1 ± 0.01 ns	0.9
Quantity/with_self/inv	3.13 ± 0.01 ns	3.11 ± 0 ns	1.01
Quantity/with_self/ustrip	2.79 ± 0 ns	3.71 ± 0.92 ns	0.754
QuantityArray/broadcasting/multi_array_of_quantities	0.143 ± 0.00069 ms	0.144 ± 0.00094 ms	0.997
QuantityArray/broadcasting/multi_normal_array	0.0529 ± 0.0002 ms	0.0558 ± 0.0029 ms	0.947
QuantityArray/broadcasting/multi_quantity_array	0.157 ± 0.00076 ms	0.157 ± 0.00062 ms	1
QuantityArray/broadcasting/x^2_array_of_quantities	22.8 ± 2.1 μs	22.9 ± 2.3 μs	0.994
QuantityArray/broadcasting/x^2_normal_array	4.56 ± 1.2 μs	4.7 ± 1.2 μs	0.97
QuantityArray/broadcasting/x^2_quantity_array	6.95 ± 0.22 μs	6.96 ± 0.26 μs	0.999
QuantityArray/broadcasting/x^4_array_of_quantities	0.0844 ± 0.00043 ms	0.0814 ± 0.0006 ms	1.04
QuantityArray/broadcasting/x^4_normal_array	0.0498 ± 0.00016 ms	0.0497 ± 0.003 ms	1
QuantityArray/broadcasting/x^4_quantity_array	0.05 ± 0.00017 ms	0.053 ± 0.0029 ms	0.943
time_to_load	0.135 ± 0.0017 s	0.144 ± 0.0014 s	0.941

Benchmark Plots

A plot of the benchmark results have been uploaded as an artifact to the workflow run for this PR.
Go to "Actions"->"Benchmark a pull request"->[the most recent run]->"Artifacts" (at the bottom).

[MilesCranmer]

Nice work!! I added some initial comments.

[ggebbie]

Failed Julia 1.6 runtest required some additional code for svd and Diagonal. I believe that the code should work for 1.9 and 1.6 (awaiting 1.6 GitHub Action to be sure). Involved adding some uglier code, however.

[MilesCranmer]

Edit: I guess we might need to add a special method for diagm as well, since it uses zero(::Type{T}) rather than zero(::T). So best to ustrip -> diagm -> QuantityArray(result, dimension(input))

[ggebbie]

Edit: I guess we might need to add a special method for diagm as well, since it uses zero(::Type{T}) rather than zero(::T). So best to ustrip -> diagm -> QuantityArray(result, dimension(input))

This is the same issue that led me to writing Diagonal. Right -- best to make a special diagm function.

[MilesCranmer]

I guess some of these issues could also be fixed if we could figure out a way to implement #76, right? Trading type instability with compatibility (though people could get type stability if they simply used zero(::T) instead... but at least it would avoid any errors).

[ggebbie]

Getting closer and hope to not have missed anything. Certainly ok if you edit this PR. I tried to test on Julia 1.6 locally but Pkg complained about Compat... package not being in registry so I quit.

[ggebbie]

Compatibility/composability with other packages is a big reason why I'm interested in this package. Unitful doesn't work with Turing, for example. It seems like there are issues with SciML and Unitful as well. Making some design choices in this package so that it is more usable in combination with other packages would be a good idea in my opinion.

[MilesCranmer]

Looking good. I took out the eigen stuff as it needs more tests but that could be added later.

There are a couple of branches that aren't tested but they are disambiguities, so I am okay with 99% coverage for the moment. We should get back to 100% later though.