[Haskell-cafe] BLAS Solve Example

[Darrin Thompson]

I'm stuck on something that I thought would be easy.

I have a matrix and a vector.

> module Main where
> import Data.Vector.Dense
> import Data.Matrix.Dense
> import BLAS.Matrix.Solve
>
> m = listMatrix (2, 3) ([1, 2, 3, 4, 5, 6]::[Double])
> v = listVector 2 ([1, 2]::[Double])
>
> main = do ???

Can I use one of the haskell libraries to solve that system? My goal
is to automate some tedium on small exercises from the Linear Algebra
text by Jim Heffron.

--
Darrin
_______________________________________________
Haskell-Cafe mailing list
Haskel...@haskell.org
http://www.haskell.org/mailman/listinfo/haskell-cafe

[Patrick Perry]

Sorry Darrin, the BLAS library only includes matrix multiplication and
solving triangular systems. To solve a general system, you would need
to use LAPACK, but there aren't any bindings for that library yet. I
would suggest you take a look at the hmatrix package, which includes a
lot more linear algebra than blas does..

Patrick

[Alberto Ruiz]

Darrin Thompson wrote:
> I'm stuck on something that I thought would be easy.
>
> I have a matrix and a vector.
>
>> module Main where
>> import Data.Vector.Dense
>> import Data.Matrix.Dense
>> import BLAS.Matrix.Solve
>>
>> m = listMatrix (2, 3) ([1, 2, 3, 4, 5, 6]::[Double])
>> v = listVector 2 ([1, 2]::[Double])
>>
>> main = do ???
>
> Can I use one of the haskell libraries to solve that system? My goal
> is to automate some tedium on small exercises from the Linear Algebra
> text by Jim Heffron.
>
> --
> Darrin

Using hmatrix:

> import Numeric.LinearAlgebra
>
> m = (2><3) [1, 2, 3,
> 4, 5, 6]
>
> v = 2 |> [1, 2 :: Double]
>
> sol = m <\> v

$ ghci solve.hs
*Main> sol
3 |> [-5.555555555555511e-2,0.11111111111111113,0.2777777777777776]

Best regards,

Alberto

[apfelmus]

Darrin Thompson wrote:

> On Wed, Jul 23, 2008 at 2:12 AM, Alberto Ruiz <ar...@um.es> wrote:
>> $ ghci solve.hs
>> *Main> sol
>> 3 |> [-5.555555555555511e-2,0.11111111111111113,0.2777777777777776]
>>
>

> I was hoping for rational solutions. If I were a true jedi master I'd
> write my own solver, which might be the right thing to do. All I know
> so far is gauss' method. Probably I'd learn something implementing the
> back substitution. hmm....

The Hugs interpreter has a few nice Haskell demos, including a simple gauss
elimination

http://darcs.haskell.org/hugs98/demos/Matrix.hs

Of course, this is unsuitable for serious floating point calculations.

Regards,
apfelmus

[Dan Weston]

A jedi master might stick with the existing double precision solver,
then convert the results to best rational approximation [1], then do a
forward solve on the rational versions of matrices, adjusting numerator
and denominator to eliminate the residual error (with a heuristic to
favor common factors). If you are very lucky, such a rational number
will exist, depending on your limits of humongous.

[1] e.g. http://www.dtashley.com/howtos/2007/01/best_rational_approximation/

Darrin Thompson wrote:
> On Wed, Jul 23, 2008 at 2:12 AM, Alberto Ruiz <ar...@um.es> wrote:

>> $ ghci solve.hs
>> *Main> sol
>> 3 |> [-5.555555555555511e-2,0.11111111111111113,0.2777777777777776]
>>
>

> I was hoping for rational solutions. If I were a true jedi master I'd
> write my own solver, which might be the right thing to do. All I know
> so far is gauss' method. Probably I'd learn something implementing the
> back substitution. hmm....
>

> Thanks.
>
> --
> Darrin