#include "auto_f2c.h"
#include <math.h>


#define PI 3.14159265358979

#define EPSILON par[0]
#define A       par[1]
#define OMEGA   par[2]
#define Y_0      par[3]
#define PHI_0    par[4]
#define T       par[10]

#define X     u[0]
#define Y     u[1]
#define PHI   u[2]

/* ---------------------------------------------------------------------- */
/* ---------------------------------------------------------------------- */
int func (integer ndim, const doublereal *u, const integer *icp,
          const doublereal *par, integer ijac,
          doublereal *f, doublereal *dfdu, doublereal *dfdp)
{
  
  f[0] = T*((Y + X - ((X*X*X)/3.))/(EPSILON));
  f[1] = T*(-X + A*sin(PHI));
  f[2] = T*(2.*PI*OMEGA);
   return 0;
}
/* ---------------------------------------------------------------------- */
/* ---------------------------------------------------------------------- */
int stpnt (integer ndim, doublereal t,
           doublereal *u, doublereal *par)
{
  /* Parameter values for the starting orbit in fvdpmap.dat : */

  EPSILON = (doublereal) 0.0001;
  A       = (doublereal) 1.1;
  OMEGA   = (doublereal) 1.55;
  Y_0      = (doublereal) -0.7;
  PHI_0    = (doublereal) 2.26194671058465;

  return 0;
}
/* ---------------------------------------------------------------------- */
/* ---------------------------------------------------------------------- */
int bcnd (integer ndim, const doublereal *par, const integer *icp,
          integer nbc, const doublereal *u0, const doublereal *u1, integer ijac,
          doublereal *fb, doublereal *dbc)
{
  /* fb[0]= u0[0]+1; fb[1]=u1[0]-1; */
  fb[0] = u0[0]+1;
  fb[1] = u1[0]-1;
  fb[2] = u0[1]-Y_0;
  fb[3] = u0[2]-PHI_0;
  return 0;
}
/* ---------------------------------------------------------------------- */
/* ---------------------------------------------------------------------- */
int icnd (integer ndim, const doublereal *par, const integer *icp,
          integer nint, const doublereal *u, const doublereal *uold,
          const doublereal *udot, const doublereal *upold, integer ijac,
          doublereal *fi, doublereal *dint)
{
  return 0;
}
/* ---------------------------------------------------------------------- */
/* ---------------------------------------------------------------------- */
int fopt (integer ndim, const doublereal *u, const integer *icp,
          const doublereal *par, integer ijac,
          doublereal *fs, doublereal *dfdu, doublereal *dfdp)
{
  return 0;
}
/* ---------------------------------------------------------------------- */
/* ---------------------------------------------------------------------- */
int pvls (integer ndim, const doublereal *u,
          doublereal *par)
{
  return 0;
}
/* ---------------------------------------------------------------------- */
/* ---------------------------------------------------------------------- */