#include <stdio.h>
#include <math.h>
#include "RK.h"

/*
	solve an ODE

		dx/dt + gain x = u(t)
		x(0)     = 0
		dx/dt(0) = 0
	
	the above can be converted as
	
	(d/dt) x = - gain x + u(t)
*/

#define	NumVars		1

#ifndef M_PI
#define	M_PI		3.14159265358979323846
#endif

static double pulse(double t)
{
	static double	period = 1.00;
	static double	duty = 0.50;
	static double	amplitude = 10.00;
	static double	offset = 0.00;
	
	if ( t - floor(t/period)*period <= period * duty ) {
		return amplitude + offset;
	} else {
		return offset;
	}
}


static void differential(int n, double x[], double t, double value[])
{
	static double	gain = 2.00;

	value[0] = - gain * x[0] + pulse(t);
}

main()
{
	int		kcnt;
	double	t, dt;
	double	x[NumVars];
	char	*filename = "oneorder.dat";
	FILE	*fd;
	
	if ((fd = fopen(filename,"w")) == NULL) {
		fprintf(stderr, "cannot open %s\n", filename);
		exit(1);
	}
	
	kcnt = 1000;
	dt = 0.01;
	/* kcnt * dt = 10 [sec] */

	/* initial value */
	t = 0.0;
	x[0] = 0.00;
	RungeKuttafprint(fd, NumVars, x, t);	// print initial values
	
	RungeKutta(fd, NumVars, x, differential, &t, dt, kcnt);
	
	fclose(fd);
}