The "star" problem, was Re: [eigen] ISO C++ working groups on Linear Algebra / Machine Learning |
[ Thread Index | Date Index | More lists.tuxfamily.org/eigen Archives ]
Summary: The right choice for what `operator*` should do depends on how many dimensions the language wants to elegantly handle. Matlab and Eigen have chosen the path of notation elegance. Numpy has chosen the path of higher dimensions.
Using operator* for matrix multiplication is concise and aligns
with how linear algebra is taught and thought. i.e. using vectors
and matrices. Eigen, which deals exclusively with rows and columns[1], feels
quite natural for linear algebra and optimization problems. It
matches how most of us (except Einstein) think about the
operations involved.
Using higher dimensions elegantly comes at a cost of making some linear algebra operations (specifically matrix multiplication) more awkward to describe and implement. If a library wants to handle the regimes {dims<=2, dims>2}, it make little sense to use operator* for matrix multiplication. The meaning is unclear in higher dimensions. There are, after all, only two sides of a '*' character: left and
right. In numpy and scipy, there is a conceptual divide between the data (numpy arrays) and the operations (e.g. scipy.linalg, numpy matrices[3]). This is easy to forget since the linear algebra functions [4] just work as expected, assuming you send in 1d and 2d data. Unfortunately, this conceptual divide makes it difficult to define `operator*` as matrix multiplication. For example, if `a` and `b` are 1d arrays of the same length, should `a*b` be an inner product or an outer product?Numpy arrays allow one to store data on meaningful hypercube axes and perform inner products and summations over those axes. In a cell tower simulation[2], this helped me consider channel responses for different users, antennas, base stations, paths, and multiple simulation realizations as axes on a hypercube. This would also be possible in Matlab, but it would be clunky: lots of `transpose`, `reshape`, and `bsxfun`. Some of you might be thinking that Matlab handles higher
dimensions just fine -- try adding a singleton axis at the end:
`size( reshape( 1:10 ,5,2,1 ) )`. I tried to be fairly objective and unbiased in my discussion
above, but I do have a preference, of course. That preference is based on how *I* use my tools in the years
I've spent with Matlab, Numpy, and Eigen. Explicit is better than implicit.
-- Mark Borgerding P.S. I do not belong to any of the working groups below. If you
do and think this worth sharing, please do so. P.P.S. Anyone learning numpy after Matlab should repeat the mantra "numpy arrays are not matrices" a thousand times a day until it sticks.
[1] There is the notable exception of the unsupported tensor
module. I've only used it indirectly from TensorFlow, so I cannot
comment further. [3] FWIW, the numpy matrix class overloads `operator*`, but few
people use the matrix class, instead opting for arrays and
`matmul`, `dot`, or `einsum` [4] https://docs.scipy.org/doc/scipy/reference/linalg.html
On 2/15/19 5:51 PM, Matthieu Brucher
wrote:
|
Mail converted by MHonArc 2.6.19+ | http://listengine.tuxfamily.org/ |