Accueil     À propos     Téléchargement     Ressources     Contactez-nous

Ce site est rarement mis à jour. Pour les informations les plus récentes, rendez-vous sur le nouveau site OCaml à l'adresse ocaml.org.

Browse thread
[Caml-list] Complex numbers in OCaml
[ Home ] [ Index: by date | by threads ]
[ Search: ]

[ Message by date: previous | next ] [ Message in thread: previous | next ] [ Thread: previous | next ]
 Date: 2001-03-28 (01:04) From: Michael Hohn Subject: Re: [Caml-list] Complex numbers in OCaml
```
>> ...
>>
>> I need complex numbers. I tried it this way:
>>
>> ...

I also needed complex numbers a while ago.
There are several options:
o   Make a FORTRAN binding.
This is tedious, and for simple low-level operations (add, etc.)
will probably not help performance.  If

o   Implement them in ocaml.
This is very simple when using user-defined operators.
Also, it is fast enough for all simple applications

>> ...
>> This works but is quite cumbersome to use. Is there any other way to do it
>> ...

Below is the module I used, based on a struct-like product type instead
of a sum type.

>> ...
>> and how can arrays of complex numbers be implemented efficiently in OCaml?
>> ...

Arrays of  (C) doubles are optimized by ocamlopt, but for complex
types, which are just structs, I do not know.

For the ultimate in efficiency of complex arrays, there is always
FORTRASH :(

Cheers,
Michael

(* \$Id: cmplx.ml,v 1.6.2.1 2000/04/05 16:29:57 hohn Exp \$
Elementary complex number support functions.
Taken from libf2c/libF77.
*)

type cmplx = { re : float; im : float}

let cx re im = { re = re ; im = im}

let zarg z = atan2 z.im z.re

let conj z = { re = z.re ; im = -. z.im}

let (+:) a b = {re =  (a.re +. b.re); im =  (a.im +. b.im) }

let (-:) a b = {re = (a.re -. b.re); im = (a.im -. b.im) }

let ( *: ) a b =
{
re = (  a.re *. b.re -. a.im *. b.im ) ;
im = (a.re *. b.im +. a.im *. b.re )
}

let f__cabs (real : float) (imag : float) =
let real = if (real < 0.0) then (-. real) else real in
let imag = if (imag < 0.0) then (-. imag) else imag in
let (real, imag) = if(imag > real) then (imag, real) else (real, imag) in
if ((real +. imag) = real) then real
else (
let temp = imag /.real in
real *. (sqrt (1.0 +. temp*.temp))
)

let zabs z = {re = f__cabs z.re z.im ; im = 0.0}

let dabs x = abs_float x

let ( /: ) a b =
let abr = if( b.re < 0.0) then (-. b.re) else (b.re)  in
let abi = if( b.im < 0.0) then (-. b.im) else b.im in
if ( abr <= abi ) then
(
if (abi = 0.0) then
raise Division_by_zero
else
let ratio = b.re /. b.im in
let den = b.im *. (1.0 +. ratio*.ratio) in
{
re = (a.re *. ratio +. a.im) /. den ;
im = (a.im *. ratio -. a.re) /. den
}
)
else
(
let ratio = b.im /. b.re in
let den = b.re *. (1.0 +. ratio*.ratio) in
{
re = (a.re +. a.im *. ratio) /. den ;
im = (a.im -. a.re *. ratio) /. den
}
)

let zexp z =
let zi = z.im in
let expx = exp z.re in
{
re = expx *. cos(zi);
im = expx *. sin(zi)
}

let zlog z =
let zi = z.im and zr = z.re in
{
im = atan2 zi zr;
re = (log ( f__cabs zr zi ) )
}

(* Integer power z**i, for small i *)
let zipow z i =
let rv z =
if (i < 0) then ((cx 1.0 0.0) /: z)
else z in
let i = if i < 0 then -i else i in
match i with
| 0 -> cx 1.0 0.0
| 1 -> rv z
| 2 -> rv (z *: z)
| 3 -> rv (z *: z *: z)
| 4 -> rv (z *: z *: z *: z)
| _ ->
let p = ref z in
for k=2 to i do p := !p *: z done ;
rv !p
-------------------
To unsubscribe, mail caml-list-request@inria.fr.  Archives: http://caml.inria.fr

```