Golf aléatoire du jour # 4: Le paradoxe de Bertrand

19

À propos de la série

Tout d'abord, vous pouvez traiter cela comme n'importe quel autre défi de golf de code et y répondre sans vous soucier de la série. Cependant, il existe un classement pour tous les défis. Vous pouvez trouver le classement avec plus d'informations sur la série dans le premier post .

Bien que j'ai un tas d'idées alignées pour la série, les défis futurs ne sont pas encore gravés dans le marbre. Si vous avez des suggestions, faites-le moi savoir sur le post sandbox correspondant .

Trou 4: Le paradoxe de Bertrand

Le paradoxe de Bertrand est un problème intéressant, qui montre comment différentes méthodes de sélection d'accords aléatoires dans un cercle peuvent produire différentes distributions d'accords, leurs points médians et leurs longueurs.

Dans ce défi, vous êtes censé générer des accords aléatoires du cercle unitaire, en utilisant la méthode "droite", c'est-à-dire qui produit une distribution des accords qui est invariante sous l'échelle et la traduction. Dans l'article Wikipédia lié, la "Méthode 2" est une telle méthode.

Voici les règles exactes:

  • Vous devez prendre un entier positifN qui spécifie le nombre d'accords à renvoyer. La sortie doit être une liste d' Naccords, chacun spécifié en deux points sur le cercle unitaire, donné par leur angle polaire en radians.
  • Votre code doit pouvoir renvoyer au moins 2 20 valeurs différentes pour chacun des deux angles . Si votre RNG disponible a une gamme plus petite, vous devez d'abord construire un RNG avec une gamme suffisamment large au-dessus de celui intégré ou vous devez implémenter votre propre RNG approprié . Cette page peut être utile pour cela.
  • La distribution des accords doit être indiscernable de celle produite par la "Méthode 2" dans l'article Wikipedia lié. Si vous implémentez un algorithme différent pour choisir les accords, veuillez inclure une preuve d'exactitude. Quel que soit l'algorithme que vous choisissez d'implémenter, il doit théoriquement être capable de générer n'importe quel accord valide dans le cercle unitaire (sauf limitations des PRNG sous-jacents ou types de données à précision limitée).
  • Votre implémentation doit utiliser et renvoyer des nombres à virgule flottante (au moins 32 bits de large) ou des nombres à virgule fixe (au moins 24 bits de large) et toutes les opérations arithmétiques doivent être précises dans un maximum de 16 ulp .

Vous pouvez écrire un programme complet ou une fonction et prendre une entrée via STDIN (ou l'alternative la plus proche), un argument de ligne de commande ou un argument de fonction et produire une sortie via STDOUT (ou l'alternative la plus proche), une valeur de retour de fonction ou un paramètre de fonction (out).

La sortie peut être dans n'importe quel format de liste ou de chaîne, à condition que les numéros individuels soient clairement distinguables et que leur nombre total soit toujours pair.

Il s'agit du code golf, donc la soumission la plus courte (en octets) l'emporte. Et bien sûr, la soumission la plus courte par utilisateur entrera également dans le classement général de la série.

Visualisation

Vous pouvez utiliser l'extrait de code suivant pour rendre les lignes générées et inspecter leur distribution. Collez simplement une liste de paires d'angles dans la zone de texte. L'extrait doit être capable de gérer presque tous les formats de liste, tant que les nombres sont de simples nombres décimaux (pas de notation scientifique). Je vous recommande d'utiliser au moins 1000 lignes pour avoir une bonne idée de la distribution. J'ai également fourni quelques exemples de données pour les différentes méthodes présentées dans l'article ci-dessous.

Exemple de données générées avec la méthode 1.

Exemple de données générées avec la méthode 2.

Exemple de données générées avec la méthode 3.

Classement

Le premier post de la série génère un classement.

Pour vous assurer que vos réponses s'affichent, veuillez commencer chaque réponse par un titre, en utilisant le modèle Markdown suivant:

# Language Name, N bytes

Nest la taille de votre soumission. Si vous améliorez votre score, vous pouvez conserver les anciens scores dans le titre, en les rayant. Par exemple:

# Ruby, <s>104</s> <s>101</s> 96 bytes

(La langue n'est pas actuellement affichée, mais l'extrait de code l'exige et l'analyse, et je pourrai ajouter un classement par langue à l'avenir.)

Martin Ender
la source

Réponses:

12

Code machine IA-32, 54 octets

Hexdump du code:

68 00 00 00 4f 0f c7 f0 50 db 04 24 58 d8 34 24
f7 d9 78 f1 d9 c0 dc c8 d9 e8 de e1 d9 fa d9 c9
d9 f3 dc c0 d9 eb de ca d8 c1 dd 1a dd 5a 08 83
c2 10 e2 d1 58 c3

Il utilise un algorithme (légèrement modifié) décrit par Wikipedia. En pseudo-code:

x = rand_uniform(-1, 1)
y = rand_uniform(-1, 1)
output2 = pi * y
output1 = output2 + 2 * acos(x)

J'utilise la plage -1...1car il est facile de faire des nombres aléatoires dans cette plage: l' rdrandinstruction génère un entier entre -2^31et 2^31-1, qui peut être facilement divisé par 2 ^ 31.

J'aurais dû utiliser la plage 0...1pour l'autre nombre aléatoire (x), qui est introduit dans acos; cependant, la partie négative est symétrique avec la partie positive - les nombres négatifs produisent des accords dont la portée est supérieure à pi radians, mais dans le but d'illustrer le paradoxe de Bertrand, cela n'a pas d'importance.

Étant donné que le jeu d'instructions 80386 (ou x87) n'a pas d' acosinstruction dédiée , j'ai dû exprimer le calcul en utilisant uniquement l' ataninstruction:

acos(x) = atan(sqrt(1-x^2)/x)

Voici le code source qui génère le code machine ci-dessus:

__declspec(naked) void __fastcall doit1(int n, std::pair<double, double>* output)
{
    _asm {
        push 0x4f000000;        // [esp] = float representation of 2^32

    myloop:
        rdrand eax;             // generate a random number, -2^31...2^31-1
        push eax;               // convert it
        fild dword ptr [esp];   // to floating-point
        pop eax;                // restore esp
        fdiv dword ptr [esp];   // convert to range -1...1
        neg ecx;
        js myloop;              // do the above 2 times

        // FPU stack contents:  // x           | y
        fld st(0);              // x           | x   | y
        fmul st(0), st;         // x^2         | x   | y
        fld1;                   // 1           | x^2 | x | y
        fsubrp st(1), st;       // 1-x^2       | x   | y
        fsqrt;                  // sqrt(1-x^2) | x   | y
        fxch;                   // x           | sqrt(1-x^2) | y
        fpatan;                 // acos(x)     | y
        fadd st, st(0);         // 2*acos(x)   | y
        fldpi;                  // pi          | 2*acos(x) | y
        fmulp st(2), st;        // 2*acos(x)   | pi*y
        fadd st, st(1);         // output1     | output2
        fstp qword ptr [edx];   // store the numbers
        fstp qword ptr [edx + 8];

        add edx, 16;            // advance the output pointer
        loop myloop;            // loop

        pop eax;                // restore stack pointer
        ret;                    // return
    }
}

Pour générer deux nombres aléatoires, le code utilise une boucle imbriquée. Pour organiser la boucle, le code profite du ecxregistre (compteur de boucle externe) positif. Il annule temporairement ecx, puis recommence pour restaurer la valeur d'origine. L' jsinstruction répète la boucle lorsqu'elle ecxest négative (cela se produit exactement une fois).

anatolyg
la source
J'aime cette réponse pour l'utilisation de l'assemblage IA32! Juste pour dire: ce n'est pas strictement un code machine 386 car 80386 n'a évidemment pas d'instruction rdrand ni besoin d'un coprocesseur FP
user5572685
Ouais, IA32 est un meilleur nom. Pas assez précis mais probablement plus correct que 80386.
anatolyg
7

Pyth, 25 23 22 octets

Un port de la réponse C ++ 11 de rcrmn. Ceci est ma première utilisation de Pyth, et je me suis bien amusé!

VQJ,*y.n0O0.tOZ4,sJ-FJ

Version 23 octets:

VQJ*y.n0O0K.tOZ4+JK-JKd

Coupez un octet en modifiant le programme pour utiliser les plis + sommes et en définissant J sur un tuple, en supprimant K.

Original:

VQJ**2.n0O0K.tO0 4+JK-JKd

Coupez 2 octets grâce à @orlp.

Explication:

VQ                         loop as many times as the input number
  J,                       set J to the following tuple expression
    *y.n0O0                2 * .n0 (pi) * O0 (a random number between 0 and 1)
            .tOZ4          .tOZ 4 (acos of OZ (a random number))
                 ,sJ-FJ    print the sum of J and folding J using subtraction in parenthesis, separated by a comma, followed by another newline
kirbyfan64sos
la source
1
Conseils Pyth: *2_est le même que y_. La variable Zest initialement 0, vous pouvez donc supprimer l'espace en .tO0 4écrivant .tOZ4.
orlp
1
Je compte 25 octets ...
Maltysen
En outre, vous pouvez mieux formater la sortie avec,+JK-JK
Maltysen
@Maltysen Both fixed. Merci!
kirbyfan64sos
@orlp Je n'avais aucune idée yet j'ai oublié Z. Fixé; Merci!
kirbyfan64sos
6

Julia, 48 octets

n->(x=2π*rand(n);y=acos(rand(n));hcat(x+y,x-y))

Il utilise l'algorithme de la méthode 2, comme la plupart des réponses jusqu'à présent. Il crée une fonction lambda qui prend une entrée entière et retourne un tableau nx 2. Pour l'appeler, donnez-lui un nom, par exemple f=n->....

Non golfé + explication:

function f(n::Int64)
    # The rand() function returns uniform random numbers using
    # the Mersenne-Twister algorithm

    # Get n random chord angles
    x = 2π*rand(n)

    # Get n random rotations
    y = acos(rand(n))

    # Bind into a 2D array
    hcat(x+y, x-y)
end

J'aime vraiment l'apparence des visualisations, j'en inclurai donc une. C'est le résultat de f(1000).

Cercle

Alex A.
la source
5

Pyth, 22 octets

Un port de la réponse C ++. J'avais une autre solution de 23 octets (maintenant 22!), Mais c'était presque une copie de la réponse pyth de @ kirbyfan64sos avec des optimisations, j'ai donc dû réfléchir un peu en dehors de la boîte et utiliser de manière créative (ab) l'opérateur de pliage.

m,-Fdsdm,y*.nZOZ.tOZ4Q

Notez que cela ne fonctionne pas pour le moment en raison d'un bogue dans l'opérateur de pliage après l'introduction de reduce2. Je fais une demande de pull.

m             Map    
 ,            Tuple of
  -Fd         Fold subtraction on input
  sd          Fold addition on input
 m      Q     Map over range input
  ,           Tuple           
   y          Double
    *         Product
     .nZ      Pi
     OZ       [0, 1) RNG
  .t  4       Acos
    OZ        [0, 1) RNG

Pour référence, c'était mon autre solution qui fonctionne de la même manière: VQKy*.nZOZJ.tOZ4,+KJ-KJ

Maltysen
la source
Vous avez mal orthographié mon nom d'utilisateur ... :(
kirbyfan64sos
@ kirbyfan64sos derp. Désolé;)
Maltysen
4

IDL, 65 octets

Évidemment, c'est le même algorithme que @rcrmn, même si je l'ai dérivé indépendamment avant de lire leur réponse.

read,n
print,[2,2]#randomu(x,n)*!pi+[-1,1]#acos(randomu(x,n))
end

La fonction randomu d'IDL utilise le Twister de Mersenne, qui a une période de 2 19937 -1.

EDIT: J'ai exécuté 1000 accords via le visualiseur ci-dessus, voici une capture d'écran du résultat:

IDL bertrand

sirpercival
la source
4

C ++ 11, 214 octets

#include<random>
#include<iostream>
#include<cmath>
int main(){using namespace std;int n;cin>>n;random_device r;uniform_real_distribution<> d;for(;n;--n){float x=2*M_PI*d(r),y=acos(d(r));cout<<x+y<<' '<<x-y<<';';}}

Il s'agit donc d'une implémentation directe du bon algorithme à partir de la page wikipedia. Le principal problème ici dans le golf est les noms oh-so-fleaking-long que les classes de générateur aléatoire ont. Mais, contrairement au bon vieux rand, il est au moins convenablement uniforme.

Explication:

#include<random>
#include<iostream>
#include<cmath>
int main()
{
    using namespace std;
    int n;
    cin>>n; // Input number
    random_device r; // Get a random number generator
    uniform_real_distribution<> d;   // Get a uniform distribution of 
                                     // floats between 0 and 1
    for(;n;--n)
    {
        float x = 2*M_PI*d(r),       // x: Chosen radius angle
              y = acos(d(r));        // y: Take the distance from the center and 
                                     // apply it an inverse cosine, to get the rotation

        cout<<x+y<<' '<<x-y<<';';    // Print the two numbers: they are the rotation
                                     // of the radius +/- the rotation extracted from
                                     // the distance to the center
    }
}
rorlork
la source
1
Ce facteur M_PI_2semble suspect. Je pense que ce devrait être 1 à la place.
anatolyg
Oui, tout à fait raison, je vais le réparer maintenant! Merci beaucoup!
rorlork
4

APL, 46 octets

f←{U←⍵ 2⍴∊{(○2×?0)(¯1○?0)}¨⍳⍵⋄⍉2⍵⍴∊(+/U)(-/U)}

Mon tout premier programme APL! Cela peut certainement être considérablement amélioré (car ma compréhension globale de l'APL fait défaut), donc toutes les suggestions seraient fantastiques. Cela crée une fonction fqui prend un entier en entrée, calcule les paires de points d'accord à l'aide de la méthode 2 et imprime chaque paire séparée par une nouvelle ligne.

Vous pouvez l' essayer en ligne !

Explication:

f←{ ⍝ Create the function f which takes an argument ⍵

    ⍝ Define U to be an ⍵ x 2 array of pairs, where the first
    ⍝ item is 2 times a random uniform float (?0) times pi (○)
    ⍝ and the second is the arccosine (¯1○) of another random
    ⍝ uniform float.

    U ← ⍵ 2 ⍴ ∊{(○2×?0)(¯1○?0)}¨⍳⍵

    ⍝ Create a 2 x ⍵ array filled with U[;1]+U[;2] (+/U) and
    ⍝ U[;1]-U[;2] (-/U). Transpose it into an ⍵ x 2 array of
    ⍝ chord point pairs and return it.

    ⍉ 2 ⍵ ⍴ ∊(+/U)(-/U)
}

Remarque: Ma précédente solution de 19 octets n'était pas valide car elle renvoyait (x, y) plutôt que (x + y, xy). La tristesse abonde.

Alex A.
la source
3

Java, 114 octets

n->{for(;n-->0;){double a=2*Math.PI*Math.random(),b=Math.acos(Math.random());System.out.println(a+b+" "+(a-b));}};

Implémentation de base en java. À utiliser comme expression lambda.

Exemple d'utilisation

Le numéro un
la source
Ne pouvez-vous pas réduire la taille en stockant Mathquelque part? Ou quelque chose? (Je ne suis pas un programmeur Java)
Ismael Miguel
@IsmaelMiguel Cela coûterait 2 caractères supplémentaires.
TheNumberOne
Désolé: / Il est tentant d'essayer de réduire le nombre de fois qui Maths'affiche. Que dit la méta sur l'utilisation d'un code pour générer un autre code pour résoudre le problème?
Ismael Miguel
2
@IsmaelMiguel C'est un jeu équitable, même si je serai surpris si vous êtes réellement meilleur en métagolf qu'en golf.
Martin Ender
3

Rubis, 72 octets

Mon premier golf ici! J'ai utilisé le même code que tout le monde, j'espère que ça va

gets.chomp.to_i.times{puts"#{x=2*Math::PI*rand},#{x+2*Math.acos(rand)}"}
rorlok
la source
2

Java, 115 123

C'est fondamentalement le même que la plupart des autres, mais j'ai besoin d'un score Java pour ce trou, alors voici:

void i(int n){for(double x;n-->0;System.out.println(x+2*Math.acos(Math.random())+" "+x))x=2*Math.PI*Math.random();}

1000 exemples d'accords peuvent être trouvés sur pastebin , voici les cinq premiers d'un même run:

8.147304676211474 3.772704020731153
8.201346559916786 3.4066194978900106
4.655131524088468 2.887965593766409
4.710707820868578 3.8493686706403984
3.3839198612642423 1.1604092552846672
Géobits
la source
1

CJam, 24 22 octets

Semblable à d'autres algorithmes, voici une version dans CJam.

{2P*mr_1dmrmC_++]p}ri*

Une entrée de 1000 produit une distribution comme:

entrez la description de l'image ici

Comment ça fonctionne

L'algorithme est simplement x = 2 * Pi * rand(); print [x, x + 2 * acos(rand())]

{                 }ri*        e# Run this loop int(input) times
 2P*mr                        e# x := 2 * Pi * rand()
      _                       e# copy x
       1dmr                   e# y := rand()
           mC                 e# z := acos(y)
             _++              e# o := x + z + z
                ]             e# Wrap x and o in an array
                 p            e# Print the array to STDOUT on a new line

Mise à jour : 2 octets économisés grâce à Martin!

Essayez-le ici

Optimiseur
la source
1

Python 3, 144 117 octets

(merci à Blckknght pour le lambdapointeur)

En utilisant la même méthode que les autres:

import math as m;from random import random as r;f=lambda n:[(x,x+2*m.acos(r()))for x in(2*m.pi*r()for _ in range(n))]

Dans la documentation Python:

Python utilise le Twister de Mersenne comme générateur de base. Il produit des flottants de précision de 53 bits et a une période de 2 19937 -1.

Production

>>> f(10)
[(4.8142617617843415, 0.3926824824852387), (3.713855302706769, 1.4014527571152318), (3.0705105305032188, 0.7693910749957577), (1.3583477245841715, 0.9120275474824304), (3.8977143863671646, 1.3309852045392736), (0.9047010644291349, 0.6884780437147916), (3.333698164797664, 1.116653229885653), (3.0027328050516493, 0.6695430795843016), (5.258167740541786, 1.1524381034989306), (4.86435124286598, 1.5676690324824722)]

Etc.

Visualisation

visualisation

Portes Zach
la source
Vous pouvez économiser environ 20 octets si vous utilisez un lambda pour la fonction et renvoyez une liste de compréhension (avec une expression de générateur interne):f=lambda n:[(x,x+2*m.acos(r()))for x in(2*m.pi*r()for _ in range(n))]
Blckknght
Ah, j'avais une lambda à l'origine. Je suppose que je n'ai pas pensé à doubler la compréhension des listes. Merci! @Blckknght
Zach Gates
Peut être raccourci à 109 octets en jouant avec les importations: tio.run/#python2
Triggernometry
1

Perl, 60

#!perl -l
use Math::Trig;print$x=2*pi*rand,$",$x+2*acos rand for 1..<>
nutki
la source
1

R, 60 56 53 49 octets

4 octets supplémentaires grâce à @JayCe et en le changeant en fonction.

En utilisant la même formule de base que les autres. R utilise la méthode Mersenne-Twister par défaut, mais d'autres peuvent être définies. Génère une liste séparée par des espaces.

function(n,y=acos(runif(n)))runif(n)*2*pi+c(y,-y)

Essayez-le en ligne!

MickyT
la source
Salut Micky, tu peux économiser 4 octets en en faisant une fonction et non en définissant x.
JayCe
@JayCe c'est beaucoup mieux merci
MickyT
0

SmileBASIC, 62 octets

INPUT N
FOR I=1TO N
X=PI()*2*RNDF()Y=ACOS(RNDF())?X+Y,X-Y
NEXT
12Me21
la source