[Numpy-discussion] type of matrixmultiply result

[Perry Greenfield]

Rick White wrote:

> On 24 Jun 2004, Todd Miller wrote:
> 
> > On Thu, 2004-06-24 at 06:14, Curzio Basso wrote:
> >
> > > I noticed that when multiplying two matrices of type Float32, 
> the result
> > > is Float64:
> >
> > I modified dot() and innerproduct() this morning to return Float32 and
> > Complex32 for like inputs.
> 
> I wonder whether it would be worth providing an option to accumulate
> the sums using Float64 and to convert to Float32 before storing them in
> an array.  I suspect that one reason this returned Float64 is that it
> is very easy to run into precision/roundoff problems in
> single-precision matrix multiplies.  You could avoid that by using
> doubles for the sum while still returning the result as a single.
> 				Rick
> 
I definitely agree. I'm pretty certain the reason it was done
with double precision floats is the sensitivity to roundoff
issues with matrix operations. I think Rick is right though that
only intermediate calculations need to be done in double precision
and that doesn't require the whole output array to be kept that way.

Perry