Monte Carlo Integral Calculator

mcint2

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

double f(double x)
{
  return exp(x) * x;
}

double random(double fs, double ls)
{
  return ((rand() / (double) RAND_MAX) *  (ls - fs)) + fs;
}

int main()
{
  srand(time(NULL));
  int N, it, b, a, i, j;
  // For a = 2 and b = 5 and f(x) = x * e^x
  double realvalue = 586.264;
  printf("Input number of iterations: ");
  scanf("%d", &it);
  printf("Input number of samples: ");
  scanf("%d", &N);
  printf("Input the boundaries (a, b): ");
  scanf("%d %d", &a, &b);

  double sum, mx = f(a), mn = 0, dif = ((double)b - a) / N;
  int inter = 0, tot = 0;

  // Maximum Value
  for (j = 0; j < N; j++)
  {
    if (f(a + j * dif) > mx)
      mx = f(a + j * dif);
    if (f(a + j * dif) < mn)
      mn = f(a + j * dif);
  }

  // Monte Carlo Method Calculation
  for (i = 0; i < it; i++)
    for (j = 0; j < N; j++)
    {
      double va = random(a, b), vb = random(mn, mx);
      inter += (f(va) > vb);
      tot++;
    }
  sum = mx * (b - a) * ((double)inter) / tot;

  printf("Real Value: %lf\n", realvalue);
  printf("Integral: %lf, error: %0.3lf%%\n", sum, (100 * fabs(sum - realvalue) / realvalue));
}