extra/math/erato/erato.factor

   1 ! Copyright (c) 2007 Samuel Tardieu.
   2 ! See http://factorcode.org/license.txt for BSD license.
   3 USING: accessors bit-arrays fry kernel lists.lazy math math.functions
   4     math.primes.list math.ranges sequences ;
   5 IN: math.erato
   6
   7 <PRIVATE
   8
   9 TUPLE: erato limit bits latest ;
  10
  11 : ind ( n -- i )
  12     2/ 1- ; inline
  13
  14 : is-prime ( n limit -- bool )
  15     [ ind ] [ bits>> ] bi* nth ; inline
  16
  17 : indices ( n erato -- range )
  18     limit>> ind over 3 * ind spin <range> ;
  19
  20 : mark-multiples ( n erato -- )
  21     2dup [ sq ] [ limit>> ] bi* <= [
  22         [ indices ] keep bits>> '[ _ f -rot set-nth ] each
  23     ] [ 2drop ] if ;
  24
  25 : <erato> ( n -- erato )
  26     dup ind 1+ <bit-array> dup set-bits 1 erato boa ;
  27
  28 : next-prime ( erato -- prime/f )
  29     [ 2 + ] change-latest [ latest>> ] keep
  30     2dup limit>> <= [
  31         2dup is-prime [ dupd mark-multiples ] [ nip next-prime ] if
  32     ] [
  33         2drop f
  34     ] if ;
  35
  36 PRIVATE>
  37
  38 : lerato ( n -- lazy-list )
  39     dup 1000003 < [
  40         0 primes-under-million seq>list swap '[ _ <= ] lwhile
  41     ] [
  42         <erato> 2 [ drop next-prime ] with lfrom-by [ ] lwhile
  43     ] if ;