core/math/intervals/intervals.factor

   1 ! Copyright (C) 2007 Slava Pestov.
   2 ! See http://factorcode.org/license.txt for BSD license.
   3 ! Based on Slate's src/unfinished/interval.slate by Brian Rice.
   4 USING: kernel sequences arrays math combinators math.order ;
   5 IN: math.intervals
   6
   7 TUPLE: interval from to ;
   8
   9 C: <interval> interval
  10
  11 : open-point ( n -- endpoint ) f 2array ;
  12
  13 : closed-point ( n -- endpoint ) t 2array ;
  14
  15 : [a,b] ( a b -- interval )
  16     >r closed-point r> closed-point <interval> ;
  17
  18 : (a,b) ( a b -- interval )
  19     >r open-point r> open-point <interval> ;
  20
  21 : [a,b) ( a b -- interval )
  22     >r closed-point r> open-point <interval> ;
  23
  24 : (a,b] ( a b -- interval )
  25     >r open-point r> closed-point <interval> ;
  26
  27 : [a,a] ( a -- interval ) closed-point dup <interval> ;
  28
  29 : [-inf,a] ( a -- interval ) -1./0. swap [a,b] ;
  30
  31 : [-inf,a) ( a -- interval ) -1./0. swap [a,b) ;
  32
  33 : [a,inf] ( a -- interval ) 1./0. [a,b] ;
  34
  35 : (a,inf] ( a -- interval ) 1./0. (a,b] ;
  36
  37 : compare-endpoints ( p1 p2 quot -- ? )
  38     >r over first over first r> call [
  39         2drop t
  40     ] [
  41         over first over first = [
  42             swap second swap second not or
  43         ] [
  44             2drop f
  45         ] if
  46     ] if ; inline
  47
  48 : endpoint< ( p1 p2 -- ? ) [ < ] compare-endpoints ;
  49
  50 : endpoint<= ( p1 p2 -- ? ) [ endpoint< ] 2keep = or ;
  51
  52 : endpoint> ( p1 p2 -- ? ) [ > ] compare-endpoints ;
  53
  54 : endpoint>= ( p1 p2 -- ? ) [ endpoint> ] 2keep = or ;
  55
  56 : endpoint-min ( p1 p2 -- p3 ) [ endpoint< ] most ;
  57
  58 : endpoint-max ( p1 p2 -- p3 ) [ endpoint> ] most ;
  59
  60 : interval>points ( int -- from to )
  61     dup interval-from swap interval-to ;
  62
  63 : points>interval ( seq -- interval )
  64     dup first
  65     [ [ endpoint-min ] reduce ] 2keep
  66     [ endpoint-max ] reduce <interval> ;
  67
  68 : (interval-op) ( p1 p2 quot -- p3 )
  69     2over >r >r
  70     >r [ first ] bi@ r> call
  71     r> r> [ second ] both? 2array ; inline
  72
  73 : interval-op ( i1 i2 quot -- i3 )
  74     pick interval-from pick interval-from pick (interval-op) >r
  75     pick interval-to pick interval-from pick (interval-op) >r
  76     pick interval-to pick interval-to pick (interval-op) >r
  77     pick interval-from pick interval-to pick (interval-op) >r
  78     3drop r> r> r> r> 4array points>interval ; inline
  79
  80 : interval+ ( i1 i2 -- i3 ) [ + ] interval-op ;
  81
  82 : interval- ( i1 i2 -- i3 ) [ - ] interval-op ;
  83
  84 : interval* ( i1 i2 -- i3 ) [ * ] interval-op ;
  85
  86 : interval-integer-op ( i1 i2 quot -- i3 )
  87     >r 2dup
  88     [ interval>points [ first integer? ] both? ] both?
  89     r> [ 2drop f ] if ; inline
  90
  91 : interval-1+ ( i1 -- i2 ) 1 [a,a] interval+ ;
  92
  93 : interval-1- ( i1 -- i2 ) -1 [a,a] interval+ ;
  94
  95 : interval-neg ( i1 -- i2 ) -1 [a,a] interval* ;
  96
  97 : interval-bitnot ( i1 -- i2 ) interval-neg interval-1- ;
  98
  99 : interval-sq ( i1 -- i2 ) dup interval* ;
 100
 101 : make-interval ( from to -- int )
 102     over first over first {
 103         { [ 2dup > ] [ 2drop 2drop f ] }
 104         { [ 2dup = ] [
 105             2drop over second over second and
 106             [ <interval> ] [ 2drop f ] if
 107         ] }
 108         [ 2drop <interval> ]
 109     } cond ;
 110
 111 : interval-intersect ( i1 i2 -- i3 )
 112     2dup and [
 113         [ interval>points ] bi@ swapd
 114         [ swap endpoint> ] most
 115         >r [ swap endpoint< ] most r>
 116         make-interval
 117     ] [
 118         or
 119     ] if ;
 120
 121 : interval-union ( i1 i2 -- i3 )
 122     2dup and [
 123         [ interval>points 2array ] bi@ append points>interval
 124     ] [
 125         2drop f
 126     ] if ;
 127
 128 : interval-subset? ( i1 i2 -- ? )
 129     dupd interval-intersect = ;
 130
 131 : interval-contains? ( x int -- ? )
 132     >r [a,a] r> interval-subset? ;
 133
 134 : interval-singleton? ( int -- ? )
 135     interval>points
 136     2dup [ second ] bi@ and
 137     [ [ first ] bi@ = ]
 138     [ 2drop f ] if ;
 139
 140 : interval-length ( int -- n )
 141     dup
 142     [ interval>points [ first ] bi@ swap - ]
 143     [ drop 0 ] if ;
 144
 145 : interval-closure ( i1 -- i2 )
 146     dup [ interval>points [ first ] bi@ [a,b] ] when ;
 147
 148 : interval-shift ( i1 i2 -- i3 )
 149     #! Inaccurate; could be tighter
 150     [ [ shift ] interval-op ] interval-integer-op interval-closure ;
 151
 152 : interval-shift-safe ( i1 i2 -- i3 )
 153     dup interval-to first 100 > [
 154         2drop f
 155     ] [
 156         interval-shift
 157     ] if ;
 158
 159 : interval-max ( i1 i2 -- i3 )
 160     #! Inaccurate; could be tighter
 161     [ max ] interval-op interval-closure ;
 162
 163 : interval-min ( i1 i2 -- i3 )
 164     #! Inaccurate; could be tighter
 165     [ min ] interval-op interval-closure ;
 166
 167 : interval-interior ( i1 -- i2 )
 168     interval>points [ first ] bi@ (a,b) ;
 169
 170 : interval-division-op ( i1 i2 quot -- i3 )
 171     >r 0 over interval-closure interval-contains?
 172     [ 2drop f ] r> if ; inline
 173
 174 : interval/ ( i1 i2 -- i3 )
 175     [ [ / ] interval-op ] interval-division-op ;
 176
 177 : interval/i ( i1 i2 -- i3 )
 178     [
 179         [ [ /i ] interval-op ] interval-integer-op
 180     ] interval-division-op interval-closure ;
 181
 182 : interval-recip ( i1 -- i2 ) 1 [a,a] swap interval/ ;
 183
 184 : interval-2/ ( i1 -- i2 ) -1 [a,a] interval-shift ;
 185
 186 SYMBOL: incomparable
 187
 188 : left-endpoint-< ( i1 i2 -- ? )
 189     [ swap interval-subset? ] 2keep
 190     [ nip interval-singleton? ] 2keep
 191     [ interval-from ] bi@ =
 192     and and ;
 193
 194 : right-endpoint-< ( i1 i2 -- ? )
 195     [ interval-subset? ] 2keep
 196     [ drop interval-singleton? ] 2keep
 197     [ interval-to ] bi@ =
 198     and and ;
 199
 200 : (interval<) ( i1 i2 -- i1 i2 ? )
 201     over interval-from over interval-from endpoint< ;
 202
 203 : interval< ( i1 i2 -- ? )
 204     {
 205         { [ 2dup interval-intersect not ] [ (interval<) ] }
 206         { [ 2dup left-endpoint-< ] [ f ] }
 207         { [ 2dup right-endpoint-< ] [ f ] }
 208         [ incomparable ]
 209     } cond 2nip ;
 210
 211 : left-endpoint-<= ( i1 i2 -- ? )
 212     >r interval-from r> interval-to = ;
 213
 214 : right-endpoint-<= ( i1 i2 -- ? )
 215     >r interval-to r> interval-from = ;
 216
 217 : interval<= ( i1 i2 -- ? )
 218     {
 219         { [ 2dup interval-intersect not ] [ (interval<) ] }
 220         { [ 2dup right-endpoint-<= ] [ t ] }
 221         [ incomparable ]
 222     } cond 2nip ;
 223
 224 : interval> ( i1 i2 -- ? )
 225     swap interval< ;
 226
 227 : interval>= ( i1 i2 -- ? )
 228     swap interval<= ;
 229
 230 : assume< ( i1 i2 -- i3 )
 231     interval-to first [-inf,a) interval-intersect ;
 232
 233 : assume<= ( i1 i2 -- i3 )
 234     interval-to first [-inf,a] interval-intersect ;
 235
 236 : assume> ( i1 i2 -- i3 )
 237     interval-from first (a,inf] interval-intersect ;
 238
 239 : assume>= ( i1 i2 -- i3 )
 240     interval-to first [a,inf] interval-intersect ;
 241
 242 : integral-closure ( i1 -- i2 )
 243     dup interval-from first2 [ 1+ ] unless
 244     swap interval-to first2 [ 1- ] unless
 245     [a,b] ;