I have written a C program to trace out the motion of a double pendulum, but am having difficulties in getting gnuplot (controlled from my c program) to trace out the paths of the masses (example). Thus far I have created the program such that it produces a number of png images at each interval (using runge kutta method), however I want to output it as a gif instead so a line traces out the path of the masses in real time.
In the code below I am assuming the problem occurs with piping out to gnuplot from within the for loop (to save you from wasting your time sifting through it)
/* Header Files */
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>
#include <assert.h>
// Definitions
#define GRAVITY 9.8
#define INCREMENT 0.0175
// Declerations of functions
double th1_xder(double t1d);
double th1_yder(double t2d);
double th2_xder(double t1d, double th1, double t2d, double th2, double l1, double l2, double m1, double m2);
double th2_yder(double t1d, double th1, double t2d, double th2, double l1, double l2, double m1, double m2);
int main (int argc, char * argv[]) {
double x1=0, y1=0, x2=0, y2=0; //coordinates
double l1=0, l2=0; //lengths
double m1=0, m2=0; //masses
double th1=0, th2=0; //angles
double t1d=0, t2d=0; //zero velocity initially
//second order Runge-Kutta equations
double k1=0, k2=0, k3=0, k4=0; //for x1
double q1=0, q2=0, q3=0, q4=0; //for y1
double r1=0, r2=0, r3=0, r4=0; //for x2
double s1=0, s2=0, s3=0, s4=0; //for y2
int i = 0;
double x0=0, y0=0;
printf("Enter '1' to input your own data or '0' to use the preset data\n");
char dummy = 'a';
scanf("%c", &dummy);
assert((dummy == '1') || (dummy == '0'));
if (dummy == '1') {
printf("Please enter a length l1:\n");
scanf("%lf", &l1);
printf("Please enter a length l2:\n");
scanf("%lf", &l2);
printf("Please enter a mass m1:\n");
scanf("%lf", &m1);
printf("Please enter a mass m2:\n");
scanf("%lf", &m2);
printf("Please enter an angle theta1:\n");
scanf("%lf", &th1);
printf("Please enter an angle theta2:\n");
scanf("%lf", &th2);
} else {
l1 = 1;
l2 = 1;
m1 = 1;
m2 = 1;
th1= 90;
th2= 0;
}
th1 = th1*(M_PI/180);
th2 = th2*(M_PI/180);
FILE *fsp;
if((fsp=fopen("origin.dat", "w"))==NULL) {
fprintf(stdout, "cannot open origin.dat\n");
exit (EXIT_FAILURE);
}
fprintf(fsp, "0\t0");
fclose(fsp);
FILE *fout;
if((fout=fopen("testout.dat", "w"))==NULL) {
fprintf(stdout, "cannot open testout.dat\n");
exit (EXIT_FAILURE);
}
printf("%f\t%f\t%f\t%f\n", x1, y1, x2, y2);
fprintf(fout, "%f\t%f\t%f\t%f\n", x1, y1, x2, y2);
for(i = 0; i < 250; i++) {
if ((fout=fopen("testout.dat", "w"))==NULL) {
fprintf(stdout, "cannot open testout.dat\n");
exit (EXIT_FAILURE);
}
k1 = th1_xder(t1d);
q1 = th1_yder(t2d);
r1 = th2_xder(t1d, th1, t2d, th2, l1, l2, m1, m2);
s1 = th2_yder(t1d, th1, t2d, th2, l1, l2, m1, m2);
k2 = th1_xder(t1d + (r1/2));
q2 = th1_yder(t2d + (s1/2));
r2 = th2_xder((t1d + (r1/2)), (th1 + (k1/2)), (t2d +(s1/2)), (th2 + (q1/2)), l1, l2, m1, m2);
s2 = th2_yder((t1d + (r1/2)), (th1 + (k1/2)), (t2d +(s1/2)), (th2 + (q1/2)), l1, l2, m1, m2);
k3 = th1_xder(t1d + (r2/2));
q3 = th1_yder(t2d + (s2/2));
r3 = th2_xder((t1d + (r2/2)), (th1 + (k2/2)), (t2d +(s2/2)), (th2 + (q2/2)), l1, l2, m1, m2);
s3 = th2_yder((t1d + (r2/2)), (th1 + (k2/2)), (t2d +(s2/2)), (th2 + (q2/2)), l1, l2, m1, m2);
k4 = th1_xder(t1d + r3);
q4 = th1_yder(t2d + s3);
r4 = th2_xder((t1d + r3), (th1 + k3), (t2d + s3), (th2 + q3), l1, l2, m1, m2);
s4 = th2_yder((t1d + r3), (th1 + k3), (t2d + s3), (th2 + q3), l1, l2, m1, m2);
t1d = t1d + (r1 + 2*r2 + 2*r3 + r4)/6;
t2d = t2d + (s1 + 2*s2 + 2*s3 + s4)/6;
th1 = th1 + (k1 + 2*k2 + 2*k3 + k4)/6;
th2 = th2 + (q1 + 2*q2 + 2*q3 + q4)/6;
x1 = l1*sin(th1);
y1 = -l1*cos(th1);
x2 = x1 + l2*sin(th2);
y2 = y1 - l2*cos(th2);
printf("%f\t%f\t%f\t%f\n", x1, y1, x2, y2);
fprintf(fout, "%f\t%f\t%f\t%f\n", x1, y1, x2, y2);
fclose(fout);
FILE *gnuplotPipe = popen("gnuplot -persist","w");
if (gnuplotPipe) {
fprintf(gnuplotPipe, "set style data lines\n");
fprintf(gnuplotPipe, "set terminal png nocrop enhanced size 1280,720; set output 'yyy%d.png'\n", i);
fprintf(gnuplotPipe, "set title 'frame%d'\n", i);
fprintf(gnuplotPipe, "set multiplot\n");
fprintf(gnuplotPipe, "set xrange [-2.5:2.5]; set yrange [-2.5:2]\n");
fprintf(gnuplotPipe, "unset key; unset ytics; unset xtics\n");
fprintf(gnuplotPipe, "plot 'testout.dat' using 3:4\n");
fprintf(gnuplotPipe, "plot '-' with lines lw 2 lc rgb 'black', 'testout.dat' u 1:2 w points pt 7 ps 2, 'testout.dat' u 3:4 w points pt 7 ps 2, 'origin.dat' u 1:2 w points pt 7 ps 2 lc 0\n");
fprintf(gnuplotPipe, "%f %f\n", x0, y0);
fprintf(gnuplotPipe, "%f %f\n", x1, y1);
fprintf(gnuplotPipe, "%f %f\n", x2, y2);
fprintf(gnuplotPipe, "e\n");
fprintf(gnuplotPipe, "\n");
fprintf(gnuplotPipe, "set nomultiplot\n");
fflush(gnuplotPipe);
fprintf(gnuplotPipe,"exit \n");
pclose(gnuplotPipe);
}
}
return EXIT_SUCCESS;
}
double th1_xder(double t1d) {
double k = t1d*INCREMENT;
return k;
}
double th1_yder(double t2d) {
double m = t2d*INCREMENT;
return m;
}
double th2_xder(double t1d, double th1, double t2d, double th2, double l1, double l2, double m1, double m2) {
double l = INCREMENT*((GRAVITY/l1)*((m2/(m1 + m2))*sin(th2)*cos(th1-th2)-sin(th1))-(m2/(m1 + m2))*sin(th1-th2)*((l2/l1)*t2d*t2d + t1d*t1d*cos(th1-th2)))/(1-((m2/(m1 + m2))*cos(th1-th2)*cos(th1-th2)));
return l;
}
double th2_yder(double t1d, double th1, double t2d, double th2, double l1, double l2, double m1, double m2) {
double p = INCREMENT*((GRAVITY/l2)*(sin(th1)*cos(th1-th2)-sin(th1)) + sin(th1-th2)*((l1/l2)*t1d*t1d + (m2/(m1 + m2))*t2d*t2d*cos(th1-th2)))/(1-((m2/(m1 + m2))*cos(th1-th2)*cos(th1-th2)));
return p;
}