
Boost : 
Subject: Re: [boost] [odeint] Iterator semantics
From: Karsten Ahnert (karsten.ahnert_at_[hidden])
Date: 20120712 02:21:04
On 07/12/2012 12:39 AM, Jeffrey Lee Hellrung, Jr. wrote:
> On Wed, Jul 11, 2012 at 2:57 PM, Karsten Ahnert
> <karsten.ahnert_at_[hidden]>wrote:
>
>> On 07/11/2012 07:44 PM, Jeffrey Lee Hellrung, Jr. wrote:
>>> On Wed, Jul 11, 2012 at 9:04 AM, Karsten Ahnert
>>> <karsten.ahnert_at_[hidden]>wrote:
>>>
>>>> Hi,
>>>>
>>>> we have posted a review request for odeint (a library for solving
>>>> ordinary differential equations) on this list some days ago. An ODE is
>>>> dx/dt = f(x(t)) and one is looking for the solution x(t) which is a
>>>> function of t. Solving ODEs is usually an iterative task, so one starts
>>>> with x(0) and iteratively applies a solver to obtain a discretized (and
>>>> approximated) representation of the solution: x(0) > x(dt) > x(2*dt)
>>>> > x(3*dt) > ...
>>>>
>>>> I played around with an iterator which is doing exactly this iterative
>>>> procedure. For example, it can be used via
>>>>
>>>> odeint::runge_kutta4< state_type > stepper; // ode solver
>>>> std::array< double , 3 > x = {{ 10.0 , 10.0 , 10.0 }}; // initial state
>>>>
>>>> double res = std::accumulate(
>>>> // parameters : solver , ode , state , time , time_step
>>>> make_const_step_iterator_begin(stepper, lorenz(), x, 0.0, 0.01) ,
>>>> make_const_step_iterator_end( stepper, lorenz(), x, 1.0, 0.01) ,
>>>> 0.0 , []( double sum, const state_type &x) {
>>>> return sum + x[0]; } );
>>>>
>>>> The first iterator increments the time until the time of the second
>>>> iterator is rearched t=0.0 > t=0.01 > t=0.02 > .. t=1.0 and the first
>>>> component of the solution is accumulated.
>>>>
>>>> In principle, this approach works quite well, but there are some
>>>> problems at least semantically. Maybe someone here has some comments or
>>>> ideas on this:
>>>>
>>>> The defintion of this iterator is:
>>>>
>>>> template< class Stepper , class System , class StepperTag >
>>>> class const_step_iterator : public boost::iterator_facade
>>>> <
>>>> const_step_iterator< Stepper , System , StepperTag > ,
>>>> typename Stepper::state_type const ,
>>>> boost::single_pass_traversal_tag
>>>>>
>>>> {
>>>> ...
>>>> };
>>>>
>>>> It is a single pass iterator, so it can only check for equality of two
>>>> iterators.
>>
> [...]
>
>>> This is problematic since equality here means that the time
>>>> (t) of the begin iterator is smaller then the time of the end iterator.
>>>>
>>>
>>> Can you elaborate? I don't understand why this is what equality means.
>>
>> One needs to implement the != operator (or equal if you use
>> Boost.Iterator). A naive implementation of operator!= (or the equal() if
>> you use boost.iterator) could check that the current time of the ODE is
>> smaller then the time of other iterator:
>>
>> bool equal( iterator other )
>> {
>> return this>t < other>t;
>> }
>>
>
> :: confused :: I would think iterator1 op iterator2 is equivalent to time1
> op time2.
No, this is not the case iterator != iterator2 is equivalaent to time1 <
time2 . The first operation is commutative the second one not.
>
> Of course, this is not a check for equality but a check for less_then.
>> It will also destroy commutativity of the !=operation and will not work
>> in general. So you need somehow introduce which iterator is the first
>> and which one is the last. I have done this with a flag which says if
>> the iterator is the first and which the last such that the operation is now
>>
>> bool equal( iterator other )
>> {
>> return ( this>first ) ?
>> ( this>t < other>t ) :
>> ( this>t > other>t );
>> }
>>
>>>
>>>
>>>> To ensure commutativity of the != or == operation I needed to flag the
>>>> begin and end iterator explicitly. This is the reason for the two
>>>> factory functions make_const_step_iterator_begin and
>>>> make_const_step_iterator_end, which are doing this flagging.
>>>>
>>>
>>> I might understand your reasoning for doing this better if I understood
>>> your rationale for equality semantics...sorry, maybe I'm being dense :/
>>>
>>>
>>>> The second semantically problem is that the end iterator in principle
>>>> does not need to know the stepper as well as the system (lorenz() in the
>>>> above example). But all algorithms from the standard library and
>>>> Boost.Range assume that the begin and end iterator are of same type.
>>>> Therefore you have to put the stepper and the system into the end
>>>> iterator too.
>>>
>>>
>>> Ugh, yes, this is an unfortunate consequence of the design of
>>> iteratorbased, STLlike algorithms. Perhaps it would be better to simply
>>> define rangeequivalents of your iterators, and encourage use of the
>>> Boost.Range algorithms (which, AFAIK, under the hood, use the STL
>>> iteratorbased algorithms)? That would eliminate a lot of the repetition
>>> you're seeing in defining begin and end iterators separately.
>>
>> Yeah, we already have a factory function returning a range. You could
>> also write the upper example as follows:
>>
>> double res = boost::range::accumulate(
>> make_const_step_range(stepper, lorenz(), x, 0.0, 1.0 , 0.01 ) ,
>> 0.0 , []( double sum, const state_type &x) {
>> return sum + x[0]; } );
>>
>> But I fear, this does not solve this problem. Boost.Range also expects
>> that begin and end are of the same type.
>>
>
> Sounds like you want your iterator to be a wrapper around optional< state
>> . The state also encapsulates the end condition, and an end iterator is an
> uninitialized optional. When incrementing an iterator, if you've reached
> the end condition, uninitialize (i.e., reset) the underlying optional. ==
> and != on the iterators are equivalent to those on the optionals.
> Equivalently, if your state is completely external (and you don't need to
> manage it within the iterator), you can use a state* and use its nullness
> to distguish end iterators from notend iterators.
Yes, maybe optional might be a good idea. I will have a look at this.
[..]
Boost list run by bdawes at acm.org, gregod at cs.rpi.edu, cpdaniel at pacbell.net, john at johnmaddock.co.uk