From: chr_rippel@...
Date: 2002-05-06T21:37:14+09:00
Subject: Fixes for the new step methods in the range.c & numeric.c

This is a multi-part message in MIME format.

------=_NextPart_000_0001_01C1F50B.6E551680
Content-Type: text/plain;
	charset="us-ascii"
Content-Transfer-Encoding: 7bit

Hi,

an included patch for the step method in range.c is a simple
bug-fix - for example

---
class Class
  alias succ superclass
end

# currently throws an exception complaining that succ for nil
# is undefined
(Fixnum..Object).step(2) {|x| p x }
---

The patch for the step method in numeric.c adds a sequence 
end point adjustment strategy - #step_fuzzyness measures the
``fuzzyness'' of step-length additions around the end-point
and adds a safety margin (there exist probably much better 
step_fuzzyness measures) - in pure Ruby:  

---
class Float
    def test_step(b,s)
      n = ((b - self)/ s).round
      last = n*s + self
      
      fuzzyness = step_fuzzyness(n,s)
      if (s > 0 && last > b + fuzzyness) || \
         (s < 0 && last < b - fuzzyness)
         n-=1
         last = n*s + self
      end
      
      if (last - b).abs < fuzzyness
        last = b
      end
        
      (n-1).times {|i| yield i*s + self }
      yield last
	self
    end
    
    EPSILON = 2.0**(-43)
    def step_fuzzyness(n, s)
      tmp = ((n-1)..(n+1)).collect do |i| 
       ((i*s + self - ((i-1)*s + self)) - s).abs
      end
      tmp.push s.abs*EPSILON
      2*tmp.max
    end
end


# with this patch on my machine ...
M = 52343
T = 12418.979
K = 20300012
L = 15673

p b = (K/T)                  # => 1634.595887
p e = (-L*M + K)/T           # => -64423.31749
p s = -M/T                   # => -4.214758717
p err = M/T * (2.0**(-40))# => 3.677147842e-12

# with this patch on my machine

last = nil
b.step(e,s) {|last| }
p e - last       # => 0.0

b.step(e+err,s) {|last| }
p (e+err) - last # => -4.214758717

b.step(e-err,s) {|last| }
p (e-err) - last # => -7.275957614e-12


# a bigger error is necessary to demonstrate
# fuzzyness in pure Ruby 

p err = M/T * (2.0**(-38.00))

b.test_step(e,s) {|last| }
p e - last       # => 0.0

b.test_step(e+err,s) {|last| }
p (e+err) - last # => -4.214758717

b.test_step(e-err,s) {|last| }
p (e-err) - last # => -1.455191523e-11
---


/Christoph

------=_NextPart_000_0001_01C1F50B.6E551680
Content-Type: application/octet-stream;
	name="range.patch"
Content-Transfer-Encoding: quoted-printable
Content-Disposition: attachment;
	filename="range.patch"

--- range.c.1.34	Fri May  3 18:44:42 2002=0A=
+++ range.c	Fri May  3 00:23:28 2002=0A=
@@ -71,7 +71,7 @@=0A=
     VALUE obj;=0A=
 {=0A=
     VALUE beg, end, flags;=0A=
-    =0A=
+=0A=
     rb_scan_args(argc, argv, "21", &beg, &end, &flags);=0A=
     /* Ranges are immutable, so that they should be initialized only =
once. */=0A=
     if (rb_ivar_defined(obj, id_beg)) {=0A=
@@ -304,26 +304,21 @@=0A=
     else {			/* generic each */=0A=
 	VALUE v =3D b;=0A=
 	long lim =3D NUM2LONG(step);=0A=
-	long i;=0A=
+	long i =3D lim;=0A=
 =0A=
 	if (lim <=3D 0) {=0A=
 	    rb_raise(rb_eArgError, "step can't be <=3D 0");=0A=
 	}=0A=
-	if (EXCL(range)) {=0A=
-	    while (r_lt(v, e)) {=0A=
-		if (r_eq(v, e)) break;=0A=
-		rb_yield(v);=0A=
-		for (i=3D0; i<lim; i++) =0A=
-		    v =3D rb_funcall(v, id_succ, 0, 0);=0A=
-	    }=0A=
-	}=0A=
-	else {=0A=
-	    while (r_le(v, e)) {=0A=
+	for(;NUM2INT(rb_funcall(v, id_cmp, 1, e)) < 0; v =3D rb_funcall(v, =
id_succ, 0, 0))=0A=
+		if (i =3D=3D lim) {=0A=
+			rb_yield(v);=0A=
+			i =3D 1;=0A=
+		} else {=0A=
+			i++;=0A=
+		}=0A=
+	if (!(EXCL(range)) && (i =3D=3D lim) && \=0A=
+	    (v =3D=3D e || NUM2INT(rb_funcall(v, id_cmp, 1, e)) =3D=3D 0)) {=0A=
 		rb_yield(v);=0A=
-		if (r_eq(v, e)) break;=0A=
-		for (i=3D0; i<lim; i++) =0A=
-		    v =3D rb_funcall(v, id_succ, 0, 0);=0A=
-	    }=0A=
 	}=0A=
     }=0A=
     return range;=0A=

------=_NextPart_000_0001_01C1F50B.6E551680
Content-Type: application/octet-stream;
	name="numeric.patch"
Content-Transfer-Encoding: quoted-printable
Content-Disposition: attachment;
	filename="numeric.patch"

--- numeric.c.1.47	Fri May  3 19:02:33 2002=0A=
+++ numeric.c	Sun May  5 12:17:11 2002=0A=
@@ -224,7 +224,7 @@=0A=
 	return rb_str_new2(value < 0 ? "-Infinity" : "Infinity");=0A=
     else if(isnan(value))=0A=
 	return rb_str_new2("NaN");=0A=
-    =0A=
+=0A=
     avalue =3D fabs(value);=0A=
     if (avalue =3D=3D 0.0) {=0A=
 	fmt =3D "%.1f";=0A=
@@ -234,13 +234,13 @@=0A=
 	while (d1 < 1.0) d1 *=3D 10.0;=0A=
 	d1 =3D modf(d1, &d2);=0A=
 	if (d1 =3D=3D 0) fmt =3D "%.1e";=0A=
-    }    =0A=
+    }=0A=
     else if (avalue >=3D 1.0e10) {=0A=
 	d1 =3D avalue;=0A=
 	while (d1 > 10.0) d1 /=3D 10.0;=0A=
 	d1 =3D modf(d1, &d2);=0A=
 	if (d1 =3D=3D 0) fmt =3D "%.1e";=0A=
-    }    =0A=
+    }=0A=
     else if ((d1 =3D modf(value, &d2)) =3D=3D 0) {=0A=
 	fmt =3D "%.1f";=0A=
     }=0A=
@@ -641,7 +641,7 @@=0A=
 =0A=
 static VALUE flo_is_nan_p(num)=0A=
      VALUE num;=0A=
-{     =0A=
+{=0A=
 =0A=
   double value =3D RFLOAT(num)->value;=0A=
 =0A=
@@ -650,7 +650,7 @@=0A=
 =0A=
 static VALUE flo_is_infinite_p(num)=0A=
      VALUE num;=0A=
-{     =0A=
+{=0A=
   double value =3D RFLOAT(num)->value;=0A=
 =0A=
   if (isinf(value)) {=0A=
@@ -662,12 +662,12 @@=0A=
 =0A=
 static VALUE flo_is_finite_p(num)=0A=
      VALUE num;=0A=
-{     =0A=
+{=0A=
   double value =3D RFLOAT(num)->value;=0A=
 =0A=
   if (isinf(value) || isnan(value))=0A=
     return Qfalse;=0A=
-  =0A=
+=0A=
   return Qtrue;=0A=
 }=0A=
 =0A=
@@ -762,6 +762,25 @@=0A=
     return flo_truncate(rb_Float(num));=0A=
 }=0A=
 =0A=
+=0A=
+static double=0A=
+step_fuzzyness(n,a,b,c)=0A=
+	long  n;=0A=
+	double a,b,c;=0A=
+{=0A=
+	/* eps =3D 2^(-43)  */=0A=
+	const double eps =3D 1.136868377e-13;=0A=
+	double n0 =3D (n-2)*c + a;=0A=
+	double n1 =3D (n-1)*c + a;=0A=
+	double n2 =3D (n)*c + a;=0A=
+	double n3 =3D ((n)+1)*c + a;=0A=
+	/* max-macro is not universal */=0A=
+	#define fuzzy_max(a,b)    (((a) > (b)) ? (a) : (b))=0A=
+	return 2*fuzzy_max(fuzzy_max(fabs(n1-n0-c),fabs(n2-n1-c)),\=0A=
+	                   fuzzy_max(fabs(n3-n2-c),fabs(c)*eps));=0A=
+}=0A=
+=0A=
+=0A=
 static VALUE=0A=
 num_step(argc, argv, from)=0A=
     int argc;=0A=
@@ -798,17 +817,32 @@=0A=
 	}=0A=
     }=0A=
     else if (TYPE(from) =3D=3D T_FLOAT || TYPE(to) =3D=3D T_FLOAT || =
TYPE(step) =3D=3D T_FLOAT) {=0A=
-	const double epsilon =3D 2.2204460492503131E-16;=0A=
 	double beg =3D NUM2DBL(from);=0A=
 	double end =3D NUM2DBL(to);=0A=
 	double unit =3D NUM2DBL(step);=0A=
-	double n =3D (end - beg)/unit;=0A=
-	long i;=0A=
+	double fuzzy,last;=0A=
+	long   n,i;=0A=
+	n =3D floor((end - beg)/unit + 0.5);=0A=
+=0A=
+	if (n < 0) {=0A=
+		return from;=0A=
+	}=0A=
+	fuzzy =3D step_fuzzyness(n,beg,end,unit);=0A=
+	last =3D n* unit + beg;=0A=
 =0A=
-	n =3D floor(n + n*epsilon) + 1;=0A=
-	for (i=3D0; i<n; i++) {=0A=
-	    rb_yield(rb_float_new(i*unit+beg));=0A=
+	if ((unit > 0.0 && last > end + fuzzy) || (unit < 0.0 && last < end - =
fuzzy)) {=0A=
+		n--;=0A=
+		last =3D n*unit + beg;=0A=
 	}=0A=
+	if (fabs(last - end) < fuzzy) {=0A=
+		last =3D end;=0A=
+	}=0A=
+=0A=
+	for (i =3D 0; i<n; i++) {=0A=
+		rb_yield(rb_float_new(i*unit + beg));=0A=
+	}=0A=
+	rb_yield(rb_float_new(last));=0A=
+=0A=
     }=0A=
     else {=0A=
 	VALUE i =3D from;=0A=

------=_NextPart_000_0001_01C1F50B.6E551680--