Benchmark MATLAB & Julia for Matrix Operations

[RoyiAvital]

I am not sure about the reason for that as I am not aware of what Julia is doing internally.
@ChrisRackauckas said it is a known open problem with fusion which doesn’t generate a code which hints the compiler there is no aliasing. Others said maybe MATLAB parallelizes this while Julia doesn’t (Though it seems like a memory bounded operation so I’m not sure parallelization can explain result).

This is the line of the function - mA = (scalarA .* mX) .+ (scalarB .* mY);.
If there is mistake on my side how it is right to write the function, please advise and I will change.

Regarding MATLAB, this is something they heavily invested in the few latest releases - Working and multi terms assignments. They might even infer fusion by themselves in this case.

I agree it doesn’t look reasonable. I will run it again.

@antoine-levitt

Have you seen the code?
It runs 7 times and the measurement is the median.

Oops I just saw that, sorry. However, it’s better for this kind of benchmarking to be done closer to the actual computation, and not do anything between two runs (in particular memory allocations). Just use BenchmarkTools, which does all this itself and avoids some common pitfalls.

I compared my technique of measuring to BenchMarksTools above. See that results are comparable within a reasonable accuracy of measuring CPU.

@antoine-levitt

Regarding the code, Please review the code and point me to any place you think I should add it and I will.
Just out of curiousity, I’d understand that for existing array it will be important to barodcast / fuse teh assignment with .= . But I’d think that for a new generated array it won’t. Not that is the reason they are missing. I was just not aware enough of that. Not very experienced with Julia.

I told you, just add dots. When you do a = b.+c, you allocate an array for a, which is slow in julia (slower than in matlab, which has an allocator optimized for that extremely common pattern). The way to get optimal performance in julia is to track down systematically allocations. a .= b.+c does not allocate anything (well, at least, no arrays).

Dots were were there to begin with (With one exception I missed mE = exp.(-(mA .^ 2)); - No . before the -).
What I wrote that it didn’t make sense to me you need to have . before the assignment operator for a new allocated array.

For instance, no reason mZ .= mX .+ mY; will be faster than mZ = mX .+ mY; for the case mZ is not defined. That was my logic. I even tried it on Julia.
It seems to be Julia’s logic as well, because you can only use mZ .= mX .+ mY; in case mZ was pre allocated.

In my code the left side hand isn’t allocated. Hence no need for . before =. Only between operations on the right hand side. The code was like that to begin with.