[SciPy-user] Calling C code (Newbie Question)
eric jones
eric at enthought.com
Tue Sep 28 01:53:22 EDT 2004
Hey Tom,
Tom Kornack wrote:
> Hi Robert:
>
>> Your code:
>> sh = sum( cn*sin(outerproduct(om,time) ), 1)
>>
>> My code:
>> sh = dot(sin(outerproduct(time, om)), cn)
>
>
> Thanks for your suggestions. Using dot() is better than sum(),
> however, the outerproduct() alone gives me malloc errors and I have 2
> GB memory. I mean, it's a huge matrix that gets created when I have a
> million points. That's why I wanted to use C.
>
> The better question for this list would be: how do I Weave this in C?
>
I inserted the following replacements into your example app:
>> for i in range(numf):
>> s2[i] = sum( sin(2.*om[i]*time) )
>> c2[i] = sum( cos(2.*om[i]*time) )
>
import weave
code = """
for (int i=0; i < Nom[0]; i++)
{
s2[i]=0.0;
c2[i]=0.0;
float om_2 = 2.0*om[i];
for (int j=0; j < Ntime[0]; j++)
{
float ang = om_2*time[j];
s2[i] += sin(ang);
c2[i] += cos(ang);
}
}
"""
weave.inline(code,['om','s2','c2','time'])
and a little later in your code:
#for i in range(numf):
# sh[i] = sum( cn*sin(om[i]*time) )
# ch[i] = sum( cn*cos(om[i]*time) )
code = """
for (int i=0; i < Nom[0]; i++)
{
sh[i] = 0.0;
ch[i] = 0.0;
for (int j=0; j < Ntime[0]; j++)
{
float ang om[i]*time[j];
sh[i] += cn[j]*sin(ang);
ch[i] += cn[j]*cos(ang);
}
}
"""
weave.inline(code,['om','sh','ch','time','cn'])
[Note: Check the accuracy of the calculations as this was just a quick
example.]
I've attached the edited example file that I used for testing. It has
two calls to your function in the main() to get rid of any caching
affects from the first run or a weave.inline() function.
I got a speedup of about 1.75 when running your function with
multiple=0. How much faster is the pure C version? I played around a
little, and it looks like the sin() and cos() are the limiting factors
here. For N=4000, the function runs in about 3.95 seconds on my
machine. It is an order of magnitude faster if I comment out all the
inner loop trig calls with something like:
sh[i] += cn[j]*ang; //sin(ang);
ch[i] += cn[j]*ang //cos(ang);
This suggest that even the pure C version that uses the same approach
isn't going to be much faster since it also will be limited byt the
speed of the sin() and cos() calls.
>
> Sorry for the rather basic question.
Not at all.
eric
-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: ex.py
URL: <http://mail.scipy.org/pipermail/scipy-user/attachments/20040928/51067b0a/attachment.ksh>
More information about the SciPy-User
mailing list