narray patch

all.rb

require "narray"
T = (RUBY_VERSION<"1.8.0") ? Time : Process
RUBY = File.join(Config::CONFIG['bindir'], Config::CONFIG['ruby_install_name'])

ruby_narray = system( "#{RUBY} dummy.rb" )
python_numeric = system( "python dummy.py numeric" )
python_numarray = system( "python dummy.py numarray" )
python_numpy = system( "python dummy.py numpy" )
octave = system( "octave -qf dummy.m" )

def array_size
  list = [
    100000, 200000, 500000,
    1000000, 2000000, 5000000,
    10000000, 20000000, 50000000
  ]
  mlist = [
    150, 200, 300, 500, 700, 1000, 1500, 2000, 3000
    #300, 400, 700, 1000, 1400, 2000, 3000, 4000, 7000
  ]

  r = 50
  n = nil
  i = nil
  list.each_with_index do |n,i|
    a = NArray.float(n)
    b = NArray.float(n)
    t = bench_time(r) { c = a+b }
    break if t>0.5
  end
  [n, mlist[i], r*2]
end

def bench_time(n)
  t1 = T.times.utime
   for i in 1..n
     yield
   end
  t = T.times.utime - t1
  puts " Time: %.2f sec\n" % [t]
  t
end

n,m,r = array_size
puts "array size = #{n}, repeat = #{r}\n\n"

system "#{RUBY}   bench.rb          float  mul #{n} #{r}" if ruby_narray
system "python bench.py numeric  float  mul #{n} #{r}" if python_numeric
system "python bench.py numarray float  mul #{n} #{r}" if python_numarray
system "python bench.py numpy    float  mul #{n} #{r}" if python_numpy
system "octave -qf bench.m       float  mul #{n} #{r}" if octave         
puts

system "#{RUBY}   bench.rb          int  add #{n} #{r}" if ruby_narray	 
system "python bench.py numeric  int  add #{n} #{r}" if python_numeric 
system "python bench.py numarray int  add #{n} #{r}" if python_numarray
system "python bench.py numpy    int  add #{n} #{r}" if python_numpy	 
system "octave -qf bench.m       int  add #{n} #{r}" if octave         
puts

system "#{RUBY}   bench.rb          complex  mul #{n} #{r}" if ruby_narray	 
system "python bench.py numeric  complex  mul #{n} #{r}" if python_numeric 
system "python bench.py numarray complex  mul #{n} #{r}" if python_numarray
system "python bench.py numpy    complex  mul #{n} #{r}" if python_numpy	 
system "octave -qf bench.m       complex  mul #{n} #{r}" if octave         
puts

system "#{RUBY}   bench.rb          float_cross mul    #{m*2} #{r}" if ruby_narray	 
system "python bench.py numeric  float_cross mul    #{m*2} #{r}" if python_numeric 
system "python bench.py numarray float_cross mul    #{m*2} #{r}" if python_numarray
system "python bench.py numpy    float_cross mul    #{m*2} #{r}" if python_numpy	 
system "octave -qf bench.m       float_cross matmul #{m*2} #{r}" if octave         
puts

system "#{RUBY}   bench.rb          float_matrix mul    #{m} 4" if ruby_narray	 
system "python bench.py numeric  float_matrix matmul #{m} 4" if python_numeric 
system "python bench.py numarray float_matrix matmul #{m} 4" if python_numarray
system "python bench.py numpy    float_matrix matmul #{m} 4" if python_numpy	 
system "octave -qf bench.m       float_matrix matmul #{m} 4" if octave         
puts

system "#{RUBY}   bench.rb          float_solve solve #{m} 2" if ruby_narray	 
system "python bench.py numeric  float_solve solve #{m} 2" if python_numeric 
system "python bench.py numarray float_solve solve #{m} 2" if python_numarray
system "python bench.py numpy    float_solve solve #{m} 2" if python_numpy	 
system "octave -qf bench.m       float_solve solve #{m} 2" if octave         
puts

exit

bench.m

#! /bin/octave -qf

arg_list = argv();
TYPE   = arg_list{1};
OP     = arg_list{2};
ARRSZ  = str2num(arg_list{3});
REPEAT = str2num(arg_list{4});
n = ARRSZ;

switch(TYPE)
  case "float"
    a = linspace(0,n-1,n);
    b = linspace(0,n-1,n);
  case "int"
    a = int32(linspace(0,n-1,n));
    b = int32(linspace(0,n-1,n));
  case "complex"
    a = complex(linspace(0,n-1,n));
    b = complex(linspace(0,n-1,n));
  case "float_cross"
    a = linspace(0,n-1,n)';
    b = linspace(0,n-1,n);
  case "float_matrix"
    a = linspace(0,n*n-1,n*n);
    a = rem(a, n+1) + 1;
    a = reshape(a,n,n);
    b = linspace(0,n*n-1,n*n);
    b = rem(b, n-1) + 1;
    b = reshape(b,n,n);
  case "float_solve"
    a = linspace(0,n*n-1,n*n);
    a = rem(a, n+1) + 1;
    a = reshape(a,n,n);
    b = reshape(linspace(1,n*n,n*n),n,n);
endswitch

[t1, u1, s1] = cputime ();
switch(OP)
  case "add"
    for i = 1:REPEAT
      c = a + b;
    endfor
  case "mul"
    for i = 1:REPEAT
      c = a .* b;
    endfor
  case "matmul"
    for i = 1:REPEAT
      c = a * b;
    endfor
    #size(c)
  case "solve"
    for i = 1:REPEAT
      c = a \ b;
    endfor
    #size(c)
endswitch
[t2, u2, s2] = cputime ();

printf ("Octave type=%s size=%d op=%s repeat=%d  Time: %.2f sec", 
	TYPE,ARRSZ,OP,REPEAT,u2 - u1);

bench.py

import time
import sys

MODULE = sys.argv[1]
TYPE   = sys.argv[2]
OP     = sys.argv[3]
ARRSZ  = int(sys.argv[4])
REPEAT = int(sys.argv[5])

if MODULE=="numeric":
    from Numeric import *
    from LinearAlgebra import *
elif MODULE=="numarray":
    from numarray import *
    from LinearAlgebra import *
elif MODULE=="numpy":
    from numpy import *
    from numpy.linalg import solve

def bench_time(func,repeat=REPEAT):
  start = time.clock()
  for i in range(repeat):
    c = func()
    stop = time.clock()
  print "Python %s type=%s size=%d op=%s repeat=%d  Time: %.2f sec" % \
      (MODULE,TYPE,ARRSZ,OP,REPEAT,stop-start)
  #print shape(c)

n = ARRSZ

if MODULE=="numpy":
    def bench_array(type=float):
        return arange(ARRSZ,dtype=type)
    
    if TYPE=="float":
        a = bench_array(float)
        b = bench_array(float)
    elif TYPE=="int":
        a = bench_array(int)
        b = bench_array(int)
    elif TYPE=="complex":
        a = bench_array(complex)
        b = bench_array(complex)
    elif TYPE=="float_cross":
        a = reshape(arange(ARRSZ,dtype=float),(ARRSZ,1))
        b = reshape(arange(ARRSZ,dtype=float),(1,ARRSZ))
    elif TYPE=="float_matrix":
        a = reshape(arange(ARRSZ**2,dtype=float),(ARRSZ,ARRSZ))
        b = reshape(arange(ARRSZ**2,dtype=float),(ARRSZ,ARRSZ))
    elif TYPE=="float_solve":
        a = reshape(arange(n**2,dtype=float)%(n+1)+1,(n,n))
        b = reshape(arange(n**2,dtype=float)+1,(n,n))
else:
    def bench_array(type=float):
        return arrayrange(ARRSZ).astype(type)
    if TYPE=="float":
        a = bench_array(Float64)
        b = bench_array(Float64)
    elif TYPE=="int":
        a = bench_array(Int32)
        b = bench_array(Int32)
    elif TYPE=="complex":
        a = bench_array(Complex64)
        b = bench_array(Complex64)
    elif TYPE=="float_cross":
        a = reshape(arrayrange(ARRSZ),(ARRSZ,1)).astype(Float64)
        b = reshape(arrayrange(ARRSZ),(1,ARRSZ)).astype(Float64)
    elif TYPE=="float_matrix":
        a = reshape(arrayrange(ARRSZ**2),(ARRSZ,ARRSZ)).astype(Float64)
        b = reshape(arrayrange(ARRSZ**2),(ARRSZ,ARRSZ)).astype(Float64)
    elif TYPE=="float_solve":
        a = reshape(arrayrange(n*n)%(n+1)+1,(n,n)).astype(Float64)
        b = reshape(arrayrange(n*n)+1,(n,n)).astype(Float64)
    dot = matrixmultiply
    solve = solve_linear_equations

def lambda_add(a=a,b=b):    c = a+b; return c;
def lambda_mul(a=a,b=b):    c = a*b; return c;
def lambda_matmul(a=a,b=b): c = dot(a,b); return c;
def lambda_solve(a=a,b=b):  c = solve(a,b); return c;

if OP=="add":
    bench_time(lambda_add)
elif OP=="mul":
    bench_time(lambda_mul)
elif OP=="matmul":
    bench_time(lambda_matmul)
elif OP=="solve":
    bench_time(lambda_solve)

bench.rb

require 'narray'
T = (RUBY_VERSION<"1.8.0") ? Time : Process

TYPE   = ARGV[0]
OP     = ARGV[1]
ARRSZ  = Integer(ARGV[2])
REPEAT = Integer(ARGV[3])

def bench_array(type=Float)
  [ NArray.new(type,ARRSZ).indgen!,
    NArray.new(type,ARRSZ).indgen!  ]
end

def bench_time(n=REPEAT)
  t1 = T.times.utime
   for i in 1..n
     yield
   end
  t = T.times.utime - t1
  printf "Ruby NArray type=%s size=%d op=%s repeat=%d  Time: %.2f sec\n",
    TYPE,ARRSZ,OP,REPEAT,t
end

n = ARRSZ

case TYPE
when "float"
  a,b = bench_array(Float)
when "int"
  a,b = bench_array(Integer)
when "complex"
  a,b = bench_array(Complex)
when "float_cross"
  a = NArray.float(n,1).indgen!
  b = NArray.float(1,n).indgen!
when "float_matrix"
  a = NArray.float(n,n).indgen!
  a = a % (n+1) + 1
  a = NMatrix.ref(a)#.transpose
  b = NArray.float(n,n).indgen!
  b = b % (n-1) + 1
  b = NMatrix.ref(b)#.transpose
when "float_solve"
  a = NMatrix.float(n,n).indgen!(1).transpose
  b = NArray.float(n,n).indgen!
  b = b % (n+1) + 1
  b = NMatrix.ref(b).transpose
end

c = 0

case OP
when "add"
  bench_time{ c = a+b }
when "mul"
  bench_time{ c = a*b }
when "solve"
  bench_time{ c = a/b }
end

dummy.m

dummy.py

import sys

MODULE = sys.argv[1]

if MODULE=="numeric":
    from Numeric import *
    from LinearAlgebra import *
elif MODULE=="numarray":
    from numarray import *
    from LinearAlgebra import *
elif MODULE=="numpy":
    from numpy import *
    from numpy.linalg import solve

dummy.rb

require 'narray'