Using Igor's sum() Command with Nested Math?

Share icon

Wed, 03/02/2016 - 01:16 pm

I'm writing some code in Igor, and porting to Python for other peoples' use. I'm also trying to optimize the code 'cause it takes minutes to run. Python has a nifty NumPy library that has a "sum" command that lets me do math within the sum which speeds things up because it avoids me creating a loop. Igor does not.

For example, in Python I can create a list [0, 1, 2, 3, 4] and b list [1, 2]. I can then say: c[0] = numpy.sum(exp(a[:]*b[0])) (I'm new at Python, so forgive some of the old-school explicit indexing.)

I wanted to do the same thing with Igor's sum command, but I get errors that it JUST wants a wave, such as c = sum(a), but I canNOT do something as simple as c = sum(a*2.0).

I'm sure I'm missing something obvious? The reason for this is it would remove one of my loops and hopefully speed things up (sped things up in my Python code by ~3x for my test dataset).
Just because it has the same name doesn't mean it does the same thing. If you check the documentation for Igor's sum command, you will find that it simply adds up all the elements in a wave.

You might be able to do what you want with a wave assignment, possibly using intermediate waves. Or perhaps MatrixOP can let you nest some of your computation. MatrixOP is designed for speed.

But no sum command like the one you want.

John Weeks
WaveMetrics, Inc.
support@wavemetrics.com

March 2, 2016 at 02:16 pm - Permalink

Yup, I realized that just 'cause it's the same name doesn't mean it's the same thing, that's why I was asking if there was some sort of equivalent. I'll look into MatrixOP. Using intermediate waves would add overhead that I'd rather avoid.

March 2, 2016 at 02:46 pm - Permalink

Quote:
Using intermediate waves would add overhead that I'd rather avoid.


I don't think it adds much overhead. You can test it with this (edited to use Duplicate instead of Make):

// Example:
//      Make/O/N=1E6 testWave = p; Demo(testWave)
Function Demo(w)
    Wave w
 
    Variable timerRefNum = StartMSTimer
 
    Duplicate/FREE w, temp          // Igor automatically kills free waves
    temp = 2 * w
    Variable result = sum(temp)
 
    Variable elapsed = StopMSTimer(timerRefNum) / 1E6
    Printf "Elapsed time: %g seconds\r", elapsed
 
    return result
End

March 2, 2016 at 03:00 pm - Permalink

astrostu wrote:

I wanted to do the same thing with Igor's sum command, but I get errors that it JUST wants a wave, such as c = sum(a), but I canNOT do something as simple as c = sum(a*2.0).


For a compound expression or to speed up execution you should look into MatrixOP.
MatrixOP/o cc=sum(a*2)


A.G.
WaveMetrics, Inc.

March 2, 2016 at 04:46 pm - Permalink

astrostu wrote:

For example, in Python I can create a list [0, 1, 2, 3, 4] and b list [1, 2]. I can then say: c[0] = numpy.sum(exp(a[:]*b[0])) ... The reason for this is it would remove one of my loops


No loop needed as far as I can tell. AG says MatrixOP for complex equations. Here's the simple example using implicit looping without MatrixOP.

make/N=5 alist={0,1,2,3,4}
make/N=2 blist={1,2}, clist
make/N=5/FREE tmp
tmp = exp(alist*blist[0])
clist[0] = sum(tmp)


A clever version might be able to do something with the second dimension in a matrix ...

make/N=(5,2)/FREE tmp
tmp = exp(alist*blist[q])
clist = sum(tmp[p][])


(though I would not swear by it)

--
J. J. Weimer
Chemistry / Chemical & Materials Engineering, UAH

March 2, 2016 at 08:32 pm - Permalink

Okay, MatrixOP definitely helped, though it took some massaging because it can't handle things like 1/wave (-> powR(wave,-1)) or scalar-wave (-> -wave+scalar). Timing test on some data took 0.413 seconds for the old method and 0.141 for the new (34% time).

When I run it once. Because this has to be bootstrapped to understand confidence intervals, running it 1000x under the old method (threaded on an 8-core machine) took 53.3 seconds, under the new method took 20.4 seconds (38% time).

Not bad.

Using a dataset 17.7x larger, old method took 16.1 seconds (39.0x longer than other dataset) to run for the main analysis, and 2223 seconds to bootstrap 1000x. New method took 6.31 seconds (39% time as old method, 44.8x longer than other dataset) and 1110 seconds to bootstrap 1000x (50.0% time as old method, 54.4x longer than other dataset).

March 3, 2016 at 01:48 pm - Permalink

Using a temporary wave to then sum() takes more time. The smaller dataset took 9% longer and the larger dataset took 15% longer, so I didn't bother trying to run the bootstrap.

Looks like MatrixOP is the way to go for this, at the moment.

March 3, 2016 at 02:06 pm - Permalink

Try this:
Make junk=p
matrixOP/O outw = rec(junk+1)
MatrixOP/O outw=powr(junk+1, -1)

rec() is MatrixOP's reciprocal function. Both of those lines produce what I think is the expected result.

John Weeks
WaveMetrics, Inc.
support@wavemetrics.com

March 3, 2016 at 02:33 pm - Permalink

johnweeks wrote:
Try this:
rec() is MatrixOP's reciprocal function. Both of those lines produce what I think is the expected result.


Cool! That sped things up another ~40-50%! So to run my code start to finish, on 6k points before today, 62 seconds, now 15 seconds (4.1x faster). On 12k points, 201 -> 43 seconds (4.7x faster). On 78k points, 1527 -> 555 (2.8x). On 107k points before today, 2239 seconds, now 821 (2.7x faster). So, it's not linear scaling, which my other version would do with speed, but it's a matrix operation so I wouldn't expect that, and it's still faster.

March 3, 2016 at 04:00 pm - Permalink