[Haskell-cafe] Re: [Haskell] Linear Equation Solver Using Arrays

Xiao-Yong Jin xj2106 at columbia.edu
Mon Sep 17 00:57:18 EDT 2007

"Andrew Wagner" <wagner.andrew at gmail.com> writes:

> (Replying on haskell-cafe)
> Let me start with a disclaimer: I haven't looked at your code
> extensively. That said, the feeling that I get from it is one of
> listening to a non-native speaker. The things in your program are
> obviously not completely inaccurate, or they wouldn't have
> type-checked. And, presumably, you've already done some testing, so it
> probably produces the right answers. Similarly, there are many things

I am ashamed.  I forgot to put an `abs' when extracting
implicit scaling factor, which would result in a failure
upon a matrix with negative elements.  I'm attaching here a
corrected version of it.  I also changed a little bit other
things, such as replacing some recursive calls with `forM_',
and using `>>=' to "sugar" it.
[2. MatrixOp, amended version --- application/octet-stream; MatrixOp.hs.gz]...

> in English that you _could_ say, which would be entirely
> grammatically correct, but which native speakers would
> never say. For example, in your code, I noticed a pattern
> which I've just recently learned: <code> do ...  tmp <-
> readArray scaleInfo j writeArray scaleInfo imax tmp
> </code>
> Now, if you were to de-sugar the last 2 lines, it would look something
> like this:
> <code>
>   readArray scaleInfo j >>= \tmp -> writeArray scaleInfo imax tmp
> </code>
> ...but as soon as you see that lambda expression, you should
> (hopefully) say "oh, that can be written more easily:
> <code>
>   readArray scaleInfo j >>= writeArray scaleInfo imax
> </code>
> Of course, a native speaker will probably not even think about the
> intermediate lambda expression, but the point is, this is an idiom.
> And, really, this is a very tiny example. Your entire program is
> written in a very imperative style. "Do this, then do that, then do
> the other thing." And this is to be expected. After all, you said that
> you were directly porting the algorithm from an imperative language.
> However, if you're going to do this idiomatically, you really need to
> think very differently.

That's the problem.  Once start using arrays, I can't stop
thinking it with indices.  Probably that's the reason I'm
using Arrays, in the first place.  I couldn't think of a way
to express the algorithm functionally, or in a Haskell way.
I was trying to write it using lists, but that failed,
because I couldn't manage so many modifications of lists and
keeping track of all the states.  And that's why I used
STRef in the function `lubksb'.  It's really hard to think 
of a functional alternative of some code in C, like

    sum += a[i][j] * b[j][k]

Possibly a `zip' or a `fold' on some lists would do in that
case.  But encountering something more complicated like the
function I wrote, I really got frustrated to think of a
idiomatic, native speakers' way of doing things.

> I suggest reading and pondering this well-known blog entry
>for more ideas about the differences between the approaches
>of functional and imperative programming:

Thanks for the link.

> Sorry, I don't have more specific suggestions, but I
>suspect this is what most people thought when they read
>your code: it's hard to make specific suggestions, when
>you're really approaching the problem in a non-idiomatic

I would appreciate it if someone could point me to certain
guide on doing numerical calculations using Haskell.

> On 9/16/07, Xiao-Yong Jin <xj2106 at columbia.edu> wrote:
>> Dear Haskellers,
>> My thanks to people on the irc channel, especially `int-e'.
>> With their help, I managed to write a linear equation solver
>> using STUArrays.
>> It is basically a direct rewrite of a C version, which
>> adopts Crout's LU decomposition method with implicit
>> pivoting, in chapter 2.3 of the book Numerical Recipes.
>> I must admit that the code is in a really bad style,
>> somewhat ugly, because of the limit of my Haskell knowledge.
>> The code is attached in this email for your amusement.  I
>> would be really thankful if you can give me any kind of
>> comments, idea, improvement, criticism, because I sincerely
>> want to make the code better.
>> Best regards,
>> Xiao-Yong
>> --
>>     c/*    __o/*
>>     <\     * (__
>>     */\      <
>> _______________________________________________
>> Haskell mailing list
>> Haskell at haskell.org
>> http://www.haskell.org/mailman/listinfo/haskell

    c/*    __o/*
    <\     * (__
    */\      <

More information about the Haskell-Cafe mailing list