Paste: spectral-norm light CSE+simd opt. (used by ATS/C in CLBG)
Author: | pruned |
Mode: | factor |
Date: | Sun, 30 May 2010 07:54:53 |
Plain Text |
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
Author: | pruned |
Mode: | factor |
Date: | Sun, 30 May 2010 08:14:52 |
Plain Text |
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 iota [| i |
n 2 / iota 0.0 dup double-2-boa
[| j |
u i j quot call v+
] reduce first2 +
] double-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 ] [ swap eval2-A' ] bi-curry bi* 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 ] [ swap eval2-A ] bi-curry bi* v* ; inline
: eval-At-times-u ( u n -- seq )
[ byte-array>double-2-array ] dip
[ (eval-At-times-u) ] inner-loop ; inline
: eval-AtA-times-u ( u n -- seq )
[ byte-array>double-2-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
New Annotation