Integración de Monte Carlo

• El método de "Monte Carlo" es una técnica numérica que usa números pseudoaleatorios.
• La Integración de Monte Carlo estima el valor de una integral.
• Toma valores de la función en puntos aleatorios.
• El área (volumen) por el valor promedio estima la integral.

• ¿Cuándo se usa la Integración de Monte Carlo?
• Cuando la función a evaluar es complicada.
• Normalmente para integrales múltiples.
• La integración de Monte Carlo siempre funciona, pero toma mucho tiempo de CPU para obtener una respuesta más exacta.
• Pero, CPU son baratas.

Es decir, los métodos de integración de Montecarlo son algoritmos para encontrar una evaluación aproximada de una integral definida, normalmente de integrales múltiples. Los algoritmos deterministas de integración numérica, para aproximar la integral, evalúan la función en un conjunto de puntos correspondientes a una parrilla regular o en un conjunto de puntos predefinidos. 

En cambio, los métodos de Montecarlo eligen de forma aleatoria los puntos en los que se evaluará la función. La integración de Montecarlo forma parte de una familia de algoritmos llamados genéricamente métodos de Montecarlo.

Método:

• Recordemos el primer año de cálculo: El teorema del valor medio.

 

• La integral de la función f, sobre el volumen V, es igual a la media de f por V.
•  La integración de Monte Carlo esta dado por:
• Estimación del valor medio de la función.
• Multiplicado por el volumen. 

•  Estimación del valor medio de la función: 
• La media de una función es calculado por la media simple (también llamado promedio)

•  La varianza de la media simple es

• El estimado de la media de la función donde generamos un gran número de valores aleatorios de xi, y los usamos para calcular la media simple de la función.

 

•  Calculando la Integral de Monte Carlo

 

• EJEMPLO: Cálculo del número PI

/* Autor: Jason Martinez
 * Calculating number PI, with Monte Carlo Integation
 *            x^2 + y^2 <= 1
 * A = PI*r^2
*/

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

const double pi = 3.14159265358979323846264338327950288419716939937510;

int main()
{
  double r, x, y;
  double avg;
  double sum, A;
  double x_max, x_min, y_max, y_min;
  int n;  //number of iterations
  int i;
  double pi_calc, err;

  scanf("%d", &n);

  sum = 0;
  A = 0;
  r = 1.0;
  x_max = y_max = 1.0;
  x_min = y_min = -1.0;
  
  for(i=0; i<n; ++i)
  {
    x = ((double)rand()/RAND_MAX)*(x_max - x_min) + x_min;
    y = ((double)rand()/RAND_MAX)*(y_max - y_min) + y_min;
    if(x*x + y*y <= r)
      sum += 1.0;
  }

  avg = sum / n; 

  A = (x_max - x_min)*(y_max - y_min) * avg;
  
  pi_calc = A/(r*r);
  printf("pi = %lf\n", pi_calc);
  err = fabs(pi - pi_calc);
  printf("error = %lf\n", err);
  return 1;
}