[erlang-questions] Calculate PI in Erlang faster than in Matlab

[Zvi]

Hi All,

I wrote "Migrating to Multi-Core" presentation, where I describe various
parallel programming models and APIs. As a practical example, I calculate PI
number using various parallel APIs, such as pthreads, Win32 threads, OpenMP,
MPI, Intel TBB and off course plain old serial C code.
So I decided to add example in my favorite parallel language - Erlang. The
code and math description is attached.

My machine hase dual-core CPU, but I have older version of Matlab, which
doesn't use all cores, so Matlab version took 125 ms:

>> tic; N=1000000; Step = 1/N; PI = Step*sum(4./(1+(((1:N)-0.5)*Step).^2));
>> toc
Elapsed time is 0.125428 seconds.

Below is Erlang results:

2> c(pi).
{ok,pi}
3> pi:test().
[{math,{1,3.14159}},
 {list_creation,{77998,3.14159}},
 {lc,{671999,3.14159}},
 {lc_no_pow,{593999,3.14159}},
 {map,{436999,3.14159}},
 {mapfoldl,{687999,3.14159}},
 {serial_decr,{170999,3.14159}},
 {serial,{171999,3.14159}},
 {parallel,{93999,3.14159}}]
4>

another run:

6> pi:test().
[{math,{1,3.14159}},
 {list_creation,{62998,3.14159}},
 {lc,{671999,3.14159}},
 {lc_no_pow,{671999,3.14159}},
 {map,{436999,3.14159}},
 {mapfoldl,{687999,3.14159}},
 {serial_decr,{170999,3.14159}},
 {serial,{171999,3.14159}},
 {parallel,{77999,3.14159}}]

NOTE: "lc" - is short for "list comprehension".
 
As you can see Erlang parallel version beats Matlab:  78-93ms vs. 123-132ms.

The results are impressive, considering, thats only list creation (i.e
list:seq(1,N) ) takes 63-78ms.

If anyone can improve parallel version not only in speed, but also in making
it clear, so it will be fitting presentation / textbook style.
Also, I want to modify HOF samples to use plists module.

If anyone can give me pointers, how modify this example to use big integer
arithmetics, i.e. to calculate first N digits of PI.

Thanks in adavnce,
Zvi

==================================

-module(pi).
-author("Zvi").
-compile([export_all]).

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% functional and parallel PI calculation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
test()->
	N = 1000000,
	Impls =
[math,list_creation,lc,lc_no_pow,map,mapfoldl,serial_decr,serial,parallel],
	[ {Impl,timer:tc(pi,calc_pi,[Impl,N])} || Impl<-Impls ].

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 1. built-in constant
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
calc_pi(math,_N) -> 
	math:pi();
	
calc_pi(list_creation,N) ->
	lists:seq(1,N),
	3.14159;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 2. serial - list comprehension - direct Matlab translation
%%	Step = 1/N;
%%	PI = Step*sum(4./(1+(((1:N)-0.5)*Step).^2))
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
calc_pi(lc,N) ->
	Step = 1/N,
	Step*lists:sum( [4.0/(1.0 + math:pow(Step*(I-0.5),2)) || I<-lists:seq(1,N)]
);

calc_pi(lc_no_pow,N) ->
	Step = 1/N,
	Xs = [Step*(I-0.5) || I<-lists:seq(1,N)],
	Step*lists:sum( [4.0/(1.0 + X*X) || X<-Xs] );

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 3. serial - using HOF: map
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
calc_pi(map,N) ->
	Step = 1/N,
	F = fun(I) -> 
		X=Step*(I-0.5),
		4.0/(1.0+X*X)
	    end,
	Step*lists:sum( lists:map( F, lists:seq(1,N) ) );

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 4. serial - using HOF: mapfoldl
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
calc_pi(mapfoldl,N) ->
	Step = 1/N,
	F = fun(I,Sum) -> 
		X=Step*(I-0.5), 
		V=4.0/(1.0+X*X), 
                 {0,Sum+V}
	    end,
	{_,Sum} = lists:mapfoldl(F, 0, lists:seq(1,N)),
	Step*Sum;

%% replace list to plist module  in tests using HOFs
	
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 5. serial - tail recursion - decrement index
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
calc_pi(serial_decr,N) ->
	calc_pi1(N,1/N,0);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 6. serial - tail recursion - increment index
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
calc_pi(serial,N) ->
	calc_pi2(1,N,1/N);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 7. parallel SMP version - process per scheduler
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
calc_pi(parallel,N) ->
	NWorkers = erlang:system_info(schedulers),
	NPerWorker = trunc(N/NWorkers+0.5),
	Step = 1/N,
	Self = self(),
	Pids = [spawn(fun() -> Self !
{sum,calc_pi2(I,lists:min([N,I+NPerWorker-1]),Step)} end) 
	        || I<-lists:seq(1,N,NPerWorker) ],
	lists:sum([receive {sum,S} -> S end || _<-Pids]).

	
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%	
calc_pi1(0,Step,Sum) -> Step*Sum;
calc_pi1(I,Step,Sum) ->
	X=Step*(I-0.5),
	V=4.0/(1.0+X*X),
	calc_pi1(I-1,Step,Sum+V).

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
calc_pi2(From,To,Step) ->
	calc_pi2(From,From,To,Step,0).

calc_pi2(To,_From,To,Step,Sum) -> Step*Sum;
calc_pi2( I,From,To,Step,Sum) ->
	X=Step*(I-0.5),
	V=4.0/(1.0+X*X),
	calc_pi2(I+1,From,To,Step,Sum+V).
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

http://www.nabble.com/file/p16440056/pi.GIF pi.GIF 
http://www.nabble.com/file/p16440056/pi.erl pi.erl 
-- 
View this message in context: http://www.nabble.com/Calculate-PI-in-Erlang-faster-than-in-Matlab-tp16440056p16440056.html
Sent from the Erlang Questions mailing list archive at Nabble.com.