[erlang-questions] Erlang and/or MPI

[Richard A. O'Keefe]

On 2/07/2013, at 2:34 AM, David Mercer wrote:

> On Sunday, June 30, 2013, Richard A. O'Keefe wrote:
> 
>> (1) How hard is it to *accurately* transcribe a correct description
>>    of an algorithm into the language in question?
>> 
>>    Here Ada, Fortran, Haskell, and Ada shine.
>>    Smalltalk and Scheme, for all their other virtues, stink.
>>    Erlang is somewhere in the middle.
> 
> Can you expand on this?  I've usually thought most programming languages make it fairly easy to transcribe an algorithm.  But even if not, why would Scheme be particularly bad at it?

Because the syntax is utterly unlike the original presentation.

Let's take the square root of a Quaternion as an example.
Let's hope I don't break it in the process!

Maths:
   sqrt (w,x,y,z) =
      (r/2,x/r,y/r,z/r) where r = sqrt((sqrt(w^2+n^2)+w)*2), if n > 0
      (0,sqrt(-w),0,0), if n = 0 and w < 0
      (sqrt(w),0,0,0),  if n = 0 and w >= 0
      where n = sqrt(x^2+y^2+z^2)

Haskell:
    sqrt q@(Quaternion w x y z)
      | n > 0  = let r = sqrt ((hypot w n + w)*2
                  in Quaternion (r/2) (x/r) (y/r) (z/r)
      | w < 0  = Quaternion 0 (sqrt (negate w)) 0 0
      | w >= 0 = Quaternion (sqrt w) 0 0 0
      where n = vector_part_abs q

This is very close to the mathematical form.  The guards and bodies
are flipped (Miranda had them the other way around).  The other
difference is that it uses hypot w n and vector_part_abs q instead
of sqrt(w^2+n^2) and sqrt(x^2+y^2+z^2) for computational reasons.

Erlang:

    qsqrt(Q = {quaternion,W,X,Y,Z}) ->
        N = vector_part_abs(Q),
        if N > 0 ->
           R = math:sqrt((math:hypot(W, N) + N)*2),
           {quaternion,R/2,X/R,Y/R,Z/R}
         ; W < 0 ->
           {quaternion,0,math:sqrt(-W),0,0}
         ; true ->
           {quaternion,math:sqrt(W),0,0,0}
        end.

Not too bad, actually, except that it doesn't work.
For some reason beyond the ability of mortals to discern,
math:sqrt/1 exists but math:hypot/2 does not.
(It's only been part of standard C for the last 14 years...)
erf/1 and erfc/1 there but hypot/2 not?

It's touch clumsy with repeating the 'math:' prefix so much,
but I don't actually have to do that.  I _could_
-import(math, [sqrt/1,hypot/2]).
if math:hypot/2 existed and then I wouldn't need the prefixes.

C99:
    quaternion qsqrt(quaternion q) {
        double const n = vector_part_abs(q);
        double const w = q.w;

        if (n > 0) {
            double const r = sqrt((hypot(w, n) + w)*2);
            return (quaternion){
                .w = r/2, .x = q.x/r, .y = q.y/r, .z = q.z/r};
        } else      
        if (w < 0) {
            return (quaternion){.w = 0, .x = sqrt(-w), .y = 0, .z = 0};
        } else {
            return (quaternion){.w = sqrt(w), .x = 0, .y = 0, .z = 0};
        }
    }

Still pretty close to the mathematics; a little clunky at making
records and a lot clunky at taking them apart.

Smalltalk:

    sqrt
      |n r|
      n := self vectorPartAbs.
      ^n > 0
         ifTrue:  [
           r := ((w hypot: n) + w) * 2) sqrt.
           Quaternion re: r/2 im: x/r jm: y/r km: z/r]
         ifFalse: [
           w < 0
             ifTrue:  [Quaternion re: 0 im: w negated sqrt jm: 0 km: 0]
             ifFalse: [Quaternion re: w sqrt im: 0 jm: 0 km: 0]]

I think you will agree that converting sqrt(-w) to w negated sqrt
is a little distracting.  What you may not know is that the
operator precedence rules for Smalltalk, well, technically Smalltalk
doesn't have operators, so 1 + 2 * 3 is 9.  Moving unary function
names from the left to the right, and adding extra parentheses
whenever there are two operators in the same expression is just
enough restructuring to make it painfully easy to mess up.

Now let's look at Scheme.
Let's suppose we have a constructor and accessors available already.
R6RS Scheme has define-record-type; R5RS Scheme did not, but there
were plenty of free equivalents.  There are also free add-ons to give
you pattern matching of a sort; I don't know how well they handle
R6RS records though.

Let's also skip over the Erlang-like nonexistence of hypot.

(define (quaternion-sqrt Q)
  (let ((N (vector-part-abs Q))
        (W (quaternion-re Q)))
    (cond
      ((> N 0)
        (let ((R (sqrt (* (+ (hypot W N) W) 2))))
          (make-quaternion (/ R 2)
                           (/ (quaternion-im Q) R)
                           (/ (quaternion-jm Q) R)
                           (/ (quaternion-km Q) R))))
      ((< W 0) (make-quaternion 0 (sqrt (- W)) 0 0))
      (else    (make-quaternion (sqrt W) 0 0 0)))))

This particular example isn't _too_ horrible.

But I invite you to consider

(1) r = sqrt ((hypot w n + w)*2
 vs (R (sqrt (* (+ (hypot W N) W) 2)))

    Transcribing that was not pleasant and I don't have a high
    degree of confidence that it's right.

(2) the fact that the arithmetic operators can be overloaded
    in Haskell, Ada, R, Smalltalk, Eiffel, Fortran,... to apply to
    Quaternions, but they cannot be in Erlang or Scheme.

Yes, I know about the possibility of enclosing a block of code
in
  (let ((outer-plus +) (outer-minus -) (outer-times *) ...)
    (define (+ x y)
      (if (quaternion? x) (quaternion-plus x y)
        (if (quaternion? y) (quaternion-plus y x)
          (outer-plus x y))))
    ...)
but I really don't want to do that, thanks.

The thing is, once Quaternions are properly set up in Haskell
or Fortran, I can write

	abs(cos(x+y) - (cos(x)*cos(y) - sin(x)*sin(y)))

and it will Just Work.  (The reason this is an interesting thing
to compute is that for real and complex numbers, the answer should
be zero, but for Quaternions in general it will not be.)
    
Again, in Haskell and Smalltalk I can set things up so that
	A*V + W
does a vector operation -- and Haskell can let me ensure that it
won't create an intermediate object, thanks to Fusion -- while in
Fortran and R it already does.  This entirely eliminates
wrong-loop-variable and off-by-one errors from the calculation
AND gives me a mathematically appropriate notation.

I don't deny that Scheme is a fine language.
I don't deny that Scheme is perfectly up to the task of
describing any reasonable algorithm.
Nor am I making any claim about the *performance* of any
particular Scheme system, such as Chez (which has some
impressive theory behind its compiler) or Racket.

All I say is that the act of *transcribing* a *numerical*
algorithm from a published paper or technical report to a
computational form is likely to have fewer errors the closer
the two notational systems are.