Simulation Toolkit (SmTk)

This example solves and make a phase plot of a damped simple 1-D pendulum. This maps the position back into the region (-n pi, m pi] or (-n,m] with the position defined by angle/pi. This will keep the plotter values in the visible region of the plot window and make plots nice numbers, given that the plotter doesn't know how to label axies nicely with pi (3.14159...).

#include <math.h>
#include <smtk.h>
#include <smtk/rk4.h>
#include <smtk/plotter_2d.h>

static const double pi = 3.1415926535897932385;
static const double num_pi = 4.0; // size of X region plus a little.
                                  // This should be a multiple of 2.
static const double xmin = -2.0,
                    xmax = xmin + num_pi; // X region (xmin, xmax]
static const double b = 0.1;           // parameter: damping constant
static const double k = sqrt(2);
static double       x[2]; // position (angular position), velocity (angular velocity)
static double       t;    // time

static void difeq(double *xdot, const double *x, double t)
{
  xdot[0] = k * x[1];
  xdot[1] = - sin(pi*x[0]) - b * x[1];
}

static smtk::rk4<double,double,double>
rk4(difeq, /* differential equations call-back function */
    2,     /* 2 equations */
    0.04   /* time step */);

static double *init(double x_in, double y_in, bool &do_draw)
{
  x[0] = x_in; // position
  x[1] = y_in; // velocity
  t = 0.0;  // time

  if(x[0] > xmax || x[0] <= xmin)
    {
      while(x[0] > xmax)
        x[0] -= 2*pi;

      while(x[0] <= xmin)
        x[0] += 2*pi;
    }

  return x;
}

static double *step(bool &do_draw, bool &show_line)
{
  rk4.go(x, &t);

  // Map the position back into the region (-pi, pi] or (-1,1] with
  // the position defined by angle/pi.  If it jumps more than 2 (2 pi)
  // then you've got bigger problems than this remapping.

  if(x[0] > xmax || x[0] <= xmin)
    {
      show_line = false;
      if(x[0] > xmax)
        x[0] -= num_pi;
      else if(x[0] <= xmin)
        x[0] += num_pi;
    }

  // usleep slows this down.  Remove/add this to speed/slow it up.
  //usleep(10000);

  return x;
}

int main(void)
{
  smtk::plotter_2d plotter;

  plotter.
    set_window_title("Damped Pendulum").
    set_callbacks(1, init, step).
    set_line_width(3).
    set_point_width(2).
    set_scale(xmin+0.2, xmax-0.2, -2.1, 2.1).
    run();

  return 0;
}

Reference	Tutorial	Example Code
smtk::plotter_2d	interactive phase space plotter	damped_SHO.cpp damped_pendulum.cpp lotka_volterra.cpp multi_plot.cpp

doc/examples/by_class/plotter_2d/damped_pendulum.cpp File Reference

Detailed Description