+ added spectralnorm benchmark

git-svn-id: trunk@4868 -

.gitattributes

@@ -5110,6 +5110,7 @@ tests/bench/shootout/io/recursive-output.txt svneol=native#text/plain
 tests/bench/shootout/io/rev.out -text
 tests/bench/shootout/io/revfile.in -text
 tests/bench/shootout/io/revfile.out -text
+tests/bench/shootout/io/spectralnorm-output.txt svneol=native#text/plain
 tests/bench/shootout/io/sumcol.in -text
 tests/bench/shootout/io/sumcol.out -text
 tests/bench/shootout/io/wc.in -text
@@ -5153,6 +5154,7 @@ tests/bench/shootout/src/binarytrees.pp svneol=native#text/plain
 tests/bench/shootout/src/hello.pp svneol=native#text/plain
 tests/bench/shootout/src/nsieve.pp svneol=native#text/plain
 tests/bench/shootout/src/recursive.pp svneol=native#text/plain
+tests/bench/shootout/src/spectralnorm.pp svneol=native#text/plain
 tests/bench/shootout/src/sumcol.pp svneol=native#text/plain
 tests/bench/shortbench.pp svneol=native#text/plain
 tests/bench/timer.pas svneol=native#text/plain

tests/bench/shootout/io/spectralnorm-output.txt

@@ -0,0 +1 @@
+1.274219991

tests/bench/shootout/src/spectralnorm.pp

@@ -0,0 +1,67 @@
+{ The Computer Language Shootout
+  http://shootout.alioth.debian.org
+
+  contributed by Ian Osgood
+  modified by Vincent Snijders
+}
+
+program spectralNorm;
+
+
+{$mode objfpc}{$inline on}
+
+var n,i : integer;
+    u,v,tmp : array of double;
+    vBv,vv : double;
+
+function A(i,j : integer): double; inline;
+begin
+  A := 1 / ((i+j)*(i+j+1) div 2 + i+1);
+end;
+
+procedure mulAv(var v, Av : array of double);
+var i,j : integer;
+begin
+  for i := low(Av) to high(Av) do
+  begin
+    Av[i] := 0.0;
+    for j := low(v) to high(v) do
+      Av[i] := Av[i] + A(i,j) * v[j];
+  end;
+end;
+
+procedure mulAtv(var v, Atv : array of double);
+var i,j : integer;
+begin
+  for i := low(Atv) to high(Atv) do
+  begin
+    Atv[i] := 0.0;
+    for j := low(v) to high(v) do
+      Atv[i] := Atv[i] + A(j,i) * v[j];
+  end;
+end;
+
+procedure mulAtAv(var v, AtAv : array of double);
+begin
+  mulAv(v, tmp);
+  mulAtv(tmp, AtAv);
+end;
+
+begin
+  Val(paramstr(1), n, i);
+  SetLength(u, n);
+  SetLength(v, n);
+  SetLength(tmp, n);
+
+  for i := low(u) to high(u) do u[i] := 1.0;
+
+  for i := 1 to 10 do begin mulAtAv(u,v); mulAtAv(v,u) end;
+
+  for i := low(u) to high(u) do
+  begin
+    vBv := vBv + u[i]*v[i];
+    vv  := vv  + v[i]*v[i];
+  end;
+
+  writeln(sqrt(vBv/vv):0:9);
+end.