/*  prog6a.c -- Solve First Order ODE by RK2 Method
    Written Apr 19, 2001 by Eric Olson for Mathematics 483

    For information see Section 4.4 in Numerical Analysis for Applied 
    Mathematics, Science, and Engineering, by Greenspan and Casulli */

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

inline double F(double x,double y){
    return y*sin(x);
}

inline double Y(double x){
    return exp(1.0-cos(x));
}

double RK2(double x0,double y,double xn,int n){
    double h=(xn-x0)/n;
    int i;
    for(i=0;i<n;i++){
        double x=x0+i*h;
        double K0=h*F(x,y);
        double K1=h*F(x+h,y+K0);
        y+=0.5*(K0+K1);
    }
    return y;
}

main(){
    double Yn=Y(10.0);
    int n;
    printf("#%7s %16s %20s\n","n","yn","err");
    for(n=16;n<=1048576;n*=2){
        double yn=RK2(0.0,1.0,10.0,n);
        double err=fabs(yn-Yn);
        printf("%8ld %16.10lg %20.10le\n",n,yn,err);
    }
}