Re: [PATCH] SEQ BUG

[Paul Eggert]

Pádraig Brady <address@hidden> writes:

> The attached patch handles this by
> only counting signficant digits from the operands.

I'd rather use the idea I proposed earlier.  Here's an implementation
of it, which works on all the test cases in your patch.  In addition,
it works on the wilder counterexamples I suggested (which alas we
can't put into the test suite since they're machine-specific).

2007-06-22  Paul Eggert  <address@hidden>

        * NEWS: seq no longer mishandles obvious cases like
        "seq 0 0.000001 0.000003" by omitting the last output number.
        * doc/coreutils.texi (seq invocation): Remove advice about workaround
        for seq off-by-one problem, since the bug is fixed now.  Replace
        it with more-generic advice about rounding error.
        * src/seq.c (long_double_format, print_numbers):
        New arg NUMERIC_FORMAT.  All uses changed.

2007-06-22  Pádraig Brady  <address@hidden>  (trivial change)

        * tests/seq/basic: Add test cases for seq off-by-one problem.

diff --git a/NEWS b/NEWS
index c587fe7..a588896 100644
--- a/NEWS
+++ b/NEWS
@@ -48,6 +48,9 @@ GNU coreutils NEWS                                    -*- 
outline -*-
   ln=target attribute) would mistakenly output the string "target"
   before the name of each symlink.  [introduced in coreutils-6.0]

+  seq no longer mishandles obvious cases like "seq 0 0.000001 0.000003",
+  so workarounds like "seq 0 0.000001 0.0000031" are no longer needed.
+
   split --line-bytes=N (-C N) no longer creates an empty file
   [this bug is present at least as far back as textutils-1.22 (Jan, 1997)]

diff --git a/doc/coreutils.texi b/doc/coreutils.texi
index 7290ab2..42558a3 100644
--- a/doc/coreutils.texi
+++ b/doc/coreutils.texi
@@ -14052,34 +14052,16 @@ $ seq 18446744073709551616 1 18446744073709551618
 18446744073709551618
 @end example

-Be careful when using @command{seq} with a fractional @var{increment};
-otherwise you may see surprising results.  Most people would expect to
-see @code{0.000003} printed as the last number in this example:
+Be careful when using @command{seq} with outlandish values: otherwise
+you may see surprising results, as @command{seq} uses floating point
+internally.  For example, on the x86 platform, where the internal
+representation uses a 64-bit fraction, the command:

 @example
-$ seq -s ' ' 0 0.000001 0.000003
-0.000000 0.000001 0.000002
+seq 1 0.0000000000000000001 1.0000000000000000009
 @end example

-But that doesn't happen on many systems because @command{seq} is
-implemented using binary floating point arithmetic (via the C
address@hidden double} type)---which means decimal fractions like 
@code{0.000001}
-cannot be represented exactly.  That in turn means some nonintuitive
-conditions like @address@hidden * 3 > 0.000003}} will end up being true.
-
-To work around that in the above example, use a slightly larger number as
-the @var{last} value:
-
address@hidden
-$ seq -s ' ' 0 0.000001 0.0000031
-0.000000 0.000001 0.000002 0.000003
address@hidden example
-
-In general, when using an @var{increment} with a fractional part, where
-(@var{last} - @var{first}) / @var{increment} is (mathematically) a whole
-number, specify a slightly larger (or smaller, if @var{increment} is negative)
-value for @var{last} to ensure that @var{last} is the final value printed
-by seq.
+outputs 1.0000000000000000007 twice and skips 1.0000000000000000008.

 @exitstatus

diff --git a/src/seq.c b/src/seq.c
index c59c6b5..de2a522 100644
--- a/src/seq.c
+++ b/src/seq.c
@@ -1,5 +1,5 @@
 /* seq - print sequence of numbers to standard output.
-   Copyright (C) 1994-2006 Free Software Foundation, Inc.
+   Copyright (C) 1994-2007 Free Software Foundation, Inc.

    This program is free software; you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
@@ -210,15 +210,53 @@ print_numbers (char const *fmt,
               long double first, long double step, long double last)
 {
   long double i;
+  long double x0 IF_LINT (= 0);

   for (i = 0; /* empty */; i++)
     {
       long double x = first + i * step;
+
       if (step < 0 ? x < last : last < x)
-       break;
+       {
+         /* If we go one past the end, but that number prints the
+            same way "last" does, and prints differently from the
+            previous number, then print "last".  This avoids problems
+            with rounding.  For example, with the x86 it causes "seq
+            0 0.000001 0.000003" to print 0.000003 instead of
+            stopping at 0.000002.  */
+
+         if (i)
+           {
+             char *x_str = NULL;
+             char *last_str = NULL;
+             if (asprintf (&x_str, fmt, x) < 0
+                 || asprintf (&last_str, fmt, last) < 0)
+               xalloc_die ();
+
+             if (strcmp (x_str, last_str) == 0)
+               {
+                 char *x0_str = NULL;
+                 if (asprintf (&x0_str, fmt, x0) < 0)
+                   xalloc_die ();
+                 if (strcmp (x0_str, x_str) != 0)
+                   {
+                     fputs (separator, stdout);
+                     fputs (x_str, stdout);
+                   }
+                 free (x0_str);
+               }
+
+             free (x_str);
+             free (last_str);
+           }
+
+         break;
+       }
+
       if (i)
        fputs (separator, stdout);
       printf (fmt, x);
+      x0 = x;
     }

   if (i)
diff --git a/tests/seq/basic b/tests/seq/basic
index 7070652..34a8c58 100755
--- a/tests/seq/basic
+++ b/tests/seq/basic
@@ -2,7 +2,8 @@
 # -*- perl -*-
 # Test "seq".

-# Copyright (C) 1999, 2000, 2003, 2005, 2006 Free Software Foundation, Inc.
+# Copyright (C) 1999, 2000, 2003, 2005, 2006, 2007 Free Software
+# Foundation, Inc.

 # This program is free software; you can redistribute it and/or modify
 # it under the terms of the GNU General Public License as published by
@@ -49,6 +50,13 @@ my @Tests =
    ['neg-3',   qw(1 -1 0),     {OUT => [qw(1 0)]}],
    ['neg-4',   qw(1 -1 -1),    {OUT => [qw(1 0 -1)]}],

+   ['float-1', qw(0.8 0.1 0.9),        {OUT => [qw(0.8 0.9)]}],
+   ['float-2', qw(0.1 0.99 1.99),      {OUT => [qw(0.10 1.09)]}],
+   ['float-3', qw(10.8 0.1 10.95),     {OUT => [qw(10.8 10.9)]}],
+   ['float-4', qw(0.1 -0.1 -0.2),      {OUT => [qw(0.1 0.0 -0.1 -0.2)]}],
+   ['float-5', qw(0.8 1e-1 0.9),       {OUT => [qw(0.8 0.9)]}],
+   ['float-6', qw(0.8 0.1 0.90000000000000000000),     {OUT => [qw(0.8 0.9)]}],
+
    ['eq-wid-1',        qw(-w 1 -1 -1), {OUT => [qw(01 00 -1)]}],

    # Prior to 2.0g, this test would fail on e.g., HPUX systems
M ChangeLog
M NEWS
M doc/ChangeLog
M doc/coreutils.texi
M src/seq.c
M tests/seq/basic
Committed as 927df10882228f17256e6f8e23b3b11d79ef6114