Re: [Numpy-discussion] Rookie problems - Why is C-code much faster?

[Zachary P. <zp...@st...>]

Mads,

The game with numpy, just as it is with Matlab or any other =20
interpreted numeric environment, is to try push as much of the =20
looping down into the C code as you can. This is because, as you now =20
know, compiled C can loop much faster than interpreted python.

A simple example for averaging 1000 (x,y,z) points:

print data.shape
(1000, 3)
# bad: explicit for loop in python
avg =3D numpy.zeros(3, numpy.float_)
for i in data: avg +=3D i
avg /=3D 1000.0

# good: implicit for loop in C
avg =3D numpy.add.reduce(data, axis =3D 0)
avg /=3D 1000.0

In your case, instead of explicitly looping through each point, why =20
not do the calculations in parallel, operating on entire vectors of =20
points at one time? Then the looping is "pushed down" into compiled C =20=

code. Or if you're really lucky, it's pushed all the way down to the =20
vector math units on your cpu if you have a good BLAS or whatever =20
installed.

Zach Pincus

Program in Biomedical Informatics and Department of Biochemistry
Stanford University School of Medicine



> 2. Here is the main loop for finding all possible pair distances,
>    which corresponds to a loop over the upper triangular part of a
>    square matrix
>
>      # Loop over all particles
>      for i in range(n-1):
>         dx =3D x[i+1:] - x[i]
>         dy =3D y[i+1:] - y[i]
>
>         dx -=3D box*rint(dx/box)
>         dy -=3D box*rint(dy/box)
>
>         r2 =3D dx**2 + dy**2      # square of dist. between points
>
> where x and y contain the positions of the particles. A naive
> implementation in C is
>
>
>     // loop over all particles
>     for (int i=3D0; i<n-1; i++) {
>         for (int j=3Di+1; j<n; j++) {
>             dx =3D x[j] - x[i];
>             dy =3D y[j] - y[i];
>
>             dx -=3D box*rint(dx/box);
>             dy -=3D box*rint(dy/box);
>
>             r2 =3D dx*dx + dy*dy;
>         }
>     }
>
> For n =3D 2500 particles, i.e. 3123750 particle pairs, the C loop is
> app. 10 times faster than the Python/Numeric counterpart. This is of
> course not satisfactory.
>
> Are there any things I am doing completely wrong here, basic
> approaches completely misunderstood, misuses etc?
>
> Any suggestions, guidelines, hints are most welcome.
>
> Best regards,
>
> Mads Ipsen
>
>
> +---------------------------------+-------------------------+
> | Mads Ipsen                      |                         |
> | Dept. of Chemistry              | phone:     +45-35320220 |
> | H.C.=D8rsted Institute            | fax:       +45-35320322 |
> | Universitetsparken 5            |                         |
> | DK-2100 Copenhagen =D8, Denmark   | mp...@os...         |
> +---------------------------------+-------------------------+
>
>
> -------------------------------------------------------
> This SF.net email is sponsored by: Splunk Inc. Do you grep through =20
> log files
> for problems?  Stop!  Download the new AJAX search engine that makes
> searching your log files as easy as surfing the  web.  DOWNLOAD =20
> SPLUNK!
> http://sel.as-us.falkag.net/sel?cmd=3Dlnk&kid=103432&bid#0486&dat=121642=

> _______________________________________________
> Numpy-discussion mailing list
> Num...@li...
> https://lists.sourceforge.net/lists/listinfo/numpy-discussion