Natural Pi # 1 - Sable

9

Objectif

Générez ( N) des segments de ligne aléatoires de longueur uniforme ( l), vérifiez s'ils croisent les tlignes parallèles équidistantes ( ).

Simulation

Que simulons-nous? L'aiguille de Buffon . Lissez le sable dans votre bac à sable, tracez un ensemble de lignes parallèles également espacées (appelez la distance entre les deux t). Prenez un bâton droit de longueur let déposez-le Nfois dans le bac à sable. Soit le nombre de fois qu'il a franchi une ligne c. Alors Pi = (2 * l * n) / (t * c)!

Comment simulons-nous cela?

  • Prendre connaissance N,t,l
  • Avec N, t, ltous des entiers positifs
  • Procédez comme suit N:
    • Générer une coordonnée entière uniformément aléatoire x,y
    • Avec 1 <= x, y <= 10^6
    • x,y est le centre d'un segment de ligne de longueur l
    • Générer un entier uniformément aléatoire a
    • Avec 1 <= a <= 180
    • Soit Ple point où le segment de ligne traverserait l'axe des x
    • Alors ac'est l'angle(x,y), P, (inf,0)
  • Compter le nombre cde segments de ligne qui traversent la ligne x = i*tpour tout entieri
  • Revenir (2 * l * N) / (t * c)

entrez la description de l'image ici

entrez la description de l'image ici

spécification

  • Contribution
    • Flexible, saisissez les données de n'importe quelle manière standard (par exemple, paramètre de fonction, STDIN) et dans n'importe quel format standard (par exemple chaîne, binaire)
  • Production
    • Flexible, donne une sortie de n'importe quelle manière standard (par exemple retour, impression)
    • Les espaces blancs, les espaces blancs arrière et avant sont acceptables
    • Précision, veuillez fournir au moins 4 décimales de précision (c.-à-d. 3.1416)
  • Notation
    • Le code le plus court gagne!

Cas de test

Votre sortie peut ne pas correspondre à ceux-ci, en raison du hasard. Mais en moyenne, vous devriez obtenir cette précision pour la valeur donnée de N, t, l.

Input (N,t,l)    ->  Output 
-----------        ------
10,10,5          -> ?.????
10,100,50        -> ?.????
1000,1000,600    -> 3.????
10000,1000,700   -> 3.1???
100000,1000,700  -> 3.14??

TL; DR

Ces défis sont des simulations d'algorithmes qui ne nécessitent que la nature et votre cerveau (et peut-être quelques ressources réutilisables) pour approximer Pi. Si vous avez vraiment besoin de Pi pendant l'apocalypse zombie, ces méthodes ne gaspillent pas de munitions ! Il y a neuf défis au total.

NonlinearFruit
la source
Je pensais que tu avais déjà fait le numéro 1?
Conor O'Brien
1
@ ConorO'Brien I zero-index it XD
NonlinearFruit
le problème est que, dans les langues sans nombres complexes, vous devez transformer le nombre 0..180 en 0..pi, ce qui va plutôt à l'encontre du but de l'expérience de l'aiguille du buffon.
Level River St,
@NonlinearFruit la direction apeut-elle également être créée par une autre méthode, si elle est uniforme? (en pensant à une bulle de Gauss 2D)
Karl Napf
1
Peut-on supposer cela t > l? Deux solutions ci-dessous font cette hypothèse, ce qui simplifie un peu la vérification de l'intersection.
primo

Réponses:

9

R, 113 100 75 70 68 67 65 59 63 57 octets

En tant que langage de programmation statistique et fonctionnel, il n'est pas surprenant que R soit assez bien adapté à ce type de tâche. Le fait que la plupart des fonctions peuvent prendre une entrée vectorisée est vraiment utile pour ce problème, car plutôt que de boucler sur les Nitérations, nous passons simplement des vecteurs de taille N. Merci à @Billywob pour quelques suggestions qui ont conduit à couper 4 octets. Un grand merci à @Primo pour m'avoir patiemment expliqué comment mon code ne fonctionnait pas pour les cas où t > l, qui est maintenant corrigé.

pryr::f(2*l*N/t/sum(floor(runif(N)+sinpi(runif(N))*l/t)))

Essayez-le en ligne!

Exemple de sortie:

N=1000, t=1000, l=500
3.037975

N=10000, t=1000, l=700
3.11943

N=100000, t=1000, l=700
3.140351

Explication

Le problème se résume à déterminer si les deux xvaleurs de l'aiguille sont de chaque côté d'une ligne parallèle. Cela a des conséquences importantes:

  1. y-les valeurs ne sont pas pertinentes
  2. L'emplacement absolu sur l' xaxe n'est pas pertinent, seule la position par rapport aux lignes parallèles les plus proches.

Il s'agit essentiellement d'une tâche dans un espace à une dimension, où nous générons une ligne de longueur en [0, l] (l'angle adétermine cette longueur), puis nous vérifions combien de fois cette longueur dépasse t. L'algorithme approximatif est alors:

  1. Exemples de x1valeurs de [0, 1000000]. Puisque des lignes parallèles se produisent à chaque tpoint sur l' xaxe -x, la position relative xest xmodulo t.
  2. Échantillonnez un angle a.
  3. Calculez la x2position en fonction de a.
  4. Vérifiez combien de fois x1+x2s'inscrit t, c'est-à-dire prenez la parole (x1+x2)/t.

L'échantillonnage des Nnombres dans [0, 1e6] modulo test équivalent à un simple échantillonnage des Nnombres dans [0, t]. Puisque (x1+x2)/test équivalent à x1/t + x2/t, la première étape devient un échantillonnage à partir de [0, t] / t, c'est-à-dire [0, 1]. Heureusement pour nous, c'est la plage par défaut pour la runiffonction de R , qui renvoie Ndes nombres réels de 0 à 1 à partir d'une distribution uniforme.

                          runif(N)

Nous répétons cette étape pour générer a, l'angle de l'aiguille.

                                         runif(N)

Ces nombres sont interprétés comme un demi-tour (c'est-à-dire .590 degrés). (L'OP demande des degrés de 1 à 180, mais dans les commentaires, il est précisé que toute méthode est autorisée si elle est aussi précise ou plus précise.) Pour un angle θ, sin(θ)nous donne la distance sur l'axe des x entre les extrémités de l'aiguille. (Normalement, vous utiliseriez le cosinus pour quelque chose comme ça; mais dans notre cas, nous considérons l'angle θcomme étant relatif à l'axe des y, pas à l'axe des x (c'est-à-dire qu'une valeur de 0 degré augmente , non à droite ), et donc nous utilisons le sinus, qui fondamentalement décale les nombres.) Multiplié par lcela nous donne l' xemplacement de l'extrémité de l'aiguille.

                                   sinpi(runif(N))*l

Maintenant, nous divisons tet ajoutons la x1valeur. Cela donne (x1+x2)/t, qui est à quelle distance l'aiguille dépasse x1, en termes de nombre de lignes parallèles. Pour obtenir l'entier du nombre de lignes franchies, nous prenons le floor.

                    floor(runif(N)+sinpi(runif(N))*l/t)

Nous calculons la somme, en nous donnant le nombre cde lignes traversées par des aiguilles.

                sum(floor(runif(N)+sinpi(runif(N))*l/t))

Le reste du code implémente simplement la formule d'approximation de pi, c'est-à-dire (2*l*N)/(t*c). Nous économisons quelques octets entre parenthèses en profitant du fait que (2*l*N)/(t*c) == 2*l*N/t/c:

        2*l*N/t/sum(floor(runif(N)+sinpi(runif(N))*l/t))

Et le tout est enveloppé dans une fonction anonyme:

pryr::f(2*l*N/t/sum(floor(runif(N)+sinpi(runif(N))*l/t)))
rturnbull
la source
@rturnbull Nice one! Ne devriez-vous pas cependant pouvoir sauter les parenthèses au début? (2*l*N) => 2*l*N?
Billywob
@Billywob Bien repéré! Merci.
rturnbull
@rturnbull Oh et au fait, (2*l*N)/(t*c) = 2*l*N/t/cvous pouvez donc économiser encore deux octets en sautant également les parenthèses sur la dernière partie.
Billywob
@Billywob Encore une fois, bien repéré! Merci encore.
rturnbull
1
@primo Merci encore, cela devrait être corrigé maintenant.
rturnbull
6

Perl, 97 octets

#!perl -p
/ \d+/;$_*=2*$'/$&/map{($x=(1+~~rand 1e6)/$&)-$a..$x+($a=$'/$&/2*sin~~rand(180)*71/4068)-1}1..$_

En comptant le shebang comme un, l'entrée est tirée de stdin, l'espace séparé. Si des valeurs aléatoires non entières étaient autorisées, cela pourrait être un peu plus court.

J'ai pris une liberté, d'approximativement π / 180 comme 71/4068 , qui est précise dans 1,48 · 10 -9 .

Exemple d'utilisation

$ echo 1000000 1000 70000 | perl pi-sand.pl
3.14115345174061

Substitutions plus ou moins mathématiquement équivalentes

En supposant que la coordonnée x représente le point le plus à gauche de l'aiguille, plutôt que son milieu, comme spécifié dans la description du problème:

89 octets

#!perl -p
/ \d+/;$_*=2*$'/$&/map{($x=(1+~~rand 1e6)/$&)..$x+($'/$&*sin~~rand(180)*71/4068)-1}1..$_

Le problème spécifie que xdoit être échantillonné comme un entier aléatoire. Si nous projetons l'espacement des lignes sur un intervalle d'un, cela nous laissera des valeurs de la forme n/tavec 0 <= n < t, pas nécessairement uniformes, si elles tne se divisent pas également 1e6. En supposant qu'une distribution uniforme soit néanmoins acceptable:

76 octets

#!perl -p
/ \d+/;$_*=2*$'/$&/map{($x=rand)..$x+($'/$&*sin~~rand(180)*71/4068)-1}1..$_

Notez que puisque randsera toujours inférieur à un (et donc tronqué à zéro), il n'est pas nécessaire au début de la plage:

70 octets

#!perl -p
/ \d+/;$_*=2*$'/$&/map{1..(rand)+($'/$&*sin~~rand(180)*71/4068)}1..$_

En supposant que l'angle de l'aiguille ne doit pas nécessairement être un degré entier, mais seulement uniformément aléatoire:

59 octets

#!perl -p
/ \d+/;$_*=2*$'/$&/map{1..(rand)+abs$'/$&*sin rand$`}1..$_

En supposant que l'angle peut être n'importe quelle distribution uniforme:

52 octets

#!perl -p
/ \d+/;$_*=2*$'/$&/map{1..(rand)+abs$'/$&*sin}1..$_

Ce qui précède est une simulation mathématiquement correcte de l'aiguille de Buffon. Cependant, à ce stade, je pense que la plupart des gens conviendraient que ce n'est pas vraiment ce que la question demandait.


Vraiment le pousser

Nous pourrions simplement jeter la moitié des cas de test, chaque fois que le deuxième point de terminaison est à gauche du premier (plutôt que de les échanger):

47 octets

#!perl -p
/ \d+/;$_*=$'/$&/map{1..(rand)+$'/$&*sin}1..$_

Notez que les valeurs de tet lsont sans conséquence pour les résultats de l'expérience. Nous pourrions simplement les ignorer (en supposant implicitement qu'elles sont égales):

28 octets

#!perl -p
$_/=map{1..(rand)+sin}1..$_

Evidemment non compétitif, mais il faut bien l'admettre, il a une certaine élégance.

primo
la source
4

Python 2, 141 octets

port sans vergogne de rtumbull, sautant déjà ycar totalement inutile.

from math import*
from random import*
lambda N,t,l:(2.*l*N)/(t*sum(randint(1,1e6)%t+abs(cos(randint(1,180)*pi/180))*l>t for _ in range(N)))

Le problème est seulement que pi est déjà connu dans le programme.

Le voici (jouable au golf) avec un pi inconnu et aucune fonction trigonométrique

def g(N,t,l):
 c=0
 for _ in range(N):
    x,y=gauss(0,1),gauss(0,1);c+=randint(1,1e6)%t+abs(x/sqrt(x*x+y*y))*l>t
 return(2.*l*N)/(t*c)

x,yen gest seulement pour la direction.

Karl Napf
la source
Requiert from random import randint;from math import cos,pi. Échoue t < l, par exemple 1000000,1000,70000.
primo