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 LAA25471; Fri, 23 Jan 2004 11:25:18 +0100 (MET) X-Authentication-Warning: pauillac.inria.fr: majordomo set sender to owner-caml-list@pauillac.inria.fr using -f Received: from concorde.inria.fr (concorde.inria.fr [192.93.2.39]) by pauillac.inria.fr (8.7.6/8.7.3) with ESMTP id LAA25765 for ; Fri, 23 Jan 2004 11:25:17 +0100 (MET) Received: from relay.felk.cvut.cz (relay.felk.cvut.cz [147.32.80.7]) by concorde.inria.fr (8.11.1/8.11.1) with ESMTP id i0NAPHP29171 for ; Fri, 23 Jan 2004 11:25:17 +0100 (MET) Received: from k333.felk.cvut.cz (k333.felk.cvut.cz [147.32.87.5]) by relay.felk.cvut.cz (8.12.10/8.12.10) with ESMTP id i0NAP6Qo089518 for ; Fri, 23 Jan 2004 11:25:06 +0100 (CET) (envelope-from kybic@fel.cvut.cz) Received: from K333/SpoolDir by k333.felk.cvut.cz (Mercury 1.48); 23 Jan 04 11:25:06 +0100 Received: from SpoolDir by K333 (Mercury 1.48); 23 Jan 04 11:24:38 +0100 Received: from cmpgw-iii-23.felk.cvut.cz (147.32.84.2) by k333.felk.cvut.cz (Mercury 1.48) with ESMTP; 23 Jan 04 11:24:32 +0100 To: caml-list@inria.fr Subject: [Caml-list] Q: Efficient operations on complex numbers From: Jan Kybic Date: 23 Jan 2004 11:24:32 +0100 Message-ID: User-Agent: Gnus/5.09 (Gnus v5.9.0) Emacs/21.2.93 MIME-Version: 1.0 Content-Type: text/plain; charset=us-ascii X-MailScanner-felk: Found to be clean X-MailScanner-SpamCheck-felk: not spam, SpamAssassin (score=-4.9, required 5, BAYES_00) X-Loop: caml-list@inria.fr X-Spam: no; 0.00; boxing:01 avoided:01 ocamlopt:01 mandrake:01 3.06:01 haskell:01 consumed:01 seq:01 implemented:01 ocamlopt:01 calc:01 calc:01 0.34:01 0.0:01 0.0:01 Sender: owner-caml-list@pauillac.inria.fr Precedence: bulk Hello, I would appreciate if you could give me same help on the following topic: Q: How to make efficient computation with complex values in OCaml? Subquestions: Q1: Can the boxing of arguments be avoided? Q2: My ocamlopt (from Linux Mandrake ocaml-3.06-12mdk.i586.rpm) does not like the '-ffast-math' switch. Why? Q2: Is OCaml the right language for me, is there any other you would recommend for me to consider? Background ========== I want to develop a fairly complex numerical algorithm involving hierarchical Boundary Element Method for the electromagnetic problem, i.e. working with surface meshes, huge data, complicated data structures, and computational kernel using complex numbers. My first attempt at coding it was in C++ but I did not achieve to program the full version of the algorithm, since I was simply getting lost in details of the implementation, handling of the data etc. My second attempt was to recode the algorithm in Haskell using ghc. This went fine and I was quite enjoying it. However, when I started to test the program on realistic size problems. I found that unwanted laziness kept creeping in that I could neither detect nor control. The program consumed a lot of memory and was quite slow, even though the computational kernel was fast enough. I improved things a little by adding `seq` all over the place but I gave up in the end. My third attempt was to code again in C++, this time using a third-party library. It went fine while I was doing simple things but then what I could not make the library do what I needed and with the limited documentation I had, I was not able to extend it. So here I am, trying to do the whole thing again, this time in OCaml, because I hope it would allow me to design the data structures I will need and at the same time the resulting algorithm should be sufficiently efficiently compiled, run fast and not used too much memory. My problem ========== So far, I have implemented a part of the computational kernel in OCaml that calculates approximations using spherical harmonics. It uses complex values. The Ocaml version is basically is straightforward rewrite of the C++ version and is almost exactly twice as slow as the C++ version. This is intriguing because the parts code that only uses real arithmetic is about as fast as the corresponding C++ version. So why is working with complex numbers so expensive in OCaml and what can be done about it? Experiments =========== 1) I compiled the following trivial function (using ocamlopt on a Linux box): let calc_sqrt _ = let rec calc_sqrt' x i = if i<=0 then x else calc_sqrt' (0.5 *. ( x +. 2.0 /. x )) (i-1) in calc_sqrt' 2.0 10000000 It takes 0.4s while the corresponding C++ version double calc_sqrt() { double x=2.0 ; for (int i=0;i<10000000;i++) x=0.5 * (x+2.0/x) ; return x ; } takes 0.34s. That's fair enough, very result good for OCaml indeed. 2) I did the same thing, only with complex numbers from the Complex module: let calc_sqrt_c _ = let two = { Complex.re=2.0 ; Complex.im=0.0 } and half = { Complex.re=0.5 ; Complex.im=0.0 } in let rec calc_sqrt_c' x i = if i<=0 then x else calc_sqrt_c' (Complex.mul half ( Complex.add x (Complex.div two x ))) (i-1) in calc_sqrt_c' two 10000000 and in C++: typedef std::complex complex_type ; complex_type calc_sqrt_c() { complex_type two(2.0,0), half(0.5,0) ; complex_type x=two ; for (int i=0;i<10000000;i++) x=half * (x+two/x) ; return x ; } Here however, the OCaml version takes 1.25s against 0.56s in C++, i.e. a 2-fold slowdown. 3) I could improve the OCaml time to 0.98s by including the code from the Complex module in my source file. (Is there a way of achieving the same result without actually copying it?) On the other hand, it did not matter whether I used tuples or records to represent the complex numbers, whether I defined the 'add' and 'mul' functions locally or at the top level, and whether I asked the compiler to inline or not. Speculation =========== I deduce, that the main factor contributing to the slowdown, is the boxing and unboxing of values. Can this be optimized away? Would MLton do it? Conclusions =========== I can live with 2-fold slowdown with respect to C++ if this permits me to actually write the program. However, it seems to be a shame that OCaml can generate code on par with C++ in some cases (real arithmetic) but not in others. So if there is a way to do it, I would like to know. Thank you, Jan -- ------------------------------------------------------------------------- Jan Kybic tel. +420 2 2435 7264 or , http://cmp.felk.cvut.cz/~kybic ------------------- To unsubscribe, mail caml-list-request@inria.fr Archives: http://caml.inria.fr Bug reports: http://caml.inria.fr/bin/caml-bugs FAQ: http://caml.inria.fr/FAQ/ Beginner's list: http://groups.yahoo.com/group/ocaml_beginners