#ifndef _LDBSODE_H
#define _LDBSODE_H

#include "f_eval.h"


typedef struct _ldbsode
{
    int N,ISTATE;
    long double EPS, HMIN, H0;
    long double *Y,T;
    LD_EVAL F;
} LDBSODE;


/*
** BSODE is a Bulirsch-Stoer method, from Press et al.
** n = number of equations.
** time = initial time.
** y = pointer to your array of n state variables.
** f = function of yours that computes ydot
** stiff_flag: 0 = no is the only option.
** eps = allowable relative error.
** zero = "problem zero"; I think it's allowable absolute error.
*/ 
LDBSODE *LDBsodeAlloc(int n, long double t, long double *y, LD_EVAL f, int stiff_flag,
    long double eps, long double zero);

/* Do the actual integration.  Return the actual TOUT */
long double LDBsodeIntegrate(LDBSODE *b, long double TOUT);

void LDBsodeFree(LDBSODE *b);

#endif  /* _BSODE_H */
