c - Improved Euler Method in Simple Harmonic Oscillator -

i have written c code using improved euler method determine position, velocity , energy of oscillator @ regular time intervals. however, run problem energy of oscillator decreasing, though there no dissipation terms. think particularly related way update position , velocity variables , on matter. code follows:

//compilation , run //gcc oscillatorimprovedeuler.c -lm -o oscillatorimprovedeuler && ./oscillatorimprovedeuler #include <stdio.h>  #include <math.h>    // global constans defined in following way (having constant value througout program #define m 1.0                 // kg #define k 1.0                 // kg/sec^2 #define h 0.1                // sec time step #define n 201 // number of time steps  int main(void) {     // avoid using arrays time     double x = 0, xint = 0;     double v = 5, vint = 0; // previous case     double t = 0;     double e = (m * v * v + k * x * x) / 2.0; // energy in units of joules      file *fp = fopen("oscillatorimprovedeuler.dat", "w+");     int = 0;     for(i = 0; < n ; i++)     {         fprintf(fp, "%f \t %f \t %f \t %f \n", x, v, e, t);         xint = x + (h) * v;         vint = v - (h) * k * x / m;         v = v - (h) * ((k * x / m) + (k * xint / m)) / 2.0;         x = x + (h) * (v + vint) / 2.0;         e = (m * v * v + k * x * x) / 2.0;         t += h;     }     fclose(fp);     return 0; } 

there may slight point miss grateful if can point out. appreciate help.

so figured out aid of math.stackexchange problem related updating position , velocity earlier time should updated , more intermediate variables needed. working code below:

//compilation , run //gcc oscillatorimprovedeuler.c -lm -o oscillatorimprovedeuler && ./oscillatorimprovedeuler #include <stdio.h>  #include <math.h>    // global constans defined in following way (having constant value througout program #define m 1.0                // kg #define k 1.0                // kg/sec^2 #define h 0.1                // sec time step #define n 200                // number of time steps  int main(void) {     // need define many variables avoid updating position , velocity     double x = 0.0, xpre = 0, xcor = 0;     double v = 5.0, vpre = 0, vcor = 0; // previous case     double t = 0;     double e = (m * v * v + k * x * x) / 2.0; // energy in units of joules      file *fp = fopen("oscillatorimprovedeuler.dat", "w+");     int = 0;     for(i = 0; < n ; i++)     {         if (i == 0)         {             fprintf(fp, "%f \t %f \t %f \t %f \n", x, v, e, t);         }         xpre = x + (h) * v;         vpre = v - (h) * k * x / m;         vcor = v - (h) * ((k * x / m) + (k * xpre / m)) / 2.0;         xcor = x + (h) * (v + vpre) / 2.0;         e = (m * vcor * vcor + k * xcor * xcor) / 2.0;         t += h;         fprintf(fp, "%f \t %f \t %f \t %f \n", xcor, vcor, e, t);         x = xcor, v = vcor;     }     fclose(fp);     return 0; }