[Numpy-discussion] additional dtype argument to numpy.dot() (Re: Numpy-discussion Digest, Vol 29, Issue 48)
David Henderson
davidh@ipac.caltech....
Wed Feb 18 11:37:54 CST 2009
Hi Paul, list:
Thanks for the reply.
numpy.sum() does indeed have a dtype for the accumulator for the
sum. numpy.sum() does not implement an inner (dot) product, just a
straight summation.
The feature I'm requesting is to add a similar accumulator type
argument for numpy.dot().
After tracing some of the code, the there is a generic dot product
algorithm for the various types. If a BLAS library is available, the
BLAS routines will supplant the library, but only for the BLAS
routines that exist.
BLAS doesn't support long double, so only a subset of the different
calls for different types call BLAS.
The place where this feature comes into play is iterative refinement
of a linear algebra solution. To implement iterative refinement, an
error term is calculated as a matrix dot product with a higher
precision than the elements of the matrix.
The floating point hardware on the i386/PPC/Sparc is ready to go for
this application, its just that BLAS doesn't implement it.
The numpy.dot() function is missing the feature to supply the type of
the accumulator for the sum.
I'm looking for a way to prototype some routines to perform iterative
refinement, and python/numpy looks like an ideal platform.
numpy has _almost_ the feature set needed to implement this. Adding
the accumulator type argument to numpy.dot() is a needed feature.
code snippet exhibiting float128 functionality that already exists:
--------------------------
Python 2.6.1 (r261:67515, Feb 12 2009, 16:56:09)
[GCC 4.0.1 (Apple Computer, Inc. build 5367)] on darwin
Type "help", "copyright", "credits" or "license" for more information.
>>> import numpy
>>> z=numpy.longdouble((1,2,3))
>>> z
array([1.0, 2.0, 3.0], dtype=float128)
>>> y=numpy.double((4,5,6))
>>> y
array([ 4., 5., 6.])
>>> numpy.dot(y,y)
77.0
>>> numpy.dot(y,z)
32.0
>>> numpy.dot(z,z)
14.0
---------------------------
I'm prepared to make the changes to numpy. It looks like they would
entail modifying the source file "multiarrayobject.c". The existing
code is:
/*NUMPY_API
Numeric.innerproduct(a,v)
*/
static PyObject *
PyArray_InnerProduct(PyObject *op1, PyObject *op2)
{
If I make changes, can I get them inserted into the numpy source base
(after appropriate testing??)
Thanks in advance,
David Henderson
On Feb 14, 2009, at 8:04 AM, numpy-discussion-request@scipy.org wrote:
> Fri, 13 Feb 2009 12:04:12 -0800, David Henderson wrote:
>> I'd like accumulate the summation in extended precision, "double" sum
>> for float inputs, "long double" sum for "double" inputs.
>>
>> A way of doing this is to add an optional third argument to dot - to
>> specify the summation type.
>
> `dot` does matrix-matrix, matrix-vector, or vector-vector products.
> These
> are usually implemented via calling the BLAS linear algebra
> library, and
> AFAIK BLAS only has routines for type-homogenous arguments. So I doubt
> this can be implemented for `dot` in any way.
>
> On the other hand, `numpy.sum` already has a `dtype` argument that
> specifies the accumulator data type. Maybe you can use that?
>
> --
> Pauli Virtanen
> -----------------------------
>
More information about the Numpy-discussion
mailing list