! double-2 hack by prunedtree: ! before: 1.271104937 seconds ( 1.274224152228617 ) ! now: 0.637893125 seconds ( 1.274224152228617 ) ! Factor port of ! http://shootout.alioth.debian.org/gp4/benchmark.php?test=spectralnorm&lang=all USING: alien.c-types specialized-arrays kernel math math.functions math.vectors sequences sequences.private prettyprint words typed locals math.vectors.simd ; SPECIALIZED-ARRAYS: double double-2 ; IN: benchmark.spectral-norm :: inner-loop ( u n quot -- seq ) n 2 / iota [| i | n iota 0.0 dup double-2-boa [| j | u i j quot call v+ ] reduce ] double-2-array{ } map-as ; inline : eval-A ( i j -- n ) [ >float ] bi@ [ drop ] [ + [ ] [ 1 + ] bi * 0.5 * ] 2bi + 1 + ; inline :: eval2-A ( i j -- n ) i 2 * j eval-A i 2 * 1 + j eval-A double-2-boa 1.0 dup double-2-boa swap v/ ; inline :: eval2-A' ( i j -- n ) i 2 * j swap eval-A i 2 * 1 + j swap eval-A double-2-boa 1.0 dup double-2-boa swap v/ ; inline : (eval-A-times-u) ( u i j -- x ) [ swap nth-unsafe ] [ eval2-A ] bi-curry bi* n*v ; inline : eval-A-times-u ( n u -- seq ) [ (eval-A-times-u) ] inner-loop ; inline : (eval-At-times-u) ( u i j -- x ) [ swap nth-unsafe ] [ eval2-A' ] bi-curry bi* n*v ; inline : eval-At-times-u ( u n -- seq ) [ byte-array>double-array ] dip [ (eval-At-times-u) ] inner-loop ; inline : eval-AtA-times-u ( u n -- seq ) [ byte-array>double-array ] dip [ eval-A-times-u ] [ eval-At-times-u ] bi ; inline : ones ( n -- seq ) 2 / [ 1.0 dup double-2-boa ] double-2-array{ } replicate-as ; inline :: u/v ( n -- u v ) n ones dup 10 [ drop n eval-AtA-times-u [ n eval-AtA-times-u ] keep ] times ; inline TYPED: spectral-norm ( n: fixnum -- norm ) u/v [ byte-array>double-array ] bi@ [ v. ] [ norm-sq ] bi /f sqrt ; : spectral-norm-main ( -- ) 2000 spectral-norm . ; MAIN: spectral-norm-main