Re: Order of matrix multiplication

[Nicholas Jankowski]

On Mon, Sep 12, 2016 at 9:05 AM, Bernardo Sulzbach <address@hidden> wrote:

I noticed today that Octave always performs multiplication from left to right.

This is, as one could expect, properly documented:


    -- Built-in Function: mtimes (X1, X2, ...)
        Return the matrix multiplication product of inputs.

        This function and x * y are equivalent.  If more arguments are
        given, the multiplication is applied cumulatively from left to
        right (...)


This is evidently not the best solution from a performance point of view if a multiplication substantially reduces the size of the product.

Locally, profiling
    m * m * m * m * x
against
    m * (m * (m * (m * x)));
when the size of m is [2048 2048] and the size of x is [2048 1] produces the following results

    Took 9.752126 seconds when using the default order.
    Took 0.035902 seconds when going from right to left.
    Maximum difference is 0.000732422.

I may be oversimplifying the analysis of numerical error here, but I would say that the second solution introduces less error as the number of floating-point operations is much smaller.

---

Understanding that trying to determine the "best" order of matrix multiplication during execution may not be a trivial task, I still wonder why Octave approach is so naive.

Additionally, would there be any function out there which I could use instead of mtimes which performed some analysis on the dimensions of the input matrices and optimized the order of multiplication without requiring me to do it manually?

I was curious so ran a similar test in Matlab R2016a

>> m = rand(2048,2048);
>> x = rand(2048,1);
>> tic;m * (m * (m * (m * x)));toc
Elapsed time is 0.021991 seconds.
>> tic;m*m*m*m*x;toc
Elapsed time is 1.339323 seconds.

ran each a few times and got no real change in the order of magnitude. my guess is both programs treat binary operators as such and follow standard order of operations.

I will point out that the octave mtimes function allows for more than two arguments, whereas Matlab gives a 'too many input arguments' error.  I've made use of this before when I'm not concerned about matlab compatibility because it means mtimes can be called with an unknown number of arguments using cell nomenclature, such as mtimes(mymatrices{:}) )  you can even reverse ordering by doing something like mtimes( mymatrices{end:-1:1} )

So mtimes definitely has the ability to be enhanced with better decision making for nargin>2. Would just require someone with the time and effort to implement:
 if (nargin >2) {do something smart} endif

nickj