#include<stdio.h>
#include<math.h>
#define M_PI 3.1415926535897932384626433
#include<gsl/gsl_odeiv2.h>
#include<gsl/gsl_errno.h>

int my_sin_prime (double t, const double y[], double dydt[], void * params) {
	dydt[0]=y[1];
	dydt[1]=-y[0];
	return GSL_SUCCESS;
}

double my_sin(double x) {
	if(x<0) return -my_sin(-x);
	int n = (int)x/2/M_PI;
	x -= n*2*M_PI;

	gsl_odeiv2_system s = { .function = my_sin_prime, .jacobian = NULL, .dimension = 2, .params = NULL };

	double hstart=1e-3, abs_prec_goal=1e-6, rel_prec_goal=1e-6;
	gsl_odeiv2_driver* d =
	gsl_odeiv2_driver_alloc_y_new( &s, gsl_odeiv2_step_rkf45, hstart, abs_prec_goal, rel_prec_goal);

	double t=0, y[]={0,1};
	int status = gsl_odeiv2_driver_apply( d, &t, x, y );
	if(status != GSL_SUCCESS) fprintf(stderr,"my_sin: status = %i\n",status);

	gsl_odeiv2_driver_free(d);
	return y[0];
}

int main(void) {
	double b=1.999999*M_PI;
	for(double x=-b; x<=b; x+=b/51) printf("%g %g\n", x, my_sin(x));
	printf("\n\n");
	for(double x=-b; x<=b; x+=b/7) printf("%g %g\n", x, sin(x));
return 0;
}
