/*
 *  pendulum_def.c
 *
 *  DSTOOL C file for the vector field named: pendulum
 *
 *  This file was generated by the program VFGEN (Version:2.4.0)
 *  Generated on 10-Jul-2008 at 13:32
 */

#include <model_headers.h>
#include <math.h>

/*
 *  The vector field.
 */

int pendulum(double *f_, double *Y_, double *p_)
    {
    const double Pi = M_PI;
    double theta, v;
    double g, b, L, m;

    theta      = Y_[0];
    v          = Y_[1];

    g          = p_[0];
    b          = p_[1];
    L          = p_[2];
    m          = p_[3];

    double t = Y_[2];


    f_[0] = v;
    f_[1] = -1.0/m*v/(L*L)*b-sin(theta)*g/L;
    }

/*
 *  The Jacobian.
 */

int pendulum_jac(double **jac_, double *Y_, double *p_)
    {
    const double Pi = M_PI;
    double theta, v;
    double g, b, L, m;

    theta      = Y_[0];
    v          = Y_[1];

    g          = p_[0];
    b          = p_[1];
    L          = p_[2];
    m          = p_[3];

    double t = Y_[2];

    jac_[0][0] = 0.0;
    jac_[0][1] = 1.0;
    jac_[1][0] = -g/L*cos(theta);
    jac_[1][1] = -1.0/m/(L*L)*b;
    }

/*
 *  The derivative with respect to the independent variable.
 */

int pendulum_dfdt(double *dfdt_, double *Y_, double *p_)
    {
    const double Pi = M_PI;
    double theta, v;
    double g, b, L, m;

    theta      = Y_[0];
    v          = Y_[1];

    g          = p_[0];
    b          = p_[1];
    L          = p_[2];
    m          = p_[3];

    double t = Y_[2];

    dfdt_[0] = 0.0;
    dfdt_[1] = 0.0;
    }

/*
 *  The Jacobian with respect to the parameters.
 */

int pendulum_dfdp(double **dfdp_, double *Y_, double *p_)
    {
    const double Pi = M_PI;
    double theta, v;
    double g, b, L, m;

    theta      = Y_[0];
    v          = Y_[1];

    g          = p_[0];
    b          = p_[1];
    L          = p_[2];
    m          = p_[3];

    dfdp_[0][0] = 0.0;
    dfdp_[0][1] = 0.0;
    dfdp_[0][2] = 0.0;
    dfdp_[0][3] = 0.0;
    dfdp_[1][0] = -sin(theta)/L;
    dfdp_[1][1] = -1.0/m*v/(L*L);
    dfdp_[1][2] =  sin(theta)*g/(L*L)+2.0*1.0/m*v/(L*L*L)*b;
    dfdp_[1][3] = 1.0/(m*m)*v/(L*L)*b;
    }

/*
 *  User-defined aux functions
 */

int pendulum_aux(double *f_, double *Y_, double *p_)
    {
    const double Pi = M_PI;
    double theta, v;
    double g, b, L, m;

    theta      = Y_[0];
    v          = Y_[1];

    g          = p_[0];
    b          = p_[1];
    L          = p_[2];
    m          = p_[3];

    double t = Y_[2];


    /*
     *  Aux function: energy
     */
    f_[0] = -m*g*L*cos(theta)+m*(v*v)*(L*L)/2.0;
    }

/*
 *  The initialization function
 */

int pendulum_init()
    {
    /* ------------ Define the dynamical system in this segment ---------- */

    const double Pi = M_PI;
    int             n_varb = 2;
    static char     *variable_names[] = {"theta","v"};
    static double   variables[] = {0.0,0.0};
    static double   variable_min[] = {-10.0,-10.0};
    static double   variable_max[] = {10.0,10.0};
    static char     *indep_varb_name = "t";
    static double   indep_varb_min = 0.0;
    static double   indep_varb_max = 10000.0;

    int             n_param = 4;
    static char     *parameter_names[] = {"g","b","L","m"};
    static double   parameters[] = {9.8100000000000005e+00,0.0000000000000000e+00,1.0000000000000000e+00,1.0000000000000000e+00};
    static double   parameter_min[] = {-10.0,-10.0,-10.0,-10.0};
    static double   parameter_max[] = {10.0,10.0,10.0,10.0};

    int             n_funct = 1;
    static char     *funct_names[] = {"energy"};
    static double   *funct_min[] = {-100.0};
    static double   *funct_max[] = {100.0};

    int             manifold_type = PERIODIC;
    static int      periodic_varb[] = {TRUE,FALSE};
    static double   period_start[] = {0,0.0};
    static double   period_end[] = {2*Pi,0.0};

    int             mapping_toggle = FALSE;
    int             inverse_flag   = FALSE;

    int             (*def_name)() = pendulum;
    int             (*jac_name)() = pendulum_jac;
    int             (*aux_func_name)() = pendulum_aux;
    int             (*inv_name)() = NULL;
    int             (*dfdt_name)() = pendulum_dfdt;
    int             (*dfdparam_name)() = pendulum_dfdp;

    c_filename = __FILE__;

    /* ------------ End of dynamical system definition          ---------- */
    #include <ds_define.c>
    return 0;
    }