[odeint] Iterator semantics

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. 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. 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. Of course, there are techniques like type erasure which might solve this problem, but they introduce some run-time overhead which we want to avoid here. Any ideas or comments about this? P.S. The code can be found in the iterator branch in odeint: boost/numeric/odeint/iterator/const_step_iterator.hpp and some examples libs/numeric/odeint/examples/const_step_iterator.cpp git clone -b iterator git@github.com:headmyshoulder/odeint-v2.git

On Wed, Jul 11, 2012 at 9:04 AM, Karsten Ahnert <karsten.ahnert@ambrosys.de>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.
It could make sense to define the remaining relational operators yourself, for convenience, but their existence in standard algorithms wouldn't be picked up due to the iterator's traversal category. 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.
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 iterator-based, STL-like algorithms. Perhaps it would be better to simply define range-equivalents of your iterators, and encourage use of the Boost.Range algorithms (which, AFAIK, under the hood, use the STL iterator-based algorithms)? That would eliminate a lot of the repetition you're seeing in defining begin and end iterators separately. Of course, there are techniques like type erasure which
might solve this problem, but they introduce some run-time overhead which we want to avoid here.
I don't see type erasure solving your problem, but regarding the overhead of type erasure only, it could be minimal relative to the cost of the actual time stepper. But that observation probably isn't relevant to your problem.
Any ideas or comments about this?
HTH, - Jeff

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@ambrosys.de>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.
It could make sense to define the remaining relational operators yourself, for convenience, but their existence in standard algorithms wouldn't be picked up due to the iterator's traversal category.
Yes one can do this, but it will not be used, only in own algorithms...
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; } 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 iterator-based, STL-like algorithms. Perhaps it would be better to simply define range-equivalents of your iterators, and encourage use of the Boost.Range algorithms (which, AFAIK, under the hood, use the STL iterator-based 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.
Of course, there are techniques like type erasure which
might solve this problem, but they introduce some run-time overhead which we want to avoid here.
I don't see type erasure solving your problem, but regarding the overhead of type erasure only, it could be minimal relative to the cost of the actual time stepper. But that observation probably isn't relevant to your problem.
With type erasure it is in principle possible to remove the dependency of iterator from the stepper and the system and put this in some member class, maybe something like template< class State , class StepperTag > class const_step_iterator : public boost::iterator_facade < const_step_iterator< State , StepperTag > , State const , boost::single_pass_traversal_tag
{ public: template< class Stepper , class System > const_step_iterator( Stepper stepper , System system , State &x , Time t_start , Time dt ) { this->data_ = new DataStepper< Stepper , System >( stepper , system , State , t , dt ); } template< class Stepper , class System > const_step_iterator( Time t_end ) { this->data_ = new DataTime( t ); } struct IData { ... }; template< class Stepper , class System > struct DataStepper : public IData { ... }; struct DataTime : public IData { ... }; ... private: std::unique_ptr< IData > data_; };

On Wed, Jul 11, 2012 at 2:57 PM, Karsten Ahnert <karsten.ahnert@ambrosys.de>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@ambrosys.de>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. 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 iterator-based, STL-like algorithms. Perhaps it would be better to simply define range-equivalents of your iterators, and encourage use of the Boost.Range algorithms (which, AFAIK, under the hood, use the STL iterator-based 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.
. The state also encapsulates the end condition, and an end iterator is an uninitialized optional. When incrementing an iterator, if you've reached
Sounds like you want your iterator to be a wrapper around optional< state 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 not-end iterators.
Of course, there are techniques like type erasure which
might solve this problem, but they introduce some run-time overhead which we want to avoid here.
I don't see type erasure solving your problem, but regarding the overhead of type erasure only, it could be minimal relative to the cost of the actual time stepper. But that observation probably isn't relevant to your problem.
With type erasure it is in principle possible to remove the dependency of iterator from the stepper and the system and put this in some member class, maybe something like
template< class State , class StepperTag > class const_step_iterator : public boost::iterator_facade < const_step_iterator< State , StepperTag > , State const , boost::single_pass_traversal_tag
{ public:
template< class Stepper , class System > const_step_iterator( Stepper stepper , System system , State &x , Time t_start , Time dt ) { this->data_ = new DataStepper< Stepper , System >( stepper , system , State , t , dt ); }
template< class Stepper , class System > const_step_iterator( Time t_end ) { this->data_ = new DataTime( t ); }
struct IData { ... };
template< class Stepper , class System > struct DataStepper : public IData { ... };
struct DataTime : public IData { ... };
...
private:
std::unique_ptr< IData > data_; };
I don't see what this buys you relative to the original problem. - Jeff

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@ambrosys.de>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@ambrosys.de>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 iterator-based, STL-like algorithms. Perhaps it would be better to simply define range-equivalents of your iterators, and encourage use of the Boost.Range algorithms (which, AFAIK, under the hood, use the STL iterator-based 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.
. The state also encapsulates the end condition, and an end iterator is an uninitialized optional. When incrementing an iterator, if you've reached
Sounds like you want your iterator to be a wrapper around optional< state 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 not-end iterators.
Yes, maybe optional might be a good idea. I will have a look at this. [..]

On Wed, Jul 11, 2012 at 11:21 PM, Karsten Ahnert <karsten.ahnert@ambrosys.de
wrote:
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@ambrosys.de>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@ambrosys.de>wrote:
[...]
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.
Well, you can define operator== and operator!= however you want, and some definitions will make sense and some won't (e.g., non-/anti-symmetry). *All* I've been trying to do here is coax out of you *why* you're defining operator== and operator!= this way :) (Especially, with my limited understanding of the application, when it seems time1 op time2 would do what you want, or perhaps a fuzzy version taking into account the time step.) [...] - Jeff

On Wed, 2012-07-11 at 23:50 -0700, Jeffrey Lee Hellrung, Jr. wrote:
On Wed, Jul 11, 2012 at 11:21 PM, Karsten Ahnert <karsten.ahnert@ambrosys.de
wrote:
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@ambrosys.de>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@ambrosys.de>wrote:
[...]
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.
Well, you can define operator== and operator!= however you want, and some definitions will make sense and some won't (e.g., non-/anti-symmetry). *All* I've been trying to do here is coax out of you *why* you're defining operator== and operator!= this way :)
(Especially, with my limited understanding of the application, when it seems time1 op time2 would do what you want, or perhaps a fuzzy version taking into account the time step.)
I agree that it1 != it2 should be implemented in the sense of t1+-dt/2 != t2+-dt/2 As far as I understand the timestep dt should be available at this point. Unfortunately, this behavior would be different from the current integrate implementations that guarantee that the endpoint is reached while in the above version we might stop in the interval t_end-dt/2. Nevertheless I think this is still better than flagging iterators. Best, Mario

On 07/12/2012 09:52 AM, Mario Mulansky wrote:
On Wed, 2012-07-11 at 23:50 -0700, Jeffrey Lee Hellrung, Jr. wrote:
On Wed, Jul 11, 2012 at 11:21 PM, Karsten Ahnert <karsten.ahnert@ambrosys.de
wrote:
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@ambrosys.de>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@ambrosys.de>wrote:
[...]
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.
Well, you can define operator== and operator!= however you want, and some definitions will make sense and some won't (e.g., non-/anti-symmetry). *All* I've been trying to do here is coax out of you *why* you're defining operator== and operator!= this way :)
(Especially, with my limited understanding of the application, when it seems time1 op time2 would do what you want, or perhaps a fuzzy version taking into account the time step.)
I agree that it1 != it2 should be implemented in the sense of t1+-dt/2 != t2+-dt/2
So you are checking for overlap. This looks good and I think it should work.
As far as I understand the timestep dt should be available at this point. Unfortunately, this behavior would be different from the current integrate implementations that guarantee that the endpoint is reached while in the above version we might stop in the interval t_end-dt/2.
Maybe this is not so problematic. The iterators can be different from the integrate routines.
Nevertheless I think this is still better than flagging iterators.

on Thu Jul 12 2012, Karsten Ahnert <karsten.ahnert-AT-ambrosys.de> wrote:
On 07/12/2012 09:52 AM, Mario Mulansky wrote:
On Wed, 2012-07-11 at 23:50 -0700, Jeffrey Lee Hellrung, Jr. wrote:
On Wed, Jul 11, 2012 at 11:21 PM, Karsten Ahnert <karsten.ahnert@ambrosys.de
wrote:
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@ambrosys.de>wrote:
Well, you can define operator== and operator!= however you want, and some definitions will make sense and some won't (e.g., non-/anti-symmetry). *All* I've been trying to do here is coax out of you *why* you're defining operator== and operator!= this way :)
(Especially, with my limited understanding of the application, when it seems time1 op time2 would do what you want, or perhaps a fuzzy version taking into account the time step.)
I agree that it1 != it2 should be implemented in the sense of t1+-dt/2 != t2+-dt/2
So you are checking for overlap. This looks good and I think it should work.
If I'm understanding correctly what you're saying, making two iterators equal if their underlying time ranges overlap would be "evil." Your equality operator would not be transitive, and thus not an equivalence relation, and therefore you would not have implemented a conforming iterator, and you'd have no right to expect an algorithm that works on iterators to work on yours. But maybe in my fast skimming I've failed to grasp what's going on. -- Dave Abrahams BoostPro Computing Software Development Training http://www.boostpro.com Clang/LLVM/EDG Compilers C++ Boost

On 07/29/2012 04:31 AM, Dave Abrahams wrote:
on Thu Jul 12 2012, Karsten Ahnert <karsten.ahnert-AT-ambrosys.de> wrote:
On 07/12/2012 09:52 AM, Mario Mulansky wrote:
On Wed, 2012-07-11 at 23:50 -0700, Jeffrey Lee Hellrung, Jr. wrote:
On Wed, Jul 11, 2012 at 11:21 PM, Karsten Ahnert <karsten.ahnert@ambrosys.de
wrote:
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@ambrosys.de>wrote:
Well, you can define operator== and operator!= however you want, and some definitions will make sense and some won't (e.g., non-/anti-symmetry). *All* I've been trying to do here is coax out of you *why* you're defining operator== and operator!= this way :)
(Especially, with my limited understanding of the application, when it seems time1 op time2 would do what you want, or perhaps a fuzzy version taking into account the time step.)
I agree that it1 != it2 should be implemented in the sense of t1+-dt/2 != t2+-dt/2
So you are checking for overlap. This looks good and I think it should work.
If I'm understanding correctly what you're saying, making two iterators equal if their underlying time ranges overlap would be "evil." Your equality operator would not be transitive, and thus not an equivalence relation, and therefore you would not have implemented a conforming iterator, and you'd have no right to expect an algorithm that works on iterators to work on yours.
Ahh, I did not thought about transitivity. You are right, checking for overlap destroys transitivity. Nevertheless, I came to point that checking for the relation (iter1.time < iter2.time) in combination with a flagging which iterator is the first is sufficient.

On 29.07.2012 6:31, Dave Abrahams wrote:
on Thu Jul 12 2012, Karsten Ahnert <karsten.ahnert-AT-ambrosys.de> wrote:
On 07/12/2012 09:52 AM, Mario Mulansky wrote:
On Wed, 2012-07-11 at 23:50 -0700, Jeffrey Lee Hellrung, Jr. wrote:
I agree that it1 != it2 should be implemented in the sense of t1+-dt/2 != t2+-dt/2
So you are checking for overlap. This looks good and I think it should work.
If I'm understanding correctly what you're saying, making two iterators equal if their underlying time ranges overlap would be "evil." Your equality operator would not be transitive, and thus not an equivalence relation, and therefore you would not have implemented a conforming iterator, and you'd have no right to expect an algorithm that works on iterators to work on yours.
IMHO the values of the iterator's underlying time variable would belong to "almost" discrete {t0 + n*dt +/- epsilon}, n from Z, where epsilon comes from floating point errors and is small. Thus, the distance between neighboring values should be no less than (dt-2*epsilon). I would say a test 'distance(t1,t2) < dt/2' would be transitive if epsilon is less than dt/4 That's of course true while t0 and dt is the same for all iterators passed to an algorithm.
But maybe in my fast skimming I've failed to grasp what's going on.
-- ------------ Sergey Mitsyn.

On Tuesday, July 31, 2012 11:47:55 am Sergey Mitsyn wrote:
On 29.07.2012 6:31, Dave Abrahams wrote:
on Thu Jul 12 2012, Karsten Ahnert <karsten.ahnert-AT-ambrosys.de> wrote:
On 07/12/2012 09:52 AM, Mario Mulansky wrote:
On Wed, 2012-07-11 at 23:50 -0700, Jeffrey Lee Hellrung, Jr. wrote:
I agree that it1 != it2 should be implemented in the sense of t1+-dt/2 != t2+-dt/2
So you are checking for overlap. This looks good and I think it should work.
If I'm understanding correctly what you're saying, making two iterators equal if their underlying time ranges overlap would be "evil." Your equality operator would not be transitive, and thus not an equivalence relation, and therefore you would not have implemented a conforming iterator, and you'd have no right to expect an algorithm that works on iterators to work on yours.
IMHO the values of the iterator's underlying time variable would belong to "almost" discrete {t0 + n*dt +/- epsilon}, n from Z, where epsilon comes from floating point errors and is small. Thus, the distance between neighboring values should be no less than (dt-2*epsilon).> I would say a test 'distance(t1,t2) < dt/2' would be transitive if epsilon is less than dt/4
That's of course true while t0 and dt is the same for all iterators passed to an algorithm.
dt even changes while iterating when methods with step-size control are used. This case is really ugly and I'm not sure if it's possible to implement a fully compatible iterator there. Maybe one should restrict to fixed-step-size methods where the overlap checking works?

On 07/31/2012 03:07 PM, Mario Mulansky wrote:
On Tuesday, July 31, 2012 11:47:55 am Sergey Mitsyn wrote:
On 29.07.2012 6:31, Dave Abrahams wrote:
on Thu Jul 12 2012, Karsten Ahnert <karsten.ahnert-AT-ambrosys.de> wrote:
On 07/12/2012 09:52 AM, Mario Mulansky wrote:
On Wed, 2012-07-11 at 23:50 -0700, Jeffrey Lee Hellrung, Jr. wrote:
I agree that it1 != it2 should be implemented in the sense of t1+-dt/2 != t2+-dt/2
So you are checking for overlap. This looks good and I think it should work.
If I'm understanding correctly what you're saying, making two iterators equal if their underlying time ranges overlap would be "evil." Your equality operator would not be transitive, and thus not an equivalence relation, and therefore you would not have implemented a conforming iterator, and you'd have no right to expect an algorithm that works on iterators to work on yours.
IMHO the values of the iterator's underlying time variable would belong to "almost" discrete {t0 + n*dt +/- epsilon}, n from Z, where epsilon comes from floating point errors and is small. Thus, the distance between neighboring values should be no less than (dt-2*epsilon).> I would say a test 'distance(t1,t2) < dt/2' would be transitive if epsilon is less than dt/4
That's of course true while t0 and dt is the same for all iterators passed to an algorithm.
dt even changes while iterating when methods with step-size control are used. This case is really ugly and I'm not sure if it's possible to implement a fully compatible iterator there. Maybe one should restrict to fixed-step-size methods where the overlap checking works?
I think already found a solution. The end iterator is marked as end (it has a flag m_end). Then, two iterators are equal if theirs flags are equal, or if the flags is not equal the time of the first iterator must be smaller then the time of the end iterator: bool equal( iterator other ) const { if( m_end == other.m_end ) return true; else { return ( m_end ) ? ( other.m_time < m_time ) : ( m_time < other.m_time ); } } This iterator is a single pass iterator, for example it is similar to the istream_iterator.

So you are checking for overlap. This looks good and I think it should work.
If I'm understanding correctly what you're saying, making two iterators equal if their underlying time ranges overlap would be "evil." Your equality operator would not be transitive, and thus not an equivalence relation, and therefore you would not have implemented a conforming iterator, and you'd have no right to expect an algorithm that works on iterators to work on yours.
[snip]
I think already found a solution. The end iterator is marked as end (it has a flag m_end). Then, two iterators are equal if theirs flags are equal, or if the flags is not equal the time of the first iterator must be smaller then the time of the end iterator:
bool equal( iterator other ) const { if( m_end == other.m_end ) return true; else { return ( m_end ) ? ( other.m_time < m_time ) : ( m_time < other.m_time ); } }
This iterator is a single pass iterator, for example it is similar to the istream_iterator.
I analyzed the problem once more and it turned out that the iterator is not transitive. The problem is the following: a range consist of a begin and an end iterator. The end iterator is flagged as end is not dereferencable and can not be incremented. When two iterator of the same type (end-end, begin-begin) are compared to each other they are equal (neither of both should be incremented). If iterators of different types are compared it is checked that the time of the begin iterator is greater then the time of the end iterator. The violation of transitivity might now occur if two end iterators and one begin iterator are compared. The two end iterators are always equal, which is not true for the comparison of the begin iterator with the end iterators. I tried to find different solutions but it is difficult since one is more or less always checking for overlap. On the other hand I wonder if the requirement of transitivity can be relaxed in this case. I can not imagine an algorithm where the above situation might occur. The iterator is a single pass iterator which limits the number of possible algorithms. Furthermore the violation can only appear if two end iterators are used. This is not the case for all algorithms in the standard and Boost.Range. Any ideas or comments? Is it really necessary to require transitivity?

on Fri Aug 10 2012, Karsten Ahnert <karsten.ahnert-AT-ambrosys.de> wrote:
I tried to find different solutions but it is difficult since one is more or less always checking for overlap. On the other hand I wonder if the requirement of transitivity can be relaxed in this case. I can not imagine an algorithm where the above situation might occur. The iterator is a single pass iterator which limits the number of possible algorithms. Furthermore the violation can only appear if two end iterators are used. This is not the case for all algorithms in the standard and Boost.Range.
Any ideas or comments? Is it really necessary to require transitivity?
Yes. Transitivity of == is fundamental. -- Dave Abrahams BoostPro Computing Software Development Training http://www.boostpro.com Clang/LLVM/EDG Compilers C++ Boost

On 08/19/2012 02:57 AM, Dave Abrahams wrote:
on Fri Aug 10 2012, Karsten Ahnert <karsten.ahnert-AT-ambrosys.de> wrote:
I tried to find different solutions but it is difficult since one is more or less always checking for overlap. On the other hand I wonder if the requirement of transitivity can be relaxed in this case. I can not imagine an algorithm where the above situation might occur. The iterator is a single pass iterator which limits the number of possible algorithms. Furthermore the violation can only appear if two end iterators are used. This is not the case for all algorithms in the standard and Boost.Range.
Any ideas or comments? Is it really necessary to require transitivity?
Yes. Transitivity of == is fundamental.
Ok. Then, the question is which kind of algorithms will fail. I think all algorithms where you pass exactly one end-iterator are ok, like algorithm( first , last , ... ) You simply can not create an iterator for which transitivity which will in such situations. This is the case for most algorithm in the STL and boost.range.

On 07/31/2012 11:47 AM, Sergey Mitsyn wrote:
On 29.07.2012 6:31, Dave Abrahams wrote:
on Thu Jul 12 2012, Karsten Ahnert <karsten.ahnert-AT-ambrosys.de> wrote:
On 07/12/2012 09:52 AM, Mario Mulansky wrote:
On Wed, 2012-07-11 at 23:50 -0700, Jeffrey Lee Hellrung, Jr. wrote:
I agree that it1 != it2 should be implemented in the sense of t1+-dt/2 != t2+-dt/2
So you are checking for overlap. This looks good and I think it should work.
If I'm understanding correctly what you're saying, making two iterators equal if their underlying time ranges overlap would be "evil." Your equality operator would not be transitive, and thus not an equivalence relation, and therefore you would not have implemented a conforming iterator, and you'd have no right to expect an algorithm that works on iterators to work on yours.
IMHO the values of the iterator's underlying time variable would belong to "almost" discrete {t0 + n*dt +/- epsilon}, n from Z, where epsilon comes from floating point errors and is small. Thus, the distance between neighboring values should be no less than (dt-2*epsilon).
Shouldn't this be t0 + n*dt +/- n*epsilon? I think the error scales with N, or even sqrt(N)?
I would say a test 'distance(t1,t2) < dt/2' would be transitive if epsilon is less than dt/4
That's of course true while t0 and dt is the same for all iterators passed to an algorithm.
But maybe in my fast skimming I've failed to grasp what's going on.

on Wed Jul 11 2012, Karsten Ahnert <karsten.ahnert-AT-ambrosys.de> wrote:
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.
Any ideas or comments about this?
Known issue. I don't know of any truly good answers other than range-based algorithms that don't require the existence of iterators. In the meantime, there's Boost.Optional. -- Dave Abrahams BoostPro Computing Software Development Training http://www.boostpro.com Clang/LLVM/EDG Compilers C++ Boost

On 07/29/2012 04:27 AM, Dave Abrahams wrote:
on Wed Jul 11 2012, Karsten Ahnert <karsten.ahnert-AT-ambrosys.de> wrote:
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.
Any ideas or comments about this?
Known issue. I don't know of any truly good answers other than range-based algorithms that don't require the existence of iterators.
Do you mean the algorithms from Boost.Range?
In the meantime, there's Boost.Optional.
I will check Boost.Optional but I fear it is not possible to use, since it requires the knowledge of the stepper type in advance.

on Mon Jul 30 2012, Karsten Ahnert <karsten.ahnert-AT-ambrosys.de> wrote:
On 07/29/2012 04:27 AM, Dave Abrahams wrote:
on Wed Jul 11 2012, Karsten Ahnert <karsten.ahnert-AT-ambrosys.de> wrote:
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.
Any ideas or comments about this?
Known issue. I don't know of any truly good answers other than range-based algorithms that don't require the existence of iterators.
Do you mean the algorithms from Boost.Range?
I mean algorithms that strictly use range operations c.f. "Iterators Must Go" by Alexandrescu (google it).
In the meantime, there's Boost.Optional.
I will check Boost.Optional but I fear it is not possible to use, since it requires the knowledge of the stepper type in advance.
I don't understand your problem deeply, but you can use type erasure within an optional if knowing types in advance is an issue. -- Dave Abrahams BoostPro Computing Software Development Training http://www.boostpro.com Clang/LLVM/EDG Compilers C++ Boost
participants (5)
-
Dave Abrahams
-
Jeffrey Lee Hellrung, Jr.
-
Karsten Ahnert
-
Mario Mulansky
-
Sergey Mitsyn