Runge Kutta in C amateur

No momento tento implementar o método de Runge-Kutta de 4º Ordem para equações diferenciais, é um método meio rustico ainda, não está bem lapidado com uma escrita em C, bem organizada e legível. Uma justificativa é o fato de ser iniciante nesta linguagem.
Uma descrição do método de Runge- Kutta de 4º ordem, consiste em determinar constantes apropriadas tais que a fórmula:

yn+1=yn+ak1+bk2+ck3dk4

Que tenha uma certa concordância com um desenvolvimento de Taylor até h⁴ , isto é até o quinto termo. O resultado do método é dado pelas seguintes equações.

yn+1=yn+(k1+k2+k3+k4)/6

k1=h f(xn,yn)

k2=h f(xn+h/2, yn+k1/2)

k3=h f(xn+h/2, yn+k2/h)

k4=h f(xn+h, yn+k3)

Onde h é o passo dado no intervalo de integração, ou seja, se meu intervalo é [0;1] então com um h=0,2, o intervalo terá a seguinte sequência numérica {0, 0,2, 0,4, 0,6, 0,8, 1,0}.
O código fonte em C que trago aqui é para uma equação em específico. Ainda estou trabalhando em um código que eu somente coloque a função sem ter que modificar o código fonte.

Segue o código em C do Método de Runge Kutta, para a resolução da equação:

dy/dx = x-y (1)
ou
y'=x-y

/*********************************************************************/
/*Bibliotecas*/
#include
#include
#include
FILE*rk4;
/*Inicio do programa*/
main(){
rk4=fopen("rk4.txt","w");
int n;
double x,yv,h,yn,x0;
double f1,f2,f3,f4;

h=0.2;
x0=0;
yv=2;

for(n=0;n<=10;n++){
/*Variaveis do método dy/dx=x+y, aqui já na equação diferencial a ser integrada*/

f1=h*(x+yv);
f2=h*((x+h/2)+(yv+f1/2));
f3=h*((x+h/2)+(yv+f2/2));
f4=h*((x+h)+(yv+f3));

/* Escreve no terminal e no arquivo*/

fprintf(rk4,"%lf%18lf\n",x,yn);
printf("%lf%18lf\n",x,yn); }

/**** Escreve no terminal e no arquivo*****/

fprintf(rk4,"%lf%18lf\n",x,yn); printf("%lf%18lf\n",x,yn);

}
fclose(rk4); }

/*Fim do programa*/
/************************************************************************/

A solução algébrica da equação (1)

y(x)=3exp(x)-x-1 (2)

A seguir os gráficos da equação resolvida pelo método e algebricamente, apesar dos erros obtidos pelos calculos computacionais, podemos notar que o seu bom resultado gráficamente.

Gráfico calculado pelo método de Runge Kutta de quarta ordem














Gráfico obtido pela função exata calculada algebricamente.