From: Tommy Nordgren Date: 2008-01-02T23:09:22+09:00 Subject: Have detected bug in file complex.rb --Apple-Mail-2--215052369 Content-Type: text/plain; charset=US-ASCII; format=flowed; delsp=yes Content-Transfer-Encoding: 7bit I have detected a bug in the implementation of Complex division. The bug gives internal underflow or overflow which causes the result to contain nans or otherwise erroneous result when the result would in fact be perfectly representable with floating point real and imaginary parts. Erroneous results will ALWAYS be returned when the absolute value of the denominator exceeds the square root of the maximum representable float. The following patch fixes this problem,and some others Index: lib/complex.rb =================================================================== --- lib/complex.rb (revision 14846) +++ lib/complex.rb (working copy) @@ -192,8 +192,23 @@ # def / (other) if other.kind_of?(Complex) - self*other.conjugate/other.abs2 + denominator = other.abs2 + numerator = self * other.conjugate + if numerator.zero? or (denominator.kind_of?(Float) and ! denominator.finite?) + magn = other.abs + self * Complex(other.real / magn, - other.image / magn) / magn + else + numerator / denominator + end elsif Complex.generic?(other) + # If dividing with an integer, we must convert to a type that can divide the left operand without truncation + if other.kind_of?(Integer) + if defined?(Rational) + other = Rational(other,1) + else + other = Float(other) + end + end Complex(@real/other, @image/other) else x, y = other.coerce(self) An incomplete unit test file follows #!/usr/bin/ruby #Delete this line to test with standard distribution complex.rb $:.unshift Dir::getwd require 'test/unit' require 'complex' class TestComplexBasics < Test::Unit::TestCase def approx_equal?(c1,c2) minerrlimit = 100.0 * Float::MIN abs1 = c1.abs abs2 = c2.abs if abs1 > abs2 maxabs = abs1 else maxabs = abs2 end rel = maxabs * Float::EPSILON * 25 if rel < minerrlimit rel = minerrlimit end d = (c1 - c2).abs return d < rel end def test_plus_integer assert_equal(Complex(8,7),Complex(3,4)+Complex(5,3)) assert_equal(Complex(5,1),Complex(4,1) + 1) assert_equal(Complex(5,1),1+ Complex(4,1)) end def test_minus_integer assert_equal(Complex(2,4),Complex(3,2)-Complex(1,-2)) assert_equal(Complex(4,5),Complex(5,5) -1) assert_equal(Complex(1,-2), 2 - Complex(1,2)) end def test_times_integer assert_equal(Complex(-1,0),Complex(0,1) * Complex(0,1)) assert_equal(Complex(3,3), Complex(1,1) * 3) assert_equal(Complex(4,-2), 2 * Complex(2,-1)) end def test_plus_float assert(approx_equal? (Complex(1.0,4.0),Complex(2.0,2.0)+Complex(-1.0,2.0))) end def test_minus_float assert(approx_equal?(Complex(3.0,0.0),Complex(2.0,2.0)- Complex(-1.0,2.0))) end def test_times_float assert(approx_equal?(Complex(-2,1),Complex(1,2)*Complex(0,1))) end def test_division_float [Complex(1.0,25.0), Complex(8,97), Complex(Float::MAX/2, Float::MAX/3), Complex(100.0 / Float::MAX, 200.0 / Float::MAX)].each { | c | assert(approx_equal?(Complex(1.0,0.0), c / c)) } assert(approx_equal?(Complex(-1.0,0.0),Complex(0.0,1.0)/ Complex(0.0,-1.0))) end def test_eps_log assert(approx_equal? (Complex(2.0,1.0),Math.log(Math.exp(Complex(2.0,1.0))))) end end --Apple-Mail-2--215052369 Content-Disposition: attachment; filename=tc_complex.rb Content-Type: text/x-ruby-script; x-mac-creator=21526368; x-unix-mode=0644; x-mac-type=54455854; name="tc_complex.rb" Content-Transfer-Encoding: 7bit #!/usr/bin/ruby #Delete this line to test with standard distribution complex.rb $:.unshift Dir::getwd require 'test/unit' require 'complex' class TestComplexBasics < Test::Unit::TestCase def approx_equal?(c1,c2) minerrlimit = 100.0 * Float::MIN abs1 = c1.abs abs2 = c2.abs if abs1 > abs2 maxabs = abs1 else maxabs = abs2 end rel = maxabs * Float::EPSILON * 25 if rel < minerrlimit rel = minerrlimit end d = (c1 - c2).abs return d < rel end def test_plus_integer assert_equal(Complex(8,7),Complex(3,4)+Complex(5,3)) assert_equal(Complex(5,1),Complex(4,1) + 1) assert_equal(Complex(5,1),1+ Complex(4,1)) end def test_minus_integer assert_equal(Complex(2,4),Complex(3,2)-Complex(1,-2)) assert_equal(Complex(4,5),Complex(5,5) -1) assert_equal(Complex(1,-2), 2 - Complex(1,2)) end def test_times_integer assert_equal(Complex(-1,0),Complex(0,1) * Complex(0,1)) assert_equal(Complex(3,3), Complex(1,1) * 3) assert_equal(Complex(4,-2), 2 * Complex(2,-1)) end def test_plus_float assert(approx_equal?(Complex(1.0,4.0),Complex(2.0,2.0)+Complex(-1.0,2.0))) end def test_minus_float assert(approx_equal?(Complex(3.0,0.0),Complex(2.0,2.0)-Complex(-1.0,2.0))) end def test_times_float assert(approx_equal?(Complex(-2,1),Complex(1,2)*Complex(0,1))) end def test_division_float [Complex(1.0,25.0), Complex(8,97), Complex(Float::MAX/2, Float::MAX/3), Complex(100.0 / Float::MAX, 200.0 / Float::MAX)].each { | c | assert(approx_equal?(Complex(1.0,0.0), c / c)) } assert(approx_equal?(Complex(-1.0,0.0),Complex(0.0,1.0)/Complex(0.0,-1.0))) end def test_eps_log assert(approx_equal?(Complex(2.0,1.0),Math.log(Math.exp(Complex(2.0,1.0))))) end end --Apple-Mail-2--215052369 Content-Disposition: attachment; filename=ruby-changes.patch Content-Type: application/octet-stream; x-unix-mode=0644; name="ruby-changes.patch" Content-Transfer-Encoding: 7bit Index: lib/complex.rb =================================================================== --- lib/complex.rb (revision 14846) +++ lib/complex.rb (working copy) @@ -192,8 +192,23 @@ # def / (other) if other.kind_of?(Complex) - self*other.conjugate/other.abs2 + denominator = other.abs2 + numerator = self * other.conjugate + if numerator.zero? or (denominator.kind_of?(Float) and !denominator.finite?) + magn = other.abs + self * Complex(other.real / magn, - other.image / magn) / magn + else + numerator / denominator + end elsif Complex.generic?(other) + # If dividing with an integer, we must convert to a type that can divide the left operand without truncation + if other.kind_of?(Integer) + if defined?(Rational) + other = Rational(other,1) + else + other = Float(other) + end + end Complex(@real/other, @image/other) else x, y = other.coerce(self) --Apple-Mail-2--215052369 Content-Type: text/plain; charset=US-ASCII; format=flowed; delsp=yes Content-Transfer-Encoding: 7bit ----------------------------------- See the amazing new SF reel: Invasion of the man eating cucumbers from outer space. On congratulations for a fantastic parody, the producer replies : "What parody?" Tommy Nordgren tommy.nordgren@comhem.se --Apple-Mail-2--215052369--