odeint (c++) - downsample observations

199 views Asked by At

Sorry if this is a simple question - but is there a "best practice" for downsampling the evolution of the state variables in odeint?

Below, I've copied a nice example for building an "observer" to log the state variables provided in this article (http://www.codeproject.com/Articles/268589/odeint-v2-Solving-ordinary-differential-equations)

struct streaming_observer
{
     std::ostream &m_out;

     streaming_observer( std::ostream &out ) : m_out( out ) {}

     void operator()( const state_type &x , double t ) const
     {
          m_out << t;
          for( size_t i=0 ; i < x.size() ; ++i )
              m_out << "\t" << x[i];
          m_out << "\n";
     }
};

// ...

integrate_const( runge_kutta4< state_type >() , lorenz , x , 0.0 , 10.0 , dt , streaming_observer( std::cout ) );

How would you alter the observer to only log the state every 10 steps (for example). I'm wondering whether there is a more elegant solution than putting in an if-statement:

struct streaming_observer
{
  std::ostream &m_out;
  int count;

  streaming_observer( std::ostream &out ) : m_out( out ) {count = 10;}

  void operator()( const state_type &x , double t ) const
  {
    if( count == 10 )  {
      count = 1;
      m_out << t;
      for( size_t i=0 ; i < x.size() ; ++i )
        m_out << "\t" << x[i];
        m_out << "\n";
      }
      else {
        count++;
      }
   }
};
2

There are 2 answers

2
mariomulansky On BEST ANSWER

I had the same issue and solved it exactly like you did. However, you could also consider using a stepper with step-size control and then use integrate_const with a dt such that the observer is called at the required intervals. When you use a stepper with step-size control (even better: dense ouput like dopri5), integrate_const adjusts the step size according to your error tolerance, but then assures that the observer is called at times t0 + n*dt.

0
headmyshoulder On

Actually I would do it exactly like you did. You can also write a small adapter for doing the striding:

template< typename Obs >
struct striding_observer {
    size_t stride , count;
    Observer obs;
    striding_observer( size_t s , Obs o ) : stride(s) , count(1) , obs(o) { }
    template< typename S , typename T >
    void operator()( State const& s , Time const& t ) {
        if( count == stride ) {
            obs( s , t );
            count = 1;
        } else {
            ++count;
        }
    }
};

template< typename Obs >
striding_observer< Obs > make_striding_observer( size_t stride , Obs o ) {
    return striding_observer< Obs >( stride , o );
}

Then the striding is optional and composable. You can then write the very first example as

integrate_const( runge_kutta4< state_type >() ,
    lorenz , x , 0.0 , 10.0 , dt ,
    make_striding_observer( 10 , streaming_observer( std::cout ) ) );