/**********************************************************************
* File: vibration3.c
* Animate the damped free vibration of
* overdamped, critical damped, underdamped systems.
* The output is for animation coordinate data.
* (1) Run this program in Ch.
* (2) Click "Go" to view the animation
* See CHHOME/docs/qanimate.pdf for detailed description of this program.
* Note: The details about this damped free vibration can be found in
* an exercise in Chapter 6 Functions in the book 
* "C for Engineers and Scientists: An Interpretive Approach"
* by Harry H. Cheng, published by McGraw-Hill, 2009,
* ISBN: 0073376051, ISBN-13: 978-0073376059.
**********************************************************************/
#include <stdio.h>
#include <math.h>

/* The amplitude of the vibration is 4.
   The initial velocity of the vibration is 0 */
double overdamped(double t) {
    return 4.12*exp(-1.57*t) - 0.12*exp(-54.2*t);
}

double criticaldamped(double t) {
    return 4*(1+6*t)*exp(-6*t);
}

double underdamped(double t) {
    return 4.06*exp(-0.5*t)*sin(3*t+1.4);
}

int main() {
    double t, t0, tf;          // time
    double y1, y2, y3;         // displacement
    double pin1x = 0, pin1y=7, // pin 1
           pin2x = 4, pin2y=7, // pin 2
           pin3x = 8, pin3y=7; // pin 3
    FILE *stream;

    stream = popen("qanimate","w");     // open qanimate pipe
    if (stream==NULL) {
        fprintf(stderr, "Error: popen() failed\n");
        exit(1);
    }

    /* A comment line starting with # */
    fprintf(stream, "# qanimate data for animation of vibration systems\n");
    /* The title displayed on the animation */
    fprintf(stream, "title \"Damped Free Vibration\"\n");
    fprintf(stream, "fixture\n");
    /* The primitives following fixture */
    fprintf(stream, "groundpin %f %f angle 180 joint %f %f\n", pin1x, pin1y, pin1x, pin1y);
    fprintf(stream, "line %f %f %f %f\n", pin1x, pin1y, pin1x, pin1y-1 );
    fprintf(stream, "text %f %f \"overdamped\"\n", pin1x+0.5, pin1y);
    fprintf(stream, "groundpin %f %f angle 180 joint %f %f\n", pin2x, pin2y, pin2x, pin2y);
    fprintf(stream, "line %f %f %f %f\n", pin2x, pin2y, pin2x, pin2y-1 );
    fprintf(stream, "text %f %f \"critical damped\"\n", pin2x+0.5, pin2y);
    fprintf(stream, "groundpin %f %f angle 180 joint %f %f\n", pin3x, pin3y, pin3x, pin3y);
    fprintf(stream, "line %f %f %f %f\n", pin3x, pin3y, pin3x, pin3y-1 );
    fprintf(stream, "text %f %f \"underdamped\"\n", pin3x+0.5, pin3y);
    fprintf(stream, "dot 11 7 pen white\n"); // to display all text corretly
    fprintf(stream, "animate restart\n");

    t0 = 0;
    tf = 10;
    for(t = t0; t<tf; t += 0.01) {
        y1 = overdamped(t);
        y2 = criticaldamped(t);
        y3 = underdamped(t);
        fprintf(stream, "rectangle %f %f %f %f fill red \\\n",-0.5, y1-1.0, 1.0, 1.0);
        fprintf(stream, "spring %f %f %f %f \\\n", pin1x, pin1y-1, pin1x, y1);
        fprintf(stream, "rectangle %f %f %f %f fill green \\\n", pin2x-0.5, y2-1.0, 1.0, 1.0);
        fprintf(stream, "spring %f %f %f %f \\\n", pin2x, pin2y-1, pin2x, y2);
        fprintf(stream, "rectangle %f %f %f %f fill blue \\\n", pin3x-0.5, y3-1.0, 1.0, 1.0);
        fprintf(stream, "spring %f %f %f %f \n", pin3x, pin3y-1, pin3x, y3);
    }
    pclose(stream);
    return 0;
}