search

haskell-numerics / hmatrix

Linear algebra and numerical computation
381 stars 104 forks source link

Matrix multiplication produces wrong result #357

Closed Javran closed 2 years ago

Javran commented 2 years ago

hmatrix version 0.20.1 with ghc 8.8.4:

> :set -XFlexibleContexts 
> import Numeric.LinearAlgebra
> let x = ((1 >< 8) (repeat 0) === ident 8) ||| (9 >< 1) ([1, 0, 1] Prelude.<> repeat 0) :: Matrix Z
> x
(9><9)
 [ 0, 0, 0, 0, 0, 0, 0, 0, 1
 , 1, 0, 0, 0, 0, 0, 0, 0, 0
 , 0, 1, 0, 0, 0, 0, 0, 0, 1
 , 0, 0, 1, 0, 0, 0, 0, 0, 0
 , 0, 0, 0, 1, 0, 0, 0, 0, 0
 , 0, 0, 0, 0, 1, 0, 0, 0, 0
 , 0, 0, 0, 0, 0, 1, 0, 0, 0
 , 0, 0, 0, 0, 0, 0, 1, 0, 0
 , 0, 0, 0, 0, 0, 0, 0, 1, 0 ]
> x * x
(9><9)
 [ 0, 0, 0, 0, 0, 0, 0, 0, 1
 , 1, 0, 0, 0, 0, 0, 0, 0, 0
 , 0, 1, 0, 0, 0, 0, 0, 0, 1
 , 0, 0, 1, 0, 0, 0, 0, 0, 0
 , 0, 0, 0, 1, 0, 0, 0, 0, 0
 , 0, 0, 0, 0, 1, 0, 0, 0, 0
 , 0, 0, 0, 0, 0, 1, 0, 0, 0
 , 0, 0, 0, 0, 0, 0, 1, 0, 0
 , 0, 0, 0, 0, 0, 0, 0, 1, 0 ]
> x ^ 80
(9><9)
 [ 0, 0, 0, 0, 0, 0, 0, 0, 1
 , 1, 0, 0, 0, 0, 0, 0, 0, 0
 , 0, 1, 0, 0, 0, 0, 0, 0, 1
 , 0, 0, 1, 0, 0, 0, 0, 0, 0
 , 0, 0, 0, 1, 0, 0, 0, 0, 0
 , 0, 0, 0, 0, 1, 0, 0, 0, 0
 , 0, 0, 0, 0, 0, 1, 0, 0, 0
 , 0, 0, 0, 0, 0, 0, 1, 0, 0
 , 0, 0, 0, 0, 0, 0, 0, 1, 0 ]

Meanwhile same computation in octave:

GNU Octave, version 6.3.0
Copyright (C) 2021 The Octave Project Developers.
This is free software; see the source code for copying conditions.
There is ABSOLUTELY NO WARRANTY; not even for MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE.  For details, type 'warranty'.

Octave was configured for "x86_64-pc-linux-gnu".

Additional information about Octave is available at https://www.octave.org.

Please contribute if you find this software useful.
For more information, visit https://www.octave.org/get-involved.html

Read https://www.octave.org/bugs.html to learn how to submit bug reports.
For information about changes from previous versions, type 'news'.

octave:1> x = [[zeros(1,8); eye(8)] [1 0 1 0 0 0 0 0 0]']
x =

   0   0   0   0   0   0   0   0   1
   1   0   0   0   0   0   0   0   0
   0   1   0   0   0   0   0   0   1
   0   0   1   0   0   0   0   0   0
   0   0   0   1   0   0   0   0   0
   0   0   0   0   1   0   0   0   0
   0   0   0   0   0   1   0   0   0
   0   0   0   0   0   0   1   0   0
   0   0   0   0   0   0   0   1   0

octave:2> x * x
ans =

   0   0   0   0   0   0   0   1   0
   0   0   0   0   0   0   0   0   1
   1   0   0   0   0   0   0   1   0
   0   1   0   0   0   0   0   0   1
   0   0   1   0   0   0   0   0   0
   0   0   0   1   0   0   0   0   0
   0   0   0   0   1   0   0   0   0
   0   0   0   0   0   1   0   0   0
   0   0   0   0   0   0   1   0   0

octave:3> x ^ 80
ans =

   126    11   126    45    84   120    37   210    20
     9   126    11   126    45    84   120    37   210
   210    20   252    56   210   165   121   330    57
    37   210    20   252    56   210   165   121   330
   120    37   210    20   252    56   210   165   121
    84   120    37   210    20   252    56   210   165
    45    84   120    37   210    20   252    56   210
   126    45    84   120    37   210    20   252    56
    11   126    45    84   120    37   210    20   252
Javran commented 2 years ago

My aplogies - turns out * are elementwise and I need <> for matrix multiplication.

idontgetoutmuch commented 2 years ago

No problem - thanks for interest :-)