Set of real numbers

From Rosetta Code
Revision as of 17:47, 4 October 2011 by Rdm (talk | contribs) (J: Optional Work)
Set of real numbers is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

All real numbers form the uncountable set ℝ. Among its subsets, relatively simple are the convex sets, each expressed as a range between two real numbers a and b where ab. There are actually four cases for the meaning of "between", depending on open or closed boundary:

  • [a, b]: {x | ax and xb }
  • (a, b): {x | a < x and x < b }
  • [a, b): {x | ax and x < b }
  • (a, b]: {x | a < x and xb }

Note that if a = b, of the four only [a, a] would be non-empty.

Task

  • Device a way to represent any set of real numbers, for the definition of 'any' in the implementation notes below.
  • Provide methods for these common set operations (x is a real number; A and B are sets):
  • xA: determine if x is an element of A
example: 1 is in [1, 2), while 2, 3, ... are not.
  • AB: union of A and B, i.e. {x | xA or xB}
example: [0, 2) ∪ (1, 3) = [0, 3); [0, 1) ∪ (2, 3] = well, [0, 1) ∪ (2, 3]
  • AB: intersection of A and B, i.e. {x | xA and xB}
example: [0, 2) ∩ (1, 3) = (1, 2); [0, 1) ∩ (2, 3] = empty set
  • A - B: difference between A and B, also written as A \ B, i.e. {x | xA and xB}
example: [0, 2) − (1, 3) = [0, 1]
  • Test your implementation by checking if numbers 0, 1, and 2 are in any of the following sets:
  • (0, 1] ∪ [0, 2)
  • [0, 2) ∩ (1, 2]
  • [0, 3) − (0, 1)
  • [0, 3) − [0, 1]

Implementation notes

  • 'Any' real set means 'sets that can be expressed as the union of a finite number of convex real sets'. Cantor's set needs not apply.
  • Infinities should be handled gracefully; indeterminate numbers (NaN) can be ignored.
  • You can use your machine's native real number representation, which is probably IEEE floating point, and assume it's good enough (it usually is).

Optional work

Define A = {x | 0 < x < 10 and |sin(π x²)| > 1/2 }, B = {x | 0 < x < 10 and |sin(π x)| > 1/2}, calculate the length of the real axis covered by the set AB. Note that |sin(π x)| > 1/2 is the same as n + 1/6 < x < n + 5/6 for all integers n; your program does not need to derive this by itself.

J

In essence, this looks like building a restricted set of statements. So we build a specialized parser and expression builder:

<lang j>has=: 1 :'(interval^:(0=L.) m)`:6' ing=: `

interval=: 3 :0

 assert. 5=#words=. ;:y
 assert. (0 { words) e. ;:'[('
 assert. (2 { words) e. ;:','
 assert. (4 { words) e. ;:'])'
 'lo hi'=.(1 3{0".L:0 words)
 'cL cH'=.0 4{words e.;:'[]'
 (lo&(<`<:@.cL) *. hi&(>`>:@.cH))ing

)


in=: 4 :'y has x' union=: 4 :'(x has +. y has)ing' intersect=: 4 :'(x has *. y has)ing' without=: 4 :'(x has *. [: -. y has)ing'</lang>

With this in place, the required examples look like this:

<lang j> ('(0,1]' union '[0,2)')has 0 1 2 1 1 0

  ('[0,2)' intersect '(1,2]')has 0 1 2

0 0 0

  ('[0,3)' without '(0,1]')has 0 1 2

1 0 1

  ('[0,3)' without '(0,1)')has 0 1 2

1 1 1

  ('[0,3)' without '[0,1]')has 0 1 2

0 0 1</lang>

Note that without the arguments these wind up being expressions. For example:

<lang j> ('(0,1]' union '[0,2)')has (0&< *. 1&>:) +. 0&<: *. 2&></lang>

In other words, this is a statement built up from inequality terminals (where each inequality is bound to a constant) and the terminals are combined with logical operations.

Optional Work

The optional work centers around expressions where the absolute value of sin pi * n is 0.5. It would be nice if J had an arcsine which gave all values within a range, but it does not have that. So:

<lang j> 1p_1 * _1 o. 0.5 0.166667</lang>

(Note on notation: 1 o. is sine in J, and 2 o. is cosine -- the mnemonic is that sine is an odd function and cosine is an even function, the practical value is that sine, cosine and sine/cosine pairs can all be generated from the same "real" valued function. Similarly, _1 o. is arcsine and _2 o. is arcsine. Also 1p_1 is the reciprocal of pi. So the above tells us that the principal value for arc sine 0.5 is one sixth.)

<lang j> (#~ 0.5 = 1 |@o. 1r6p1&*) i. 30 1 5 7 11 13 17 19 23 25 29

  2 -~/\ (#~ 0.5 = 1 |@o. 1r6p1&*) i. 30

4 2 4 2 4 2 4 2 4</lang>

Here we see the integers which when multiplied by pi/6 give 0.5 for the absolute value of the sine, and their first difference. Thus:

<lang j>zeros0toN=: ((>: # ])[:+/\1,$&4 2@<.)&.(6&*)</lang>

is a function to generate the values which correspond to the boundaries of the intervals we want:

<lang j>zB=: zeros0toN 10 zA=: zeros0toN&.*: 10

  zA

0.408248 0.912871 1.08012 1.35401 1.47196 1.68325 1.77951 1.95789 2.04124 2.1984...

  zB

0.166667 0.833333 1.16667 1.83333 2.16667 2.83333 3.16667 3.83333 4.16667 4.8333...

  #zA

200

  #zB

20</lang>

And, here are the edges of the sets of intervals we need to consider.

To find the length of the the set A-B we can find the length of set A and subtract the length of the set A-B:

<lang j> (+/_2 -~/\zA) - +/,0>.zA (<.&{: - >.&{.)"1/&(_2 ]\ ]) zB 2.07586</lang>

Here, we have paired adjacent elements from the zero bounding list (non-overlapping infixes of length 2). For set A's length we sum the results of subtracting the smaller number of the pair from the larger. For set A-B's length we consider each combination of pairs from A and B and subtract the larger of the beginning values from the smaller of the ending values (and ignore any negative results).

Perl

<lang perl>use strict; use utf8;

  1. numbers used as boundaries to real sets. Each has 3 components:
  2. the real value x;
  3. a +/-1 indicating if it's x + ϵ or x - ϵ
  4. a 0/1 indicating if it's the left border or right border
  5. e.g. "[1.5, ..." is written "1.5, -1, 0", while "..., 2)" is "2, -1, 1"

package BNum;

use overload ( '""' => \&_str, '<=>' => \&_cmp, );

sub new { my $self = shift; bless [@_], ref $self || $self }

sub flip { my @a = @{+shift}; $a[2] = !$a[2]; bless \@a }

my $brackets = qw/ [ ( ) ] /; sub _str { my $v = sprintf "%.2f", $_[0][0]; $_[0][2] ? $v . ($_[0][1] == 1 ? "]" : ")") : ($_[0][1] == 1 ? "(" : "[" ) . $v; }

sub _cmp { my ($a, $b, $swap) = @_;

# if one of the argument is a normal number if ($swap) { return -_ncmp($a, $b) } if (!ref $b || !$b->isa(__PACKAGE__)) { return _ncmp($a, $b) }

$a->[0] <=> $b->[0] || $a->[1] <=> $b->[1] }

sub _ncmp { # $a is a BNum, $b is something comparable to a real my ($a, $b) = @_; $a->[0] <=> $b || $a->[1] <=> 0 }

package RealSet; use Carp; use overload ( '""' => \&_str, '|' => \&_or, '&' => \&_and, '~' => \&_neg, '-' => \&_diff, );

my %pm = qw/ [ -1 ( 1 ) -1 ] 1 /; sub range { my ($cls, $a, $b, $spec) = @_; $spec =~ /^( \[ | \( )( \) | \] )$/x or croak "bad spec $spec";

$a = BNum->new($a, $pm{$1}, 0); $b = BNum->new($b, $pm{$2}, 1); normalize($a < $b ? [$a, $b] : []) }

sub normalize { my @a = @{+shift}; # remove invalid or duplicate borders, such as "[2, 1]" or "3) [3" # note that "(a" == "a]" and "a)" == "[a", but "a)" < "(a" and # "[a" < "a]" for (my $i = $#a; $i > 0; $i --) { splice @a, $i - 1, 2 if $a[$i] <= $a[$i - 1] } bless \@a }

sub _str { my (@a, @s) = @{+shift} or return '()'; join " ∪ ", map { shift(@a).", ".shift(@a) } 0 .. $#a/2 }

sub _or { # we may have nested ranges now; let only outmost ones survive my $d = 0; normalize [ map { $_->[2] ? --$d ? () : ($_) : $d++ ? () : ($_) } sort{ $a <=> $b } @{+shift}, @{+shift} ]; }

sub _neg { normalize [ BNum->new('-inf', 1, 0), map($_->flip, @{+shift}), BNum->new('inf', -1, 1), ] }

sub _and { my $d = 0; normalize [ map { $_->[2] ? --$d ? ($_) : () : $d++ ? ($_) : () } sort{ $a <=> $b } @{+shift}, @{+shift} ]; }

sub _diff { shift() & ~shift() }

sub has { my ($a, $b) = @_; for (my $i = 0; $i < $#$a; $i += 2) { return 1 if $a->[$i] <= $b && $b <= $a->[$i + 1] } return 0 }

sub len { my ($a, $l) = shift; for (my $i = 0; $i < $#$a; $i += 2) { $l += $a->[$i+1][0] - $a->[$i][0] } return $l }

package main; use List::Util 'reduce';

sub rng { RealSet->range(@_) } my @sets = ( rng(0, 1, '(]') | rng(0, 2, '[)'), rng(0, 2, '[)') & rng(0, 2, '(]'), rng(0, 3, '[)') - rng(0, 1, '()'), rng(0, 3, '[)') - rng(0, 1, '[]'), );

for my $i (0 .. $#sets) { print "Set $i = ", $sets[$i], ": "; for (0 .. 2) { print "has $_; " if $sets[$i]->has($_); } print "\n"; }

  1. optional task

print "\n####\n"; sub brev { # show only head and tail if string too long my $x = shift; return $x if length $x < 60; substr($x, 0, 30)." ... ".substr($x, -30, 30) }

  1. "|sin(x)| > 1/2" means (n + 1/6) pi < x < (n + 5/6) pi

my $x = reduce { $a | $b } map(rng(sqrt($_ + 1./6), sqrt($_ + 5./6), '()'), 0 .. 101); $x &= rng(0, 10, '()');

print "A\t", '= {x | 0 < x < 10 and |sin(π x²)| > 1/2 }', "\n\t= ", brev($x), "\n";

my $y = reduce { $a | $b } map { rng($_ + 1./6, $_ + 5./6, '()') } 0 .. 11; $y &= rng(0, 10, '()');

print "B\t", '= {x | 0 < x < 10 and |sin(π x)| > 1/2 }', "\n\t= ", brev($y), "\n";

my $z = $x - $y; print "A - B\t= ", brev($z), "\n\tlength = ", $z->len, "\n";</lang>output<lang>Set 0 = [0.00, 2.00): has 0; has 1; Set 1 = (0.00, 2.00): has 1; Set 2 = [0.00, 0.00] ∪ [1.00, 3.00): has 0; has 1; has 2; Set 3 = (1.00, 3.00): has 2;

A = {x | 0 < x < 10 and |sin(π x²)| > 1/2 }

       = (0.41, 0.91) ∪ (1.08, 1.35) ∪  ...  ∪ (9.91, 9.94) ∪ (9.96, 9.99)

B = {x | 0 < x < 10 and |sin(π x)| > 1/2 }

       = (0.17, 0.83) ∪ (1.17, 1.83) ∪  ...  ∪ (8.17, 8.83) ∪ (9.17, 9.83)

A - B = [0.83, 0.91) ∪ (1.08, 1.17] ∪ ... ∪ (9.91, 9.94) ∪ (9.96, 9.99)

       length = 2.07586484118467</lang>

Python

<lang python>class Setr():

   def __init__(self, lo, hi, includelo=True, includehi=False):
       self.eqn = "(%i<%sX<%s%i)" % (lo,
                                     '=' if includelo else ,
                                     '=' if includehi else ,
                                     hi)
   def includes(self, X):
       return eval(self.eqn, locals())
   def union(self, b):
       ans = Setr(0,0)
       ans.eqn = "(%sor%s)" % (self.eqn, b.eqn)
       return ans
   def intersection(self, b):
       ans = Setr(0,0)
       ans.eqn = "(%sand%s)" % (self.eqn, b.eqn)
       return ans
   def difference(self, b):
       ans = Setr(0,0)
       ans.eqn = "(%sand not%s)" % (self.eqn, b.eqn)
       return ans
   def union(self, b):
       ans = Setr(0,0)
       ans.eqn = "(%sor%s)" % (self.eqn, b.eqn)
       return ans
   def __repr__(self):
       return "Setr%s" % self.eqn


sets = [Setr(0,1, 0,1).union(Setr(0,2, 1,0)),

       Setr(0,2, 1,0).intersection(Setr(0,2, 0,1)),
       Setr(0,3, 1,0).difference(Setr(0,1, 0,0)),
       Setr(0,3, 1,0).difference(Setr(0,1, 1,1)),
       ]

settexts = '(0, 1] ∪ [0, 2);[0, 2) ∩ (1, 2];[0, 3) − (0, 1);[0, 3) − [0, 1]'.split(';')

for s,t in zip(sets, settexts):

   print("Set %s %s. %s" % (t,
                            ', '.join("%scludes %i"
                                    % ('in' if s.includes(v) else 'ex', v)
                                    for v in range(3)),
                            s.eqn))</lang>
Output
Set (0, 1] ∪ [0, 2) includes 0, includes 1, excludes 2. ((0<X<=1)or(0<=X<2))
Set [0, 2) ∩ (1, 2] excludes 0, includes 1, excludes 2. ((0<=X<2)and(0<X<=2))
Set [0, 3) − (0, 1) includes 0, includes 1, includes 2. ((0<=X<3)and not(0<X<1))
Set [0, 3) − [0, 1] excludes 0, excludes 1, includes 2. ((0<=X<3)and not(0<=X<=1))