Je ne peux pas trouver mon erreur dans ce programme Scheme pour calculer Pi

Question

Je suis en train de faire une expérience de Monte Carlo pour calculer une approximation de PI. De SICP:

  

La méthode de Monte Carlo est constitué de   choisir des expériences d'échantillons au hasard   d'un grand ensemble et puis faire   déductions sur la base du   probabilités estimée à partir de   tabuler les résultats de ces   expériences. Par exemple, nous pouvons   approximation en utilisant le fait que   6 / pi ^ 2 est la probabilité que deux   entiers choisis au hasard auront pas   facteurs en commun; c'est que leur   plus grand commun diviseur sera 1.   obtenir le rapprochement de nous   effectuer un grand nombre d'expériences.   Dans chaque expérience, nous choisissons deux   des nombres entiers au hasard et effectuer un test   pour voir si leur GCD est 1. La fraction   de fois que le test est passé donne   nous notre estimation de 6 / pi ^ 2, et de   ce que nous obtenons notre approximation   pi.

Mais quand je lance mon programme j'obtiens des valeurs telles que 3,9 ...

Voici mon programme:

(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)))

Mon interprète est MIT / GNU 07/07/90

Merci pour toute aide.

Answer 1

Eh bien, pour répondre à votre question directement, vous avez le instruction if arrière; il devrait en être ainsi.

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

Vous calculez: 1 - 6 / π ² dans la limite que le nombre de procès tend vers l'infini. On obtient "pi" = sqrt (6 / (1 - 6 / π ² )) = sqrt (6π ² / (π ² - 6)) = 3,911).

---

Prenons un peu de recul ici, cependant, et regarder ce que la méthode de Monte Carlo fait pour nous avec ce calcul (indice:.? Attendre la convergence très lent Combien de fois sont en cours d'exécution, il vous) ...

Chaque essai soit nous donne 0 ou 1, avec la probabilité p = 6 / π ² . Ceci est un exemple d'un processus de Bernoulli , qui a, pour le nombre de 1 s m dans un certain nombre d'essais n , .

Considérons ρ = m / n, la fraction de temps qui passent le test common-diviseurs. Ceci a une valeur moyenne de p et une variance de p (1-p) / n, ou std dev σ _ρ = sqrt (p (1-p) / n). Pour n = 10000, vous devriez vous attendre un std dev de 0,00488. moyennes et 5% du temps, vous serez en dehors de 2 std devs, ou entre 0,5982 et 0,6177. Ainsi, une estimation de π à partir de ce procédé, étant donné n = 10 000, sera comprise entre 3,117 et 3,167 95% du temps, et en dehors de cette plage de 5% du temps.

Si vous voulez augmenter le nombre d'essais de 100, qui réduit std dev par un facteur de 10 et l'estimation de π entre 3,1391 obtient rétrécies et 3,1441 avec 95% de confiance.

méthodes de Monte Carlo sont bien pour une estimation approximative, mais ils ont besoin de beaucoup et beaucoup d'essais pour une réponse précise, et atteignent habituellement un point de rendements décroissants.

Non que ce soit une façon inutile de se rapprocher pi, juste être conscient de la question.

Answer 2

Je trouve mon erreur. Merci à tous. J'incrémentant le cas avec succès au mauvais endroit.

Le code corrigé est:

(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)))