Parallel Ant Colony optimization

[Jon Harrop]

Eric Rollins posted an interesting parallel numerical benchmark with
implementations in Erlang, Haskell and MLton SML:

http://proliferationofniches.blogspot.com/2008/06/functional-languages-on-quad-core_03.html

I translated the programs to OCaml and obtained the following results for
1000 iterations on a dual core 4400+ Athlon64:

OCaml: 0.49s (parallel)
MLton: 2.9s (serial)
GHC 6.8.2: 52s (parallel)

I am interested in benchmarking parallel programs written in various
functional languages but, with the benefit of hindsight, this is not a
great benchmark of parallelism because it is embarrassingly parallel.
However, it raises an interesting question of what to do with PRNGs in the
context of parallelism. Specifically, whether or not to skip the
synchronization that is only required to make a random algorithm
deterministic and reproducible.

Anyway, here's my code:

open Printf

let arg i d = try int_of_string Sys.argv.(i) with _ -> d

let seed = arg 1 1
let boost = arg 2 5
let iter = arg 3 100
let n = arg 4 200
let decr = float boost /. float iter

module Set = struct
type t = {mutable count: int; array: bool array}

let make n = {count = 0; array = Array.make n false}

let add m x =
if not(m.array.(x)) then begin
m.count <- m.count + 1;
m.array.(x) <- true
end

let mem m x = m.array.(x)

let count m = m.count
end

module Path = struct
let do_sum_weight used current cities pher =
let sum = ref 0.0 in
for c=0 to n-1 do
if not(Set.mem used c) then
sum := !sum +. cities.(current).(c) *. (1.0 +. pher.(current).(c))
done;
!sum +. 0.0

let find_sum_weight used current cities pher sought =
let next = ref 0 in
let sum = ref 0.0 in
try
for c=0 to n-1 do
if not(Set.mem used c) then begin
if !sum >= sought then raise Exit;
next := c;
sum := !sum +. cities.(current).(c) *. (1.0 +. pher.(current).(c))
end
done;
!next
with Exit -> !next

let make cities pher rng =
let used = Set.make n in
let start = Random.State.int rng n in
let rec lp(path, current) =
if Set.count used = n then path else
let sum = do_sum_weight used current cities pher in
let next = find_sum_weight used current cities pher
(Random.State.float rng sum) in
Set.add used next;
lp(next::path, next) in
Set.add used start;
lp([start], start)

let rec iter_aux first f = function
| [] -> f first first
| [last] -> f first last
| c1::(c2::_ as t) -> f c1 c2; iter_aux first f t

let iter f = function
| [] -> ()
| first::_ as t -> iter_aux first f t

let length cities path =
let len = ref 0.0 in
iter (fun r c -> len := !len +. cities.(r).(c)) path;
!len
end

module Pher = struct
let make n =
Array.make_matrix n n 0.0

let update pher path =
Path.iter (fun r c -> pher.(r).(c) <- pher.(r).(c) +. float boost) path

let evaporate pher =
for i=0 to n-1 do
for j=0 to n-1 do
let x = pher.(i).(j) in
pher.(i).(j) <- if x <= decr then 0. else x -. decr
done
done
end

let maxd = 100.0

let cities =
let rng = Random.State.make [|seed|] in
Array.init n (fun _ -> Array.init n (fun _ -> Random.State.float rng
(maxd -. 1.0) +. 1.0))

let invoke (f : 'a -> 'b) x : unit -> 'b =
let input, output = Unix.pipe() in
match Unix.fork() with
| -1 -> (let v = f x in fun () -> v)
| 0 ->
Unix.close input;
let output = Unix.out_channel_of_descr output in
Marshal.to_channel output (try `Res(f x) with e -> `Exn e) [];
close_out output;
exit 0
| pid ->
Unix.close output;
let input = Unix.in_channel_of_descr input in
fun () ->
let v = Marshal.from_channel input in
ignore (Unix.waitpid [] pid);
close_in input;
match v with
| `Res x -> x
| `Exn e -> raise e

let work seed =
let pher = Pher.make n in
let rng = Random.State.make [|seed|] in
let path = ref (Path.make cities pher rng) in
let len = ref (Path.length cities !path) in
for i=1 to iter do
let path' = Path.make cities pher rng in
let len' = Path.length cities path' in
if len' <= !len then begin
printf "%d %f\n" i len';
Pher.update pher path';
len := len';
path := path';
end;
Pher.evaporate pher;
done;
!len, !path

let () =
printf "seed=%d boost=%d iter=%d n=%d\n%!" seed boost iter n;
let futures = List.map (invoke work) [1; 2] in
let results = List.map (fun f -> f()) futures in
match List.sort compare results with
| (len, path)::_ ->
List.iter (printf "%d ") path;
printf "= %f\n" len
| _ -> ()

--
Dr Jon D Harrop, Flying Frog Consultancy
http://www.ffconsultancy.com/products/?u