On Tuesday, 24 December 2013 at 19:08:40 UTC, H. S. Teoh wrote:
> On Tue, Dec 24, 2013 at 06:57:50PM +0000, Craig Dillabaugh wrote:
>> On Tuesday, 24 December 2013 at 18:56:02 UTC, Craig Dillabaugh
>> wrote:
>> >On Tuesday, 24 December 2013 at 17:10:53 UTC, Andrei Alexandrescu
>> >wrote:
>> >>On 12/24/13 5:09 AM, Joseph Rushton Wakeling wrote:
>> clip
>> >
>> >Doesn't think work, or am I missing something?
>> >
>> >low + floor( (up-low)/step ) * step
>> 
>> I meant "Doesn't this work". Why can't I ever make a post without a
>> typo!
>
> You mean, a tpyo? :-P
>
>
> T

Exactly!

On 12/24/13 10:56 AM, Craig Dillabaugh wrote:
> On Tuesday, 24 December 2013 at 17:10:53 UTC, Andrei Alexandrescu
> wrote:
>> On 12/24/13 5:09 AM, Joseph Rushton Wakeling wrote:
>>> On 24/12/13 13:58, monarch_dodra wrote:
>>>> I think you are missing the point of what happens if the step is not 1
>>>> (or if
>>>> the passed in type can have fractional input). EG:
>>>>
>>>> iota(0, 105, 10);
>>>> or
>>>> iota(0, 10.5);
>>>>
>>>> In this case, "back" should be 100, and not 95. To compute back, you
>>>> need to be
>>>> able to evaluate length, and to add length*inc to front.
>>>
>>> Oh, snap.  Have we been working on the same problems for too long? :-)
>>
>> The integral cases are easy. We need to crack the floating point case:
>> given numbers low, up, and step, what's the closest number smaller
>> than up that's reached by repeated adds of step to low?
>>
>> Andrei
>
> Doesn't think work, or am I missing something?
>
> low + floor( (up-low)/step ) * step

I doubt it's anything as simple as that. The magnitudes of up, low, and step must be taken into account.

Andrei

On 12/24/13 11:17 AM, Craig Dillabaugh wrote:
> On Tuesday, 24 December 2013 at 19:08:40 UTC, H. S. Teoh wrote:
>> On Tue, Dec 24, 2013 at 06:57:50PM +0000, Craig Dillabaugh wrote:
>>> On Tuesday, 24 December 2013 at 18:56:02 UTC, Craig Dillabaugh
>>> wrote:
>>> >On Tuesday, 24 December 2013 at 17:10:53 UTC, Andrei >Alexandrescu
>>> >wrote:
>>> >>On 12/24/13 5:09 AM, Joseph Rushton Wakeling wrote:
>>> clip
>>> >
>>> >Doesn't think work, or am I missing something?
>>> >
>>> >low + floor( (up-low)/step ) * step
>>>
>>> I meant "Doesn't this work". Why can't I ever make a post without a
>>> typo!
>>
>> You mean, a tpyo? :-P
>>
>>
>> T
>
> Exactly!

Use a spellchequer.

Andrei

On Tue, Dec 24, 2013 at 11:35:49AM -0800, Andrei Alexandrescu wrote:
> On 12/24/13 10:56 AM, Craig Dillabaugh wrote:
> >On Tuesday, 24 December 2013 at 17:10:53 UTC, Andrei Alexandrescu wrote:
[.[..]
> >>The integral cases are easy. We need to crack the floating point case: given numbers low, up, and step, what's the closest number smaller than up that's reached by repeated adds of step to low?
> >>
> >>Andrei
> >
> >Doesn't think work, or am I missing something?
> >
> >low + floor( (up-low)/step ) * step
> 
> I doubt it's anything as simple as that. The magnitudes of up, low, and step must be taken into account.
[...]

What about low + fmod(up-low, step)*step? Is that better? (Assuming that
fmod would take relative magnitudes into account, of course.)


T

-- 
Only boring people get bored. -- JM

On 12/24/13 11:59 AM, H. S. Teoh wrote:
> On Tue, Dec 24, 2013 at 11:35:49AM -0800, Andrei Alexandrescu wrote:
>> On 12/24/13 10:56 AM, Craig Dillabaugh wrote:
>>> On Tuesday, 24 December 2013 at 17:10:53 UTC, Andrei Alexandrescu
>>> wrote:
> [.[..]
>>>> The integral cases are easy. We need to crack the floating point
>>>> case: given numbers low, up, and step, what's the closest number
>>>> smaller than up that's reached by repeated adds of step to low?
>>>>
>>>> Andrei
>>>
>>> Doesn't think work, or am I missing something?
>>>
>>> low + floor( (up-low)/step ) * step
>>
>> I doubt it's anything as simple as that. The magnitudes of up, low,
>> and step must be taken into account.
> [...]
>
> What about low + fmod(up-low, step)*step? Is that better? (Assuming that
> fmod would take relative magnitudes into account, of course.)

No, again, I think nothing like that would work. It's hopelessly naive (in addition to being plain wrong for simple inputs like low=1, high=10, step=1).

There are combinations of low, high, and step that don't even make progress, i.e. step is sufficiently small compared to low to effectively make low + step == low.

Try this to verify approaches:

#!/Users/aalexandre/bin/rdmd
import std.stdio, std.conv, std.math;

void main(string[] args)
{
    auto low = to!double(args[1]);
    auto up = to!double(args[2]);
    auto step = to!double(args[3]);
    writeln("Predicted: ", low + fmod(up-low, step)*step);
    for (;;)
    {
        auto next = low + step;
        if (next >= up) break;
        low = next;
    }
    writeln("Actual: ", low);
}


Andrei

On Tue, Dec 24, 2013 at 12:57:02PM -0800, Andrei Alexandrescu wrote:
> On 12/24/13 11:59 AM, H. S. Teoh wrote:
> >On Tue, Dec 24, 2013 at 11:35:49AM -0800, Andrei Alexandrescu wrote:
> >>On 12/24/13 10:56 AM, Craig Dillabaugh wrote:
> >>>On Tuesday, 24 December 2013 at 17:10:53 UTC, Andrei Alexandrescu wrote:
> >[...]
> >>>>The integral cases are easy. We need to crack the floating point case: given numbers low, up, and step, what's the closest number smaller than up that's reached by repeated adds of step to low?
> >>>>
> >>>>Andrei
> >>>
> >>>Doesn't think work, or am I missing something?
> >>>
> >>>low + floor( (up-low)/step ) * step
> >>
> >>I doubt it's anything as simple as that. The magnitudes of up, low, and step must be taken into account.
> >[...]
> >
> >What about low + fmod(up-low, step)*step? Is that better? (Assuming
> >that fmod would take relative magnitudes into account, of course.)
> 
> No, again, I think nothing like that would work. It's hopelessly naive (in addition to being plain wrong for simple inputs like low=1, high=10, step=1).
> 
> There are combinations of low, high, and step that don't even make progress, i.e. step is sufficiently small compared to low to effectively make low + step == low.

Then what should be returned in that case?


> Try this to verify approaches:
> 
> #!/Users/aalexandre/bin/rdmd
> import std.stdio, std.conv, std.math;
> 
> void main(string[] args)
> {
>     auto low = to!double(args[1]);
>     auto up = to!double(args[2]);
>     auto step = to!double(args[3]);
>     writeln("Predicted: ", low + fmod(up-low, step)*step);
>     for (;;)
>     {
>         auto next = low + step;

Aren't you accumulating roundoff errors here? Since + may introduce an error of 1/2 ulp, if your loop runs in n steps, wouldn't it accumulate an error of n/2 ulps? For large n, this would significantly skew your results.


>         if (next >= up) break;
>         low = next;
>     }
>     writeln("Actual: ", low);
> }
[...]


T

-- 
All problems are easy in retrospect.

On 12/24/13 1:09 PM, H. S. Teoh wrote:
> On Tue, Dec 24, 2013 at 12:57:02PM -0800, Andrei Alexandrescu wrote:
>> On 12/24/13 11:59 AM, H. S. Teoh wrote:
>>> On Tue, Dec 24, 2013 at 11:35:49AM -0800, Andrei Alexandrescu wrote:
>>>> On 12/24/13 10:56 AM, Craig Dillabaugh wrote:
>>>>> On Tuesday, 24 December 2013 at 17:10:53 UTC, Andrei Alexandrescu
>>>>> wrote:
>>> [...]
>>>>>> The integral cases are easy. We need to crack the floating point
>>>>>> case: given numbers low, up, and step, what's the closest number
>>>>>> smaller than up that's reached by repeated adds of step to low?
>>>>>>
>>>>>> Andrei
>>>>>
>>>>> Doesn't think work, or am I missing something?
>>>>>
>>>>> low + floor( (up-low)/step ) * step
>>>>
>>>> I doubt it's anything as simple as that. The magnitudes of up, low,
>>>> and step must be taken into account.
>>> [...]
>>>
>>> What about low + fmod(up-low, step)*step? Is that better? (Assuming
>>> that fmod would take relative magnitudes into account, of course.)
>>
>> No, again, I think nothing like that would work. It's hopelessly
>> naive (in addition to being plain wrong for simple inputs like
>> low=1, high=10, step=1).
>>
>> There are combinations of low, high, and step that don't even make
>> progress, i.e. step is sufficiently small compared to low to
>> effectively make low + step == low.
>
> Then what should be returned in that case?

The same as an iota with step zero.

>>          auto next = low + step;
>
> Aren't you accumulating roundoff errors here?

Ideally one addition is all needed for popBack.


Andrei

On Tue, Dec 24, 2013 at 01:18:34PM -0800, Andrei Alexandrescu wrote:
> On 12/24/13 1:09 PM, H. S. Teoh wrote:
> >On Tue, Dec 24, 2013 at 12:57:02PM -0800, Andrei Alexandrescu wrote:
> >>On 12/24/13 11:59 AM, H. S. Teoh wrote:
> >>>On Tue, Dec 24, 2013 at 11:35:49AM -0800, Andrei Alexandrescu wrote:
> >>>>On 12/24/13 10:56 AM, Craig Dillabaugh wrote:
> >>>>>On Tuesday, 24 December 2013 at 17:10:53 UTC, Andrei Alexandrescu wrote:
> >>>[...]
> >>>>>>The integral cases are easy. We need to crack the floating point case: given numbers low, up, and step, what's the closest number smaller than up that's reached by repeated adds of step to low?
> >>>>>>
> >>>>>>Andrei
> >>>>>
> >>>>>Doesn't think work, or am I missing something?
> >>>>>
> >>>>>low + floor( (up-low)/step ) * step
> >>>>
> >>>>I doubt it's anything as simple as that. The magnitudes of up, low, and step must be taken into account.
> >>>[...]
> >>>
> >>>What about low + fmod(up-low, step)*step? Is that better? (Assuming
> >>>that fmod would take relative magnitudes into account, of course.)
> >>
> >>No, again, I think nothing like that would work. It's hopelessly naive (in addition to being plain wrong for simple inputs like low=1, high=10, step=1).
> >>
> >>There are combinations of low, high, and step that don't even make progress, i.e. step is sufficiently small compared to low to effectively make low + step == low.
> >
> >Then what should be returned in that case?
> 
> The same as an iota with step zero.
> 
> >>         auto next = low + step;
> >
> >Aren't you accumulating roundoff errors here?
> 
> Ideally one addition is all needed for popBack.
[...]

You're missing my point. I'm not talking about popBack specifically here. I'm talking about the problem of accumulated roundoff error in the implementation of iota. You're suggesting that iota should keep a running accumulator representing the current value of .front, and popFront should simply add step to this accumulator.

But this kind of implementation is problematic, because it accumulates roundoff errors. Suppose evaluating current + step introduces a roundoff error of E. Then the next time we call popFront, the error becomes 2*E. Then the third time the error becomes 3*E.  So after calling popFront n times, the accumulator is now off from the correct value by n*E. Since n is increasing, eventually n*E > step, which will cause iota to return the wrong number of elements in the range (not to mention that the last few elements will be completely inaccurate).

Yet you're suggesting to use this method of incrementally adding step to an accumulator as the standard by which to verify the correctness of .back? That makes no sense.

The only floating-point implementation of iota that guarantees the correct number of iterations (up to floating-point precision, of course) is one that *directly* calculates the value of front_i = low + i*step, instead of consecutively adding step to an accumulator, no matter whether you're calling popFront or popBack or opIndex.


T

-- 
The fact that anyone still uses AOL shows that even the presence of options doesn't stop some people from picking the pessimal one. - Mike Ellis

On 12/24/13 3:11 PM, H. S. Teoh wrote:
> You're missing my point.

It's Christmas, let's be nice to one another :o). The only problem here is that I didn't explain my point enough, which fostered confusion.

> I'm not talking about popBack specifically
> here. I'm talking about the problem of accumulated roundoff error in the
> implementation of iota. You're suggesting that iota should keep a
> running accumulator representing the current value of .front, and
> popFront should simply add step to this accumulator.
>
> But this kind of implementation is problematic, because it accumulates
> roundoff errors.

No need to explain. I understand all of this already.

> Yet you're suggesting to use this method of incrementally adding step to
> an accumulator as the standard by which to verify the correctness of
> .back? That makes no sense.

Of course it does. It's just that it's quite a difficult problem. The current implementation is conservative at the expense of speed.

> The only floating-point implementation of iota that guarantees the
> correct number of iterations (up to floating-point precision, of course)
> is one that *directly* calculates the value of front_i = low + i*step,
> instead of consecutively adding step to an accumulator, no matter
> whether you're calling popFront or popBack or opIndex.

No.

Consider:

struct IotaDouble
{
private:
    double begin, end, step, current;
    size_t index;
    enum UpdateMethod { indexed, addition, increment }
    UpdateMethod m;
    UpdateMethod initUpdate()
    {
        ... this is the difficult part ...
    }
    ...
    void popFront()
    {
        final switch (m)
        {
        case UpdateMethod.indexed:
            current = begin + ++index * step;
            break;
        case UpdateMethod.addition:
            current += step;
            break;
        case UpdateMethod.increment:
            ++current;
            break;
        }
    }
}

I suspect increment is not much faster than addition, so they could be merged subject to measurement. But I think the addition is quite a bit faster than indexing, which makes it worth special casing for. Also, the update method will be invariant through the lifetime of the object, which will play well with the branch predictor (the test may come at virtually or literally no cost).

The difficulty, of course, is choosing the appropriate update method depending on the limits and step.


Andrei

Francesco Cattoglio:

> Sorry for the amount of text, but I tried to explain everything in a clear and simple way.
> I'm willing to volunteer for adding support for more types in iota.

Probably directly addressed by issue 10762:

Problem with iota(long)
https://d.puremagic.com/issues/show_bug.cgi?id=6446

iota(BigInt) too
https://d.puremagic.com/issues/show_bug.cgi?id=6447

assertion failure in std.range.iota
https://d.puremagic.com/issues/show_bug.cgi?id=6531

iota should generate char intervals too
https://d.puremagic.com/issues/show_bug.cgi?id=9447


But if you rewrite iota() also take in account the following:

Purity:

Make std.range.iota strongly pure
https://d.puremagic.com/issues/show_bug.cgi?id=5746

map, filter, iota and zip in pure (and nothrow) functions
https://d.puremagic.com/issues/show_bug.cgi?id=8882


An important enhancement:

Optional "[]" syntax for std.range.iota too
https://d.puremagic.com/issues/show_bug.cgi?id=10466

(A less important enhancement is issue 11252).

Bye,
bearophile