[Numpy-discussion] Can this function by Numpy-ized?
John J. Lee
jjl at pobox.com
Sun Aug 26 13:27:08 EDT 2001
On Sun, 26 Aug 2001, Rob wrote:
> I finally got my FEM EM code working. I profiled it and this function
> uses up a big hunk of time. It performs gaussian integration over a
> triangle. I am trying to figure out how to slice the arrays so as to
> push it down into the C level. Does anyone have any ideas? Thanks,
> Rob.
>
> ps. it looks to be intractible to me. Maybe I need to look at writing a
> C extension. I've never done that before.
[...]
I'm not a good enough ufunc / array function hacker to come up with
anything, but I doubt it, and I wonder how much time it would save anyway,
given the size of the arrays involved. From the look of your comments
before the function, it looks like a) you're a C programmer, not quite
comfortable with Python yet, or b) you wrote it in C first, then moved it
into Python. If the latter, you might want to try wrapping that C
function with SWIG, though to be honest I'm not sure the overhead of
learning to use SWIG is any less than for writing a C extension manually
(but less error-prone I expect, and less work if you have a lot of C to
wrap). I think the Scipy project has some typemaps for passing Numeric
arrays; if not, I've seen some collections of SWIG typemaps for Numeric
somewhere...
BTW, the best place to put those comments is in a docstring. Here is a
slightly more Pythonically-formatted version (I must be bored today):
def ComputeGaussQuadPoint(QuadPoint, a):
"""Return coordinates of 7-point Gauss nodes of a triangular patch.
QuadPoint: node index, from 0 to 6
a: triangular patch?
"""
SrcPointCol=zeros((3),Float)
tn = a.TrngleNode # the three nodes of the triangular patch
SrcPointCol[0] = (a.Qpnt[QuadPoint,0]*a.NodeCord[tn[0],0] +
a.Qpnt[QuadPoint,1]*a.NodeCord[tn[1],0] +
a.Qpnt[QuadPoint,2]*a.NodeCord[tn[2],0])
SrcPointCol[1] = (a.Qpnt[QuadPoint,0]*a.NodeCord[tn[0],1] +
a.Qpnt[QuadPoint,1]*a.NodeCord[tn[1],1] +
a.Qpnt[QuadPoint,2]*a.NodeCord[tn[2],1])
SrcPointCol[2] = (a.Qpnt[QuadPoint,0]*a.NodeCord[tn[0],2] +
a.Qpnt[QuadPoint,1]*a.NodeCord[tn[1],2] +
a.Qpnt[QuadPoint,2]*a.NodeCord[tn[2],2])
return SrcPointCol
John
More information about the NumPy-Discussion
mailing list