This comes following Michael Croucher’s post MATLAB: Repeated multiplication is faster than integer powers.

He noticed that in MATLAB doing element-wise exponentiation of matrices (using .^ operator) was slower than explicitly (element-wisely) multiplying the matrix by itself. I’m pretty sure this is only true for small exponents, anyway. He also wondered how was it in other languages.

Well, I did my job and tried it in Julia:

function testpow() #% function to compare integer powers with repeated multiplication a=rand(1,10000000)*10; println("speed of ^4 using pow") tic(); test4p = a.^4; toc() println("speed of ^4 using multiply") tic(); test4m = a.*a.*a.*a; toc() println("maximum difference in results"); max_diff = maximum(abs(test4p-test4m)) end

Calling testpow() , on my Intel(R) Core(TM) i3 M380 @ 2.53GHz laptop, it outputs (after warm-up)

speed of ^4 using pow elapsed time: 1.074429726 seconds speed of ^4 using multiply elapsed time: 0.40391358 seconds maximum difference in results 1.8189894035458565e-12

So it is the same, except for a different slowness factor and for the absolute values.

But then I wondered — what’s if I just (re)define the .^ operator for integer exponents to perform chain multiplications? So I wrote this “Cristóvão’s style” Julia code:

import Base.(.^) function .^{T<:FloatingPoint, N}(A::Array{T, N}, x::Integer) if abs(x) > 42 # the "Answer to the Ultimate Question of Life, the Universe, and Everything" A.^float(x) elseif x > 1 B = similar(A) @inbounds for i in 1:length(A) B[i] = A[i] for k in 1:x-1 B[i] *= A[i] end end B elseif x < 0 B = similar(A) @inbounds for i in 1:length(A) B[i] = one(T) for k in 1:abs(x) B[i] *= A[i] end B[i] \= one(T) end B elseif x == 1 copy(A) else # x == 0 ones(A) end end

By running this before the testpow() function definition (and not touching it), I get the following output

speed of ^4 using pow elapsed time: 0.100119984 seconds <— notice this speed of ^4 using multiply elapsed time: 0.40944625 seconds maximum difference in results 0.0

Cool! A 10x speed-up comparing to default .^ .

I think it is really nice to be able to write improved different algorithms in the very same language and then get improvements changes on other code without having to touch it.

Although I didn’t tried it, I guess one cannot achieve this in MATLAB, or Python 😛

UPDATE: Using real world real numbers (floating point numbers), precision can be lost in chain multiplication exponential implementation.