[Numpy-discussion] [C++-sig] Overloading sqrt(5.5)*myvector

Bruce Sherwood Bruce_Sherwood@ncsu....
Sat Dec 29 17:00:53 CST 2007

Okay, I've implemented the scheme below that was proposed by Scott 
Daniels on the VPython mailing list, and it solves my problem. It's also 
much faster than using numpy directly: even with the "def "and "if" 
overhead: sqrt(scalar) is over 3 times faster than the numpy sqrt, and 
sqrt(array) is very nearly as fast as the numpy sqrt.

Thanks to those who made suggestions. There remains the question of why 
operator overloading of the kind I've described worked with Numeric and 
Boost but not with numpy and Boost. There is also the question of 
whether it would pay for numpy to make what is probably an exceedingly 
fast check and do much faster calculations of sqrt(scalar) and other 
such mathematical functions.

Bruce Sherwood

Bruce Sherwood wrote:
> I found by timing measurements that a faster scheme with less penalty 
> for the case of sqrt(array) looks like this:
> def sqrt(x):
>      if type(x) is float: return mathsqrt(x)
>      return numpysqrt(x)
> Bruce Sherwood wrote:
>> Roman Yakovenko wrote:
>>> On Dec 29, 2007 7:47 AM, Bruce Sherwood <Bruce_Sherwood@ncsu.edu> wrote:
>>>> I realized belatedly that I should upgrade from Boost 1.33 to 1.34.
>>>> Alas, that didn't cure my problem.
>>> Can you post small and complete example of what you are trying to 
>>> achieve?
>> I don't have a "small and complete" example available, but I'll 
>> summarize from earlier posts. VPython (vpython.org) has its own vector 
>> class to mimic the properties of 3D vectors used in physics, in the 
>> service of easy creation of 3D animations. There is a beta version 
>> which imports numpy and uses it internally; the older production 
>> version uses Numeric. Boost python and thread libraries are used to 
>> connect the C++ VPython code to Python.
>> There is operator overloading that includes scalar*vector and 
>> vector*scalar, both producing vector. With Numeric, sqrt produced a 
>> float, which was a scalar for the operator overloading. With numpy, 
>> sqrt produces a numpy.float64 which is caught by vector*scalar but not 
>> by scalar*vector, which means that scalar*vector produces an ndarray 
>> rather than a vector, which leads to a big performance hit in existing 
>> VPython programs. The overloading and Boost code is the same in the 
>> VPython/Numeric and VPython/numpy versions. I don't know whether the 
>> problem is with numpy or with Boost or with the combination of the two.
>> Here is the relevant part of the vector class:
>> inline vector
>> operator*( const double s) const throw()
>> { return vector( s*x, s*y, s*z); }
>> and here is the free function for right multiplication:
>> inline vector
>> operator*( const double& s, const vector& v)
>> { return vector( s*v.x, s*v.y, s*v.z); }
>> Maybe the problem is in the Boost definitions:
>> py::class_<vector>("vector", py::init< py::optional<double, double, 
>> double> >())
>>    .def( self * double())
>>    .def( double() * self)
>> Left multiplication is fine, but right multiplication isn't.
>> A colleague suggested the following Boost declarations but cautioned 
>> that he wasn't sure of the syntax for referring to operator, and 
>> indeed this doesn't compile:
>> .def( "__mul__", &vector::operator*(double), "Multiply vector times 
>> scalar")
>> .def( "__rmul__", &operator*(const double&, const vector&), "Multiply 
>> scalar times vector")
>> I would really appreciate a Boost or numpy expert being able to tell 
>> me what's wrong (if anything) with these forms. However, I may have a 
>> useful workaround as I described in a post to the numpy discussion 
>> list. A colleague suggested that I do something like this for sqrt and 
>> other such mathematical functions:
>> def sqrt(x):
>>   try: return mathsqrt(x)
>>   except TypeError: return numpysqrt(x)
>> That is, first try the simple case of a scalar argument, handled by 
>> the math module sqrt, and only use the numpy sqrt routine in the case 
>> of an array argument. Even with the overhead of the try/except 
>> machinery, one gets must faster square roots for scalar arguments this 
>> way than with the numpy sqrt.
>> Bruce Sherwood
> _______________________________________________
> Numpy-discussion mailing list
> Numpy-discussion@scipy.org
> http://projects.scipy.org/mailman/listinfo/numpy-discussion

More information about the Numpy-discussion mailing list