extra/benchmark/spectral-norm/spectral-norm.factor

   1 ! Factor port of
   2 ! http://shootout.alioth.debian.org/gp4/benchmark.php?test=spectralnorm&lang=all
   3 USING: specialized-arrays kernel math math.functions
   4 math.vectors sequences sequences.private prettyprint words hints
   5 locals ;
   6 SPECIALIZED-ARRAY: double
   7 IN: benchmark.spectral-norm
   8
   9 :: inner-loop ( u n quot -- seq )
  10     n [| i |
  11         n 0.0 [| j |
  12             u i j quot call +
  13         ] reduce
  14     ] double-array{ } map-as ; inline
  15
  16 : eval-A ( i j -- n )
  17     [ >float ] bi@
  18     [ drop ] [ + [ ] [ 1 + ] bi * 0.5 * ] 2bi
  19     + 1 + recip ; inline
  20
  21 : (eval-A-times-u) ( u i j -- x )
  22     tuck [ swap nth-unsafe ] [ eval-A ] 2bi* * ; inline
  23
  24 : eval-A-times-u ( n u -- seq )
  25     [ (eval-A-times-u) ] inner-loop ; inline
  26
  27 : (eval-At-times-u) ( u i j -- x )
  28     tuck [ swap nth-unsafe ] [ swap eval-A ] 2bi* * ; inline
  29
  30 : eval-At-times-u ( u n -- seq )
  31     [ (eval-At-times-u) ] inner-loop ; inline
  32
  33 : eval-AtA-times-u ( u n -- seq )
  34     [ eval-A-times-u ] [ eval-At-times-u ] bi ; inline
  35
  36 : ones ( n -- seq ) [ 1.0 ] double-array{ } replicate-as ; inline
  37
  38 :: u/v ( n -- u v )
  39     n ones dup
  40     10 [
  41         drop
  42         n eval-AtA-times-u
  43         [ n eval-AtA-times-u ] keep
  44     ] times ; inline
  45
  46 : spectral-norm ( n -- norm )
  47     u/v [ v. ] [ norm-sq ] bi /f sqrt ;
  48
  49 HINTS: spectral-norm fixnum ;
  50
  51 : spectral-norm-main ( -- )
  52     2000 spectral-norm . ;
  53
  54 MAIN: spectral-norm-main