Linear combinations in python / numpy

Question

Hi,

I'm not sure if this is a stupid question or not.

Let's say I have 3 arrays of numpy, A1, A2, A3 and 3 floats, c1, c2, c3

and I would like to evaluate B = A1 * c1 + A2 * c2 + A3 * c3

will numpy calculate this for example

E1 = A1*c1 E2 = A2*c2 E3 = A3*c3 D1 = E1+E2 B = D1+E3 

or is he smarter? In C ++, I had a neat way to abstract this operation.

I have defined a number of common LC, LC template functions for a linear combination, such as:

 template<class T,class D> void LC( T & R, T & L0,D C0, T & L1,D C1, T & L2,D C2) { R = L0*C0 +L1*C1 +L2*C2; } 

and then specializes for different types,

for example, for an array, the code looked like

 for (int i=0; i<L0.length; i++) R.array[i] = L0.array[i]*C0 + L1.array[i]*C1 + L2.array[i]*C2; 

which avoids creating new intermediate arrays.

This may seem messy, but it works very well.

I could do something similar in python, but I'm not sure if its nescesary.

Thanks in advance for your understanding. -nick

Answer 1

While numpy , theoretically, can always update its internal functions at any time to perform miraculous optimizations, this is currently not the case: B = A1*c1 + A2*c2 + A3*c3 really produces and then discards the intermediate temporary arrays ("spends" some auxiliary memory, of course - nothing more).

B = A1 * c1 followed by B += A2 * c2; B += A3 * c3 B += A2 * c2; B += A3 * c3 , again at this time, so you don’t have to spend part of this temporary memory.

Of course, you can tell the difference only if you work in an environment with scarce real memory (where some of this auxiliary memory is simply virtual and leads to page errors) and for large enough arrays, “waste” all real memory, and then some. In such extreme conditions, however, a little refactoring may buy you some performance.

Answer 2

This is the idea of numexpr (a quick expression parser for numeric arrays for Python and NumPy). You can try this package before compiling your routines.