Zufällige 3D Vektoren innerhalb einer Kugel

[SeppJ]

Heaviside-Stufenfunktion .

Die Verteilungsfunktion ergibt sich dadurch, dass die Kugelflächen wie R^2 gehen, die Verteilung der Radien die wir ziehen, muss daher auch proportional zu R^2 sein, damit man hinterher eine gleichmäßige Verteilung im Raum bekommt. Die Werte sollen durch R1 und R2 eingeschränkt sein, daher die Stufenfunktionen. Und die Gesamtwahrscheinlichkeit muss 1 sein, das ist die Normierung unter dem Bruch.

[*Rewind*]

SeppJ schrieb:

Nun, und das habe ich dann umgekehrt . Und während draußen Gewitterstürme tobten und der strenge Geruch von Formaldehyd von den eingemachten Organen im Regal aufstieg, erwachte so in meinem Labor beim kalten Lichte der Gasüberschlagslampen das MONSTER . Muahuahuahuahua!!! Sie haben mich wahnsinnig genannt!! Ich werde ihnen zeigen, wer hier wahnsinnig ist!!!!!

Zu viel Sinn City geguckt?

[volkard]

SeppJ schrieb:

Ob sich die Rechnung nun lohnt sei dahingestellt. Bei der Wegwerfmethode kann man höchstens (4/3 r^3)/(2r)3 = 1/6 der Zufallszahlen benutzen, wenn man eine Vollkugel hat.

Da ist Π verlorengegangen.
(4/3 Π r^3)/(2r)3 = 1/6 Π = 52%

Im Fall, daß die innere Kugel den halben Radius hat, sind's noch 46%, wie pumuckl auf Seite 2 festgestellt hatte.

[otze]

Werner Salomon schrieb:

was dann kommt habe ich zwar nicht verstanden - insbesondere die Ausgangsgleichung

Okay, ich geb mir mal mehr Muehe

Eine vereinfachende Definition zuerst, ich bezeichne die Oberflaeche einer Kugel mit Radius r mit:
A(r) = 4* pi * r^2

was ich nun machen will ist, wie du bereits sagtest: die Auftrittswahrscheinlichkeit eines r muss proportional zu A(r) sein. Da wir zugleich die Nebenbedingung haben, dass das \int_a^1 p(r) dr = 1 sein muss, ergibt sich der Proportionalitaetsfaktor als die Normalisierungskonstante.

p(r) = A(r) / \int_a^1 A(r) dr = A(r)/ (4/3*pi - 4/3*pi*a^3)

und dann nurnoch A(r) einsetzen.

Jetzt bin ich aber auf deinen Gedankengang gespannt. Du scheinst irgendwie straight forward gewusst zu haben, wie die r Verteilt sein muessen.

[brotbernd]

volkard schrieb:

SeppJ schrieb:

Ob sich die Rechnung nun lohnt sei dahingestellt. Bei der Wegwerfmethode kann man höchstens (4/3 r^3)/(2r)3 = 1/6 der Zufallszahlen benutzen, wenn man eine Vollkugel hat.

Da ist Π verlorengegangen.
(4/3 Π r^3)/(2r)3 = 1/6 Π = 52%

Im Fall, daß die innere Kugel den halben Radius hat, sind's noch 46%, wie pumuckl auf Seite 2 festgestellt hatte.

Zieht man in einer Würfelschale, die die Kugelschale enthält sind es auch wieder 52%. Man kann also auch bei sehr dünnen Schalen genauso viele Zufallszahlen nutzen.

[Dobi]

Wenn die Kugelschale sehr dünn wird, wird man auch bei der Würfelschalenmethode ziemlich viele Zahlen wegwerfen müssen. Man könnte so eine dünne Kugelschale dann mit Würfeln approximieren. Ist dann dann aber vermutlich kaum noch einfacher als die hier diskutierte elegante/komplizierte Lösung.

[brotbernd]

Dobi schrieb:

Wenn die Kugelschale sehr dünn wird, wird man auch bei der Würfelschalenmethode ziemlich viele Zahlen wegwerfen müssen. Man könnte so eine dünne Kugelschale dann mit Würfeln approximieren. Ist dann dann aber vermutlich kaum noch einfacher als die hier diskutierte elegante/komplizierte Lösung.

Es ist immer 0,52.
4/3 * pi * (r1^3 - r2^3) / (8 * (r1^3 - r2^3) = 1/6 pi

Edit: Falscher Fehler. Der innere Würfel ist natürlich innerhalb der inneren Kugel. Das war also quatsch...

[brotbernd]

...

[otze]

Was dem Thread jetzt fehlt ist die Verallgemeinerung auf Hyperkugeln.

[krümelkacker]

Der Trick mit den Gaußverteilungen für eine Gleichverteilung auf einer Kugeloberfläche klappt immer -- unabhängig von der Dimension n. Für den Radius r benötigt man dann aber eine W'keitsdichte P(r), die von n abhängt:
$P(r) \propto r^{n-1}$

[krümelkacker]

Ich sehe gerade, dass Boost sogar schon eine Gleichverteilung über einer Kugeloberfläche anbietet mit demselben Trick, den ich auch jetzt schon mehrfach vorschlug: boost::random::uniform_on_sphere.

Das mit dem Radius ist auch ganz einfach: Für den n-dimensionalen Raum zwischen zwei Kugeloberflächen mit Radius a und b: Einfach eine Zahl aus dem Intervall [a^n, b^n] per Gleichverteilung auswählen und daraus die n-te Wurzel ziehen. Fertig.

[Dobi]

Kann man sich den ganzen Aufwand mit der Gleichverteilung bei der Polarkoordinatendarstellung nicht einfach sparen, indem mal einen Zufallspunkt im Würfel erzeugt, und dessen Abstand zum Mittelpunkt dann mit der von krümelkacker beschriebenen Methode für den Radius korrigiert?

[krümelkacker]

Dobi schrieb:

Kann man sich den ganzen Aufwand mit der Gleichverteilung bei der Polarkoordinatendarstellung nicht einfach sparen, indem mal einen Zufallspunkt im Würfel erzeugt, und dessen Abstand zum Mittelpunkt dann mit der von krümelkacker beschriebenen Methode für den Radius korrigiert?

Nein. Das, was Du dafür brauchst ist eine Normalverteilung für jede Achse, keine Gleichverteilung im Würfel. Sonst bekommst in den Richtungen, wo der Würfel seine Ecken hat, Häufungen. Aber scheinbar muss ich das 5mal sagen, dass das mit der Normalverteilung elegant einfach wird.... :p

[Dobi]

Ah, OK, stimmt.
Man könnste es dann nur so machen, dass man alle Punkte wegwirft, die ausserhalb der Kugel liegen. Damit würde man dann auch bei dünnen Kugelschalen nicht mehr Verlust haben als bei einer Vollkugel.

Vermutlich musst du es für mich auch noch ein 6stes Mal sagen, weil ich noch nicht verstanden hab, wie das funktioniert, und deshalb nach was mathematisch einfacherem suche.

[Werner Salomon]

otze schrieb:

Jetzt bin ich aber auf deinen Gedankengang gespannt. Du scheinst irgendwie straight forward gewusst zu haben, wie die r Verteilt sein muessen.

Hallo miteinander,

entschuldigt bitte, dass meine Antwort so lange auf sich warten lies, aber ich bin z.Zt. offline stark beschäftigt.

@otze Ich habe Deine Herleitung inzwischen verstanden - mit funktionierendem Latex wäre es auch ein wenig einfacher gewesen.

Mein Gedankengang war einfach ein zufälliges Volumen zwischen zwei gegebenen Volumina - also dem der inneren Kugel und dem der äußeren Kugel - zu wählen und dann daraus den Radius zu bestimmen. Das ist schon alles. Finde ich eher noch naheliegender, als Deine Herleitung.
Bei der Herleitung des Breitenwinkels (in meinem Code der Winkel 'phi' in class Sphere) bin ich genauso vorgegangen. Ich wähle eine zufällige Oberfläche und bestimmen daraus den Winkel (bzw.den sin des Winkels).

Was SeppJ dort gemacht hat, habe ich nicht entschlüsseln können.

krümelkacker schrieb:

Ich sehe gerade, dass Boost sogar schon eine Gleichverteilung über einer Kugeloberfläche anbietet mit demselben Trick, den ich auch jetzt schon mehrfach vorschlug: boost::random::uniform_on_sphere.

Das mit dem Radius ist auch ganz einfach: Für den n-dimensionalen Raum zwischen zwei Kugeloberflächen mit Radius a und b: Einfach eine Zahl aus dem Intervall [a^n, b^n] per Gleichverteilung auswählen und daraus die n-te Wurzel ziehen. Fertig.

@krümelkacker: Danke für diesen Hinweis. Ich habe in Deinem Beitrag zu spät das Wort normalverteilt gelesen, und hatte ihn daher daher zunächst ignoriert.

Kann denn jemand (nachvollziehbar!) beweisen, dass dieser Algorithmus tatsächlich zu gleich verteilten Punkten auf einer Kugelschale führt. Ich habe ein kleines Testprogramm geschrieben - rein numerisch ist das so.

Faszinierend ist, dass SeppJ, ich und auch boost::uniform_sphere zu identischen Ergebnissen kommen; obwohl in jedem Fall der gewählte Algorithmus völlig unterschiedlich ist.

Gruß
Werner

[krümelkacker]

Werner Salomon schrieb:

Kann denn jemand (nachvollziehbar!) beweisen, dass dieser Algorithmus tatsächlich zu gleich verteilten Punkten auf einer Kugelschale führt.

Ja. Sei N die Anzahl der Dimensionen und X_k mit k=1...N unabhängige Zufallsvariablen mit derselben Verteilung: normalverteilt mit Erwartungswert 0 und Streuung 1. Wenn wir X als Vektor [X_1 X_2 X_3 ...]^T betrachten, lässt sich aufgrund der Unabhängigkeit der einzelnen Komponenten die Wahrscheinlichkeitsdichte für X per Multiplikation angeben:
P(x) = (1/sqrt(2*pi))^N * exp(-1/2 x_1^2) * exp(-1/2 x_2^2) * ... * exp(-1/2 x_N^2)
Eines der Potenzgesetze sagt nun, dass wie das so umschreiben können
P(x) = (1/sqrt(2*pi))^N * exp(-1/2 (x_1^2+x_22+x_3^2+...+x_N2) )
Und x_1^2+x_22+x_3^2+...+x_N2 ist einfach das Quadrat der Läng des Vektors. Man sieht also, dass die Wahrscheinlichkeitsdichte nicht mehr von der Richtung des Vektors x abhängt, sondern nur noch von seiner 2-Norm. Die Wahrscheinlichkeitsverteilung auf der Kugeloberfläche nach einem Normalisierungsschritt ist also uniform.