Skip to content

modular odeint

Stephen Berry edited this page Jul 31, 2017 · 4 revisions

Running modular systems of ODEs with boost's odeint as the integration engine

Necessary changes to use odeint are commented with stars (***)

#include "Damper.h"
#include "Spring.h"

#include <boost/numeric/odeint.hpp>

using namespace asc;

int main()
{
   state_t x;
   x.reserve(100); // We reserve more space than necessary, but Ascent will only allocate what is needed
   double t = 0.0;
   double dt = 0.01;
   double t_end = 1.5;

   Body b0(x);
   Body b1(x);
   b1.m = 1.0;
   b1.s = 1.0;
   b1.v = 40.0;

   Spring spring(b0, b1);
   spring.k = 2000.0;

   Damper damper(b0, b1);
   damper.c = 5.0;

   boost::numeric::odeint::runge_kutta4<asc::state_t> integrator;
   Recorder recorder;

   auto system = [&](const asc::state_t& xs, asc::state_t& D, const double t)
   {
      // We must run the spring and damper before the body in order to accumulate forces
      x = xs; // *** We must copy the internal integrator's state (xs) to the global state (x) so that parameters such as body position and velocity are updated properly
      spring(x, D, t);
      damper(x, D, t);
      b1(x, D, t);
   };

   asc::state_t x0;
   while (t < t_end)
   {
      recorder({ t, b1.s });
      x0 = x; // *** We must copy to an initialization state because the integrator will change the passed in state
      integrator.do_step(system, x0, t, dt);
      t += dt; // *** odeint requires the time to be updated here
   }

   recorder.csv("spring-damper", { "t", "b1 position" });

   return 0;
}