No puedo encontrar mi error en este programa Esquema para calcular Pi

Question

Estoy haciendo un experimento de Monte Carlo para calcular una aproximación del PI. De SICP:

  

El método Monte Carlo consiste en   la elección de los experimentos de muestras al azar   de un conjunto grande y luego hacer   deducciones sobre la base de la   probabilidades estimadas a partir   la tabulación de los resultados de las   experimentos. Por ejemplo, podemos   aproximar utilizando el hecho de que   6 / pi ^ 2 es la probabilidad de que dos   números enteros escogidos al azar no tendrán   factores en común; es decir, que su   máximo común divisor será 1. Para   obtener la aproximación a, nos   realizar un gran número de experimentos.   En cada experimento elegimos dos   enteros al azar y realizar una prueba   para ver si su MCD es 1. La fracción   de veces que la prueba se pasa da   nosotros nuestra estimación de 6 / pi ^ 2, y desde   esto obtenemos nuestra aproximación a   pi.

Pero cuando ejecuto mi programa obtengo valores como 3.9 ...

Este es mi programa:

(define (calculate-pi trials)
  (define (this-time-have-common-factors?)
    (define (get-rand)
      (+ (random 9999999999999999999999999999999) 1))
    (= (gcd (get-rand) (get-rand)) 1))
  (define (execute-experiment n-times acc)
    (if (> n-times 0)
        (if (this-time-have-common-factors?)
            (execute-experiment (- n-times 1) acc)
            (execute-experiment (- n-times 1) (+ acc 1)))
        acc))
  (define n-success (execute-experiment trials 0))
  (define prob (/ n-success trials))
  (sqrt (/ 6 prob)))

Mi intérprete es el MIT / GNU 07/07/90

Gracias por cualquier ayuda.

Answer 1

Bueno, para responder a su pregunta directamente, usted tiene la sentencia if hacia atrás; debe ser de esta manera.

    (if (this-time-have-common-factors?)
        (execute-experiment (- n-times 1) (+ acc 1)
        (execute-experiment (- n-times 1) acc))

Así que estás cálculo de 1 - 6 / π ² en el límite cuando el # de ensayos tiende a infinito. Esto produce "pi" = sqrt (6 / (1 - 6 / π ² )) = sqrt (6π ² / (π ² - 6)) = 3.911).

---

Vamos a dar un paso atrás aquí, sin embargo, y ver lo que hace el método de Monte Carlo para nosotros con este cálculo (pista:.? Esperar una convergencia muy lenta ¿Cuántas veces usted está funcionando con ella) ...

Cada ensayo ya sea nos da 0 o 1, con una probabilidad de p = 6 / π ² . Este es un ejemplo de un proceso de Bernoulli , que tiene, para el número de de 1 m en un número de ensayos n , un binomial distribución .

Considere ρ = m / n, la fracción de veces que pasan la prueba de divisores comunes. Esto una tiene un valor de p y una varianza de p (1-p) / n, o std dev σ _ρ = sqrt (p (1-p) / n) significa. Para n = 10.000, que debe esperar un std dev de 0,00488. 95% del tiempo que será dentro de 2 std desarrolladores de los medias, y 5% del tiempo que va a estar fuera de desarrolladores 2 std, o entre 0,5982 y 0,6177. Así que una estimación de π de este método, dado n = 10 000, será de entre 3.117 y 3.167 95% del tiempo y fuera de este rango 5% del tiempo.

Si desea aumentar el número de ensayos por 100, que reduce std dev por un factor de 10 y la estimación de π consigue acortado entre 3,1391 y 3,1441, con una confianza del 95%.

métodos de Monte Carlo son agradables para una estimación aproximada pero necesitan montones y montones de ensayos para una respuesta precisa, y por lo general llegan a un punto de rendimientos decrecientes.

No es que esto es una forma infructuosa para aproximar pi, acaba de ser conscientes del problema.

Answer 2

Encuentro mi error. Gracias a todos. Yo estaba incrementando los casos exitosos en el lugar equivocado.

El código corregido es:

(define (calculate-pi trials)
  (define (this-time-have-common-factors?)
    (define (get-rand)
      (+ (random 9999999) 1))
    (= (gcd (get-rand) (get-rand)) 1))
  (define (execute-experiment n-times acc)
    (if (> n-times 0)
        (if (this-time-have-common-factors?)
            (execute-experiment (- n-times 1) (+ acc 1))
            (execute-experiment (- n-times 1) acc))
        acc))
  (define n-success (execute-experiment trials 0))
  (define prob (/ n-success trials))
  (sqrt (/ 6 prob)))