From mboxrd@z Thu Jan  1 00:00:00 1970
Received: (from majordomo@localhost) by pauillac.inria.fr (8.7.6/8.7.3) id EAA19572; Wed, 28 Mar 2001 04:19:05 +0200 (MET DST)
X-Authentication-Warning: pauillac.inria.fr: majordomo set sender to owner-caml-list@pauillac.inria.fr using -f
Received: from nez-perce.inria.fr (nez-perce.inria.fr [192.93.2.78]) by pauillac.inria.fr (8.7.6/8.7.3) with ESMTP id EAA19568 for <caml-list@pauillac.inria.fr>; Wed, 28 Mar 2001 04:19:04 +0200 (MET DST)
Received: from smtp10.atl.mindspring.net (smtp10.atl.mindspring.net [207.69.200.246])
	by nez-perce.inria.fr (8.11.1/8.10.0) with ESMTP id f2S2J3T14726
	for <caml-list@inria.fr>; Wed, 28 Mar 2001 04:19:03 +0200 (MET DST)
Received: from dylan (1Cust65.tnt2.tucson.az.da.uu.net [63.11.142.65])
	by smtp10.atl.mindspring.net (8.9.3/8.8.5) with SMTP id VAA14335
	for <caml-list@inria.fr>; Tue, 27 Mar 2001 21:18:59 -0500 (EST)
Message-ID: <001801c0b72d$7fd17a60$210148bf@dylan>
From: "David McClain" <dmcclain1@mindspring.com>
To: <caml-list@inria.fr>
Subject: [Caml-list] Complex Numbers
Date: Tue, 27 Mar 2001 19:19:15 -0700
MIME-Version: 1.0
Content-Type: text/plain;
	charset="iso-8859-1"
Content-Transfer-Encoding: 7bit
X-Priority: 3
X-MSMail-Priority: Normal
X-Mailer: Microsoft Outlook Express 5.00.2919.6600
X-MimeOLE: Produced By Microsoft MimeOLE V5.00.2919.6600
Sender: owner-caml-list@pauillac.inria.fr
Precedence: bulk

"One of the articles by William Kahan (www.cs.berkeley.edu/~wkahan/)
"How JAVA's Floating-Point Hurts Everyone Everywhere" discusses
(among other things) some traps in complex number implementations,
including those in Fortran.
"

My NML, written in OCaml, paid particular attention to these issues after I
read everything I could from Kahan. To my knowledge, aside from Common Lisp
(perhaps), it is the only proper implementation of complex arithmetic
available.

In particular, Kahan's demo of "Borda's Pipe", shows that stream lines are
incorrectly computed in every language I have tried: Mathematica, RSI/IDL,
Fortran, C/C++, Basic, with the possible exception of Common Lisp. Only NML
has passed demonstrably passed this test.

The difficulty stems from the fact that, to paraphrase Kahan, "a tuple
representation of complex numbers is insufficient, given the Riemann
surfaces of common transcendental functions". There are two zeros defined in
the IEEE Floating Point Standard: 0+, and 0-.  Arithmetic in the floating
point domain do not obey classical algrebraic relations with respect to
conventional identity elements.

I sweated out two solid weeks in gaining a full understanding of his rather
cryptic statement. NML represents my final victory over this issue and I can
thank OCaml and the INRIA team for their wonderful implementation language
that made it all so readily possible.

- DM


Borda NML code follows:

--------------------------------------------------
(* borda.nml -- test of complex arithmetic in NML *)
(* DM/MCFA  09/99 *)

let g z = z^2 + z*sqrt(z^2+1)
let f z = 1 + g(z) + log(g z)

let ls =
  let ptor r theta =
    complex (r*cos theta) (r*sin theta)
  in
  let t = tan(pi/2*[0, 0.001, .. 1]) in
    map (fun ang ->
    ptor t ang)
      (lis (pi * [-0.5, -0.45, .. 0.5]))

let rec pshow (l, ... as args) =
  let clip x =
    x #< 400 #> -400
  in
  let z = f(l) in
    oplot(clip(re z), clip(im z) | args @ [clip: true])

let borda() =
  let s = 5 in
    axes(xrange: [-s,s], yrange: [-s,s],
 title: "The Correct Borda!");
    pshow(hd ls, color: red, thick: 2);
    app pshow (tl ls)

let _ = borda()



-------------------
To unsubscribe, mail caml-list-request@inria.fr.  Archives: http://caml.inria.fr