diffkk

[bruceravel]

This may be the unwelcome feature request of the week.

It would be useful to have the functionality of diffkk implemented in Larch.

[newville]

Agreed!

[maurov]

Good idea. If it can helps, here the one I got from you (Matt?) some years ago ;)

https://gist.github.com/33997083e96ab4036fe7.git

Mauro

[bruceravel]

I found myself with a bit of time this week (which is what happens when you are between synchrotrons!), so I started taking a stab at this. I am working on it in the "diffkk" branch.

It's pretty fun. The funnest part was vectorizing the MacLaurin series calculation. Doing so reduced the time it takes to transform copper foil from about 15 seconds to about 1.6 seconds, with no difference in the output array. Neat-o!

My motivation for all this is needing f'(E) for Ce L edge data. I am not 100% pleased with how the L series data matches with the CL f2(E) function, so I have more playing around to do before I make a pull request.

Here's a teaser, copper foil:

[newville]

Cool! If speed of the KK transform is really an issue, we could just use Fortran. We're going to need Fortran-compiled code for lots of small tasks anyway, and it's not really that hard to support (well, except that Mac OS X is painful, even as I write this from a Macbook....). Of course it's easiest to use Python until it's obviously a problem, but I have to admit I'm surprised at how little compiled code we have.

[bruceravel]

The following is a flawed comparison, but ...

Running the diffkk program from ifeffit like this:

for i in $(seq 1 10); do
    diffkk
done

Takes 9.1 seconds of user+sys time. A Fortran program that just did diffkk without using the ifeffit library would be faster.

Starting larch, then doing run 'dkk.lar' from the command line:

import time
a=read_ascii('examples/xafsdata/cu_10k.xmu')
b=diffkk(a.energy, a.mu, e0=8979, z=29, order=4, form='mback')
start = time.clock()
for i in range(10):
    b.kktrans()
endfor
finish = time.clock()
print finish - start

takes 15.5 seconds. Here, I have excluded larch's start-up time.

So, you are right that using an external library would be wiser for an application that does hundreds or thousands of KK transforms. You'd get less than order of magnitude, but probably a factor of 3 or 4 improvement. My immediate need is to do one :smile: It seems to be working pretty well.

The more interesting comparison is with the algorithm that @maurov linked to above -- the one from CORE. That's purely scalar -- a direct translation of Fortran into python -- and is a factor of 10 slower than using numpy's vectorization.

My take-away message is that these vectorized math systems like numpy are super-cool! But we already knew that.

[bruceravel]

See #110