Benutzer:SiriusB/Formelsammlung

Auf dieser Seite werden Hintergrundinformationen, insbesondere Formeln, Tabellen und ähnliches aufgeführt werden, die den Rahmen der entsprechenden Wikipedia-Artikeln sprengen bzw. die Artikel unübersichtlicher machen würden. Da Interessierte sich ggf. doch für diese Details interessieren könnten, will ich sie nicht unter dem Scheffel behalten.

Ergänzungen zum Artikel Atombombenexplosion

Freiluftexplosion

Druckwellen-Standardkurve

Abb. 1: DNA-Standardkurve für Überdruck OP und davon abgeleitetem dynamischen Druck DP

Die frühere US-Behörde Defense Nuclear Agency (DNA) hat um 1984 ein Modell zur rechnerischen Abschaetzung von Explosionsdruckwellen entwickelt, welches die Abhängkeit der Überdruckspitze als Funktion der Distanz sowie der Detonationshöhe und der Sprengenergie liefert. Obgleich das Modell ursprünglich auf Nuklearexplosionen orientiert ist, kann es über Skalierungsregeln auch auf nahezu beliebige andere Explosionstypen angewendet werden, sofern die Explosion von einer Punktquelle ausgeht. Das Modell, das in Gestalt eines DOS-Programms BLAST vorliegt, unterliegt keiner Geheimhaltung. Es baut im Wesentlichen auf der Rankine-Hugoniot-Gleichung sowie empirischer Fits auf der Basis von Atomtest-Daten.

Der in Abb. 1 dargestellten Standardkurve für die Druckwelle einer 1-Kilotonnen-Atombombenexplosion liegt die folgende, von der entwickelte Beziehung zwischen dem Abstand R vom Explosionstentrum und dem Druckpegel OP zugrunde, wobei von einer Freiluftexplosion in einer homogenen unbegrenzten Atmosphäre unter Meeresniveaubedingungen (P = P₀ = 101325 Pa und T = T₀ = 288{,}15 K) ausgegangen wird:

{\mathit {OP}}_{\mathrm {DNA} }={\frac {3.04\times 10^{11}}{R^{3}}}+{\frac {1.13\times 10^{9}}{R^{2}}}+{\frac {7.9\times 10^{6}}{R{\sqrt {\ln \left({\frac {R}{445.42}}+3\exp \left(-{\frac {1}{3}}{\sqrt {\frac {R}{445.42}}}\right)\right)}}}}

Dabei ist R in m einzusetzen (der Übersichtlichkeit wegen wurde auf die Einheiten-Divisoren in der Formel verzichtet), und das Ergebnis ist OP in Pa. Der dynamische Druck ergibt sich aus

{\mathit {DP}}={\frac {n-1}{2}}\,{\mathit {OP}}\,,

wobei n das Dichteverhältnis vor und hinter der Stoßfront und P der Druck der ungestörten Atmosphäre ist. Für Luft ist

n={\frac {1+\mu _{\mathrm {s} }\left({\frac {\mathrm {OP} }{P}}+1\right)}{\mu +\left({\frac {\mathrm {OP} }{P}}+1\right)}}\quad {\mbox{mit}}\quad \mu ={\frac {\kappa +1}{\kappa -1}}\quad {\mbox{und}}\quad \mu _{\mathrm {s} }={\frac {\kappa _{\mathrm {s} }+1}{\kappa _{\mathrm {s} }-1}}\,.

κ ist der Adiabatenexponent, und der Index s deutet an. dass κ infolge der Stoßerhitzung hinter der Stoßfront nicht mehr den klassischen Wert für Luft von 1,402 besitzt, da durch Ionisation mehr Freiheitsgrade hinzukommen. Für die weitere Berechnung werden zunächst einige temporäre Variablen definiert, über die dann die Korrektur für κ resultiert:

x={\frac {\mathit {OP}}{P}}+1\,;\quad y=10^{-12}\,x^{6}\,;\quad z=\ln x-{\frac {0{,}47y}{100+y}}

Dann ist

\kappa _{\mathrm {s} }=1{,}402-{\frac {3{,}4\times 10^{-4}\,z^{4}}{1+2{,}22\times 10^{-5}\,z^{6}}}\,.

Unterhalb von etwa 1000 kPa ist der Korrekturterm vernachlässigbar. Aus den hier berechneten Größen folgt auch der Normalreflexionsfaktor R_n, der die Druckerhöhung bei senkrechter Reflexion wiedergibt (trivialerweise 2 bei gewöhnlichen Schallwellen):

R_{\mathrm {n} }=2+{\frac {(\kappa _{\mathrm {s} }+1)(n-1)}{2}}\quad {\mbox{und}}\quad {\mathit {OP}}_{\mathrm {GZ} }=R_{\mathrm {n} }\cdot {\mathit {OP}}_{\mathrm {DNA} }(R=H)

Mit diesen Resultaten ergibt sich aus den Rankine-Hugoniot-Gleichungen die Ausbreitungsgeschwindigkeit und die Windgeschwindigkeit der Stoßfront:

v_{\mathrm {s} }=c_{\mathrm {s} }\cdot {\sqrt {{\frac {(\kappa _{\mathrm {s} }+1){\mathit {OP}}_{\mathrm {DNA} }}{2\kappa _{\mathrm {s} }P}}+1}}\ ;\quad c_{\mathrm {s} }={\mbox{Schallgeschwindigkeit}}

w_{\mathrm {s} }=v_{\mathrm {s} }\cdot \left(1-{\frac {1}{n}}\right)

Daraus ergibt sich durch Integration der reziproken Geschwindigkeit auch die Laufzeit der Druckwelle zu einem bestimmten Radius.

Skalierungsfaktoren

Der wichtigste Skalierungsfaktor ist die Sprengkraft- oder Sachs-Skalierung. Die Standardkurven sind für 1 kT definiert; für beliebige Energien W (gleicher mechanischer Anteil vorausgesetzt) ist für beliebige Längen- und Zeitgrößen L die Gleichung

L(W)=L(1\,\mathrm {kT} )\cdot \left({W \over \mathrm {kT} }\right)^{1/3}

anzuwenden. Für gegebenen Luftdruck P und Temperatur T gehen zudem

S_{P}={\frac {P}{P_{0}}}\ ;\quad S_{T}={\frac {T}{T_{0}}}\ ;

{\mbox{Länge:}}\ S_{L}=S_{P}^{-1/3}\ ;\quad {\mbox{Dauer:}}\ S_{D}={\frac {S_{L}}{\sqrt {S_{T}}}}

in die korrigierte 1-kT-Druck- und Zeitkurve ein:

{\mathit {OP}}(R,P)=S_{P}\cdot {\mathit {OP}}(R/S_{L},P_{0})\ ;\quad t(R,P,T)=S_{D}\cdot t(R/S_{L},P_{0},T_{0})\ .

Der Zeitskalierungsfaktor ist für die Laufzeit der Druckwelle von Bedeutung. Über Druck bzw. Temperatur wird auch die weiter unten benötigte Schallgeschwindigkeit skaliert.

Modell für Luftexplosion

Abb. 2: OP-Zonen für Luftexplosionen nach dem Modell des Programms BLAST.

Das DNA-Modell beschreibt auch die Druckwelle bei Luftexplosionen, also unter Berücksichtigung der Reflexion an der Oberfläche.

Wendet man auf diese Formeln einen numerischen Solver an, so erhält man für diskrete Werte von OP Kurven, wie sie in Abb. 2 dargestellt sind.

Bei der Berechnung sind zwei Regime zu unterscheiden, nämlich das Regime der regulären Reflexion und das der Mach-Reflexion. Für das erstere benötigt man neben dem oben berechneten Normalreflexionsfaktor R_n noch einige temporäre Variablen und Koeffizienten, um dann den Überdruck für reguläre Reflexion, OP_reg berechnen zu können. Dabei wurde gelegentlich von Druckgrößen nur der Betrag verwendet (z.B. in Exponenten oder bei nicht-ganzzahliger Potenzierung usw.), der Übersichtlichkeit wegen aber auf zusätzliche Symbole verzichtet.

Höhenkoeffizienten

Seien H und GR für eine 1-kT-Explosion gegeben. Dann seien der Einfallswinkel α der Primärfront sowie einige temporäre Variablen definiert:

\alpha =\arctan {\frac {H}{\mathit {GR}}}

T={\frac {340}{OP_{\mathrm {DNA} }^{0{,}55}}}\ ;\quad U=\left({\frac {7782}{OP_{\mathrm {DNA} }^{0{,}7}}}+0{,}9\right)^{-1}

{\tilde {U}}=\left({\frac {7473}{OP_{\mathrm {DNA} }^{0{,}5}}}+6{,}6\right)^{-1}\ ;\quad V=\left({\frac {647}{OP_{\mathrm {DNA} }^{0{,}8}}}+{\tilde {U}}\right)^{-1}

Grenzwinkel α_M zwischen regulärer und Mach-Region und (Winkel-)Breite β der Zone, wo die Wellen verschmelzen, Hilfsvariablen und Schaltparameter σ:

\alpha _{\mathrm {M} }=\arctan {\frac {1}{T+U}}\ ;\quad \beta =\arctan {\frac {1}{T+V}}\ ;\quad s={\frac {\alpha -\alpha _{\mathrm {M} }}{\beta }}

s_{0}=\max(\min(s,1),-1)\ ;\quad \sigma =0{,}5\left(\sin \left({\frac {\pi }{2}}\,s_{0}\right)+1\right)

Reguläre Reflexion

f={\frac {OP_{\mathrm {DNA} }}{75842\,\mathrm {Pa} }}\ ;\quad g={\frac {f^{6}\left(1{,}2+0{,}07\,f^{0{,}5}\right)}{f^{6}+1}}\ ,

und daraus

OP_{\mathrm {reg} }=OP_{\mathrm {DNA} }\left((R_{\mathrm {n} }-2)\sin ^{g}\alpha +2\right)\ .

Mach-Reflexion

Zunächst wieder ein paar Zwischenvariablen:

A=\min(3{,}7-0{,}94\,\ln({\mathit {GR}}),\ 0{,}7)\ ;\quad B=0{,}77\,\ln({\mathit {GR}})-{\frac {18}{\mathit {GR}}}-3{,}8

C=\max(A,B)

Nun setze 2^(-1/3)*GR anstelle von R in die DNA-Formel OP(R) ein und setze das daraus erhaltene Ergebnis OP₁ hier ein:

{\mathit {OP}}_{\mathrm {mach} }={\frac {{\mathit {OP}}_{1}}{1-C\sin \alpha }}

Der Abstand GR_M, an dem die Verschmelzung von direkter und reflektierter Welle einsetzt, lässt sich auch annähern durch

GR_{\mathrm {M} }={\frac {H^{2{,}5}}{5822}}+2{,}09\,H^{0{,}75}\ .

Gesamtamplitude

Der Gesamt-Überdruck OP_air ist dann

OP_{\mathrm {air} }=\sigma OP_{\mathrm {reg} }+(1-\sigma )OP_{\mathrm {mach} }\ .

Der dynamische Druck folgt dann mit dem Dichteverhältnis in der resultierenden Druckfront (anstelle der Standarddruckwelle) n_a=n(OP_air):

DP_{\mathrm {air} }=OP_{\mathrm {air} }(n_{a}-1)\left(1-\sigma \sin ^{2}\alpha \right)/2

Die Ausbreitungsgeschwindigkeit der Welle kann nun mit OP_air statt OP_DNA geschehen; die Windgeschwindigkeit w_s berechnet man sinnvollerweise aus DP_air, weil dort die Restriktion auf die Horizontale bereits berücksichtigt worden ist:

w_{s}={\sqrt {2{\frac {{\mathit {DP}}_{\mathrm {air} }}{\rho _{2}}}}}\ ,

wobei ρ₂ = n*ρ_0 die Dichte hinter der Stoßfront ist.

Der Fit

Die in dem Artikel ebenfalls angegebenen Gleichungen für die optimale Explosionshöhe sowie der Grundradius GR für den gewählten Überdruck

{H_{\mathit {OP}} \over \mathrm {m} }\approx 16100\,\left({{\mathit {OP}} \over \mathrm {Pa} }\right)^{-0{,}38}=560\,\left({{\mathit {OP}} \over \mathrm {PSI} }\right)^{-0{,}38}\,.

{{\mathit {GR}}_{\mathit {OP}} \over \mathrm {m} }\approx \left(\left(9{,}54\cdot 10^{-3}\,{H_{\mathit {OP}} \over \mathrm {m} }^{1{,}95}\right)^{4}+\left(3{,}01\,{H_{\mathit {OP}} \over \mathrm {m} }^{0{,}75}\right)^{4}\right)^{1/4}\,.

{{\mathit {GR}}_{{\mathit {OP}},H=0} \over \mathrm {m} }\approx \left(\left(2{,}738\cdot 10^{6}\,{{\mathit {OP}} \over \mathrm {Pa} }^{-0{,}879}\right)^{1{,}37}+\left(1{,}012\cdot 10^{4}\,{{\mathit {OP}} \over \mathrm {Pa} }^{-0{,}342}\right)^{1{,}37}\right)^{1/1{,}37}\,.

entstammen nicht BLAST, sondern sind auf ähnliche Weise wie der GR(OP)-Fit für Bodenexplosionen via Least-Squares-Methode zustande gekommen. Allerdings wurden nicht allein Resultate aus BLAST berücksichtigt, sondern ebenso Druckkurven von Glasstone und Dolan, The Effects of Nuclear Weapons, 1977, welche geringfügig von ersteren abweichen. Der Fit stellt einen Kompromiss zwischen beiden dar. Sämtliche Genauigkeitsangaben im Artikel beziehen sich auf Abweichungen zu diesen Quellen.

Die zweite und dritte Formel sind Näherungen auf der Basis einer Potenzsumme mit reellen Exponenten, die in doppelt logarithmischer Darstellung für sehr kleine und sehr große OP lineare Asymptoten mit fließendem Übergang ergibt. Die Krümmung der Übergangskurve wurde durch die Potenz einer p-Norm angepasst. Für die Bodenexplosion wurde einfach die Umkehrbeziehung der DNA-Formel gefittet, allerdings für 2 kT, was einer 1-kT-Bodenexplosion unter Berücksichtigung der Reflexion der Druckwelle an der Oberfläche entspricht. Alle diese Gleichungen beziehen sich auf 1 kT; für andere Energien skalieren sie mit der Kubikwurzel des Betrages in kT.

Laufzeit der Druckwelle

Um die Laufzeit der Druckwelle zu errechnen, ist die reziproke Geschwindigkeit der Stoßfront über die (auf die tatsächliche Sprengkraft skalierte) Entfernung zu integrieren. Bei Luftexplosionen sind hier zwei Wegabschnitte zu unterscheiden:

Innerhalb der Mach-Zone OP_M ist die Laufzeit mit der einer Freiluftexplosion (skalierte DNA-Standardexplosion) identisch. R ist die Raumdiagonale vom Explosionszentrum zum Messpunkt. Die Refraktion der Welle in der inhomogenen Atmosphäre kann normalerweise vernachlässigt und von einer geradlinigen Ausbreitung ausgegangen werden.
Außerhalb der Mach-Zone bewegt sich die Stoßfront horizontal, daher ist erst die Strecke vom Zentrum zum Rand des Mach-Radius als Freiluftexplosion und von dort zum Messpunkt als Bodenexplosion mit der reflexionsverstärkten Druckwelle zu betrachten, welche sich schneller fortpflanzt als die unverstärkte Welle.

Das DNA-Modell verwendet einen Näherungsfit, der ohne rechenaufwendige Integration auskommt. Zunächst wird eine Laufzeitfunktion für eine 1-kT-Freiluftexplosion definiert:

t_{\mathrm {free} }(R)={\frac {R^{2}(R+6{,}7)}{340{,}5\,R^{2}+73200\,R+7{,}12\cdot 10^{6}}}

Für Freiluftexplosionen müssen jetzt nur R und t_free um die Kubikwurzel von W/kT skaliert werden. Für Luftexplosionen wird ein weiterer Korrekturfaktor benötigt:

\nu =\left\{{\begin{matrix}1&{\mbox{wenn}}\ GR\leq GR_{\mathrm {M} }\\{\sqrt[{3}]{2}}-\left({\sqrt[{3}]{2}}-1\right){\frac {GR_{\mathrm {M} }}{GR}}&{\mbox{sonst}}\end{matrix}}\right.

Nun wird in die Laufzeitformel einfach R_ν=R/ν statt R eingesetzt und das Ergebnis mit W^(1/3) skaliert. Für von den Standardbedingungen stark abweichende Bedingungen können die obigen Skalierungsfaktoren auch für Luftexplosionen angewendet werden. Für große Detonationshöhen muss die numerische Integration mit anhand lokaler Längen- und Zeitskalierungen korrigierter Ausbreitungsgeschwindigkeit verwendet werden.

Hinweis: Die Variablen haben in BLAST teilweise andere Bezeichnungen.

Feuerball-Entwicklung

Abb 3.: Effektivtemperatur und Durchmesser des Feuerballs einer 20-kT-Explosion als Funktion der Zeit.

Zur näheren Illustration werde ich noch Diagramme zur Feuerballentwicklung in geeignete Artikel (Atombombenexplosion oder Feuerball) verlinken. Die Ausdehnung des Feuerballs einer 20-Kilotonnen-Luftexplosion verläuft gemäß Abb. 3 unten, während die Effektivtemperatur dem oberen Graphen in Abb. 3 folgt. Dabei ist die beobachtete Effektivtemperatur (orange Teilkurve) niedriger als die theoretische (dunkle Kurve), da ein Teil der Strahlung durch den Feuerball umgebenden Nebel aus Ozon und Stickoxiden absorbiert wird. Hier wird der Unterschied zur wahren Temperatur deutlich. Die Leuchtkraft ergibt sich dann durch das Produkt der Feuerballoberfläche und der Intensität nach dem Stefan-Boltzmann-Gesetz, wobei für die beobachtete Temperatur eher der Feuerballdurchmesser zum Zeitpunkt des Übergangs und darüber (ca. 0,005 s nach Zündung, entsprechend ca. 140 m) zu verwenden ist, da dies etwa dem Radius der absorbierenden Nebelschicht entspricht, die schließlich vom expandierenden Feuerball eingeholt wird.

Der Durchmesser des Feuerballs zwischen 0,0001 und 10 Sekunden lässt sich genähert durch die Beziehung

D(t)=\left(\left(1010\cdot t^{0.384}\right)^{-3}+475^{-3}\right)^{-1/3}

darstellen. Die Skalierung auf andere Sprengenergien ist nicht ganz so einfach wie bei der Druckwelle; Gründe und Näherungen hierfür sind im Artikel Atombombenexplosion angegeben.

Pilzwolke

Endhöhe und -radius der Pilzwolke

Die Endhöhe der Pilzwolke ist stark von der Energie abhängig; folgende Kurven wurden mit folgender Funktion an das Diagramm der Originalquelle angepasst:

H(W)\approx 3.0\,\mathrm {km} \cdot 10^{0.006941{\mathit {LW}}^{4}-0.06216{\mathit {LW}}^{3}+0.1526{\mathit {LW}}^{2}+0.1878{\mathit {LW}}}

,

R(W)\approx 0.6\,\mathrm {km} \cdot 10^{0.0137{\mathit {LW}}^{3}-0.0358{\mathit {LW}}^{2}+0.37{\mathit {LW}}}

.

Dabei ist ${\mathit {LW}}=\log _{10}W$ . Für große Energien (> 1 MT) tritt "overshoot" auf, d.h. die Wolke steigt zunächst noch höher, stabilisiert sich dann aber auf einer etwas geringeren Endhöhe. Das Bild zeigt die Endhöhen nach ca. 10 min.

Tabelle der typischen Auswirkungen

Sichtweite: 20 km
Explosionshöhe: Optimiert für 15 PSI (Fit)
Druck im GZ: 42 PSI

Wirkung bis GR / km	Sprengenergie / Explosionshöhe
Wirkung bis GR / km	1 kT / 200 m	20 kT / 540 m	150 kT / 1,1 km	1 MT / 2,0 km	20 MT / 5,4 km
Druckwirkung
Totale Zerstörung (20 PSI)	0,24	0,64	1,26	2,4	6,4
Fast totale Zerstörung (10 PSI)	0,40	1,07	1,8	4,0	10,7
Weitgehende Zerstörung (5 PSI)	0,61	1,66	3,2	6,2	16,6
Schwere zivile Schäden (3 PSI)	0,84	2,3	4,5	8,4	23
Mäßige zivile Schäden (1 PSI)	1,73	4,7	9,2	17,3	47
Themische Wirkung
Starke Brandwirkung	0,5	2,0	5	10	30
Verbrennungen 3. Grades	0,6	2,5	5,9	12,5	38
Verbrennungen 2. Grades	0,8	3,2	7,3	15,0	44
Verbrennungen 1. Grades	1,1	4,2	9,5	19,2	53
Wirkung der ionisierenden Direktstrahlung (Raumdiagonale¹ SR / km)
Tödliche² Neutronendosis	0,7	1,2	1,6	2,1	2,8
Akut schädliche² Neutronendosis	1,1	1,6	2,0	2,5	3,4
Tödliche γ-Dosis	0,6	1,1	1,6	2,2	4,7
Akut schädliche γ-Dosis	1,0	1,6	2,3	2,9	5,4
Tödliche Gesamtdosis	0,85	1,4	1,8	2,3	4,7
Akut schädliche Gesamtosis	1,2	1,8	2,3	2,9	5,4

¹) Die Angaben der Strahlendosen beziehen sich auf die raumdiagonale Luftstrecke (slant range, SR) bei 90% der Luftdichte in Meereshöhe (also 0,9·1,225 kg/m³ = 1,1 kg/m³), was einer Höhe von 1100 m entspricht. Die Grunddistanz GR ist entsprechend kleiner oder besitzt gar keine reelle Lösung, falls SR < H. Die tatsächliche mittlere Luftdichte weicht um nicht mehr als etwa 15% von dem angenommenen Wert ab.

²) "Akut schädlich" bezieht sich auf eine Dosis von 1 Gy, "Tödlich" auf 10 Gy. Das Programm WE (Weapon Effects; gleiche Quelle wie BLAST) gibt di Werte zwar in rad an (1 rad = 0.01 Gy), berechnet aber Korrekturfaktoren für "Tissue" (Gewebe), bei denen ich zumindest vermute, dass es sich um eine Art biologische Gewichtung handelt (sollte jemand aus den in der eingebauten Online-Hilfe des DOS-Programms etwas anderes herauslesen, möge er/sie mich bitte berichtigen). Die Werte sollten also zumindest in der Größenordnung denen in Sv entsprechen. Die entsprechenden Wirkungen sind unter Strahlenkrankheit nachzulesen.

Quellen

Programm BLAST und andere

Ergänzungen zum Artikel Sonne

Mitte-Rand-Variation

Mitte-Rand-Variations in Abhängigkeit des Austrittswinkels θ bzw. des Abstands ρ vom Zentrum der projizierten Sonnenscheibe (1 = Rand):

{I(\theta ) \over I(0)}={\frac {1+\beta \cos \theta }{1+\beta }}\ ;\quad {I(\rho ) \over I(0)}={\frac {1+\beta {\sqrt {1-\rho ^{2}}}}{1+\beta }}

Näherungsformel für β (gefittet anhand von Tabellendaten)

\ln \beta =c_{0}+c_{1}\mu +c_{2}\mu ^{2}+c_{3}\mu ^{3}+c_{4}\mu ^{4}\ ;\quad \mu =\ln {\frac {\lambda }{10^{-10}\,\mathrm {m} }}

Tabelle der Polynomkoeffozienten c_i

c₀  = 6352,1488616
c₁  =-2701,6346146
c₂  =  431,30312867
c₃  =  -30,621513866
c₄  =    0,81542004069

Mittlere Sonnenintensität

Die gemittelte Intensität <I> erhält man durch (analytisch mögliche) Integration über die Scheibenfläche und Normierung. Der Normierungsfaktor f definiert und durch β darstellbar durch

{\frac {\langle I\rangle }{I(0)}}={\frac {1}{f}}\quad \Rightarrow \quad f={\frac {1+\beta }{1+{\frac {2}{3}}\beta }}\ .

Die Helligkeit B der Korona in Einheiten der Sonnenhelligkeit ist dann

B=2f\cdot \int _{\rho }^{\infty }I(r)r\,\mathrm {d} r

=2\cdot 10^{-6}f\cdot \left({\frac {3{,}670}{16}}\rho ^{-16}+{\frac {1{,}939}{5{,}8}}\rho ^{-5{,}8}+{\frac {0{,}0551}{0{,}5}}\rho ^{-0{,}5}\right)\cdot B_{\mathrm {Sun} }

Tabelle der MRV-Koeffizienten

            nm      beta     f        I1/I0    I1/Imittel
Grenze UV   380,0   9,5411   1,4321   0,0949   0,1359
Violett     422,5   5,2556   1,3890   0,1599   0,2220
Blau        473,5   3,1658   1,3393   0,2401   0,3215
Grün        534,5   2,0850   1,2908   0,3242   0,4184
V-Filter    555,0   1,8697   1,2774   0,3485   0,4451
Gelb        587,0   1,6141   1,2592   0,3825   0,4817
Orange      609,5   1,4760   1,2480   0,4039   0,5040
Rot         701,0   1,1143   1,2131   0,4730   0,5738
Grenze IR   780,0   0,9364   1,1922   0,5164   0,6157

Die Werte wurden mit der Formel

$\ln \beta =c_{4}\mu ^{4}+c_{3}\mu ^{3}+c_{2}\mu ^{2}+c_{1}\mu +c_{0}\ ;\mu =\ln {\frac {\lambda }{10^{-10}\,\mathrm {m} }}$

berechnet. Dabei ist

c₀ =  6,3521488616·10³
c₁ = -2.7016346146·10³
c₂ =  4.3130312867·10²
c₃ = -3.0621513866·10¹
c₄ =  8.1542004069·10^-1

Wechselwirkende Galaxien

Passage zweier Galaxien auf parabolischen Bahnen um den gemeinsamen Schwerpunkt

Diese z.T. im Artikel Galaxie verlinkten Bilder entstanden mit Hilfe einer reduzierten N-Körper-Simulation, d.h. die Galaxien werden durch Punkmassen ersetzt, die von Testteilchen umkreist werden. Das Verfahren erhebt keinerlei Anspruch auf wissenschaftlichen Wert (heutzutage wird mit hydrodynamischen Codes wie SPH gearbeitet), liefert aber eindrucksvolle Bilder, die zumindest die grundlegenden Effekte wie Gezeitendeformation, Gezeitenarme usw. zeigt. Zudem ist das Verfahren sehr schnell, da die Rechenzeit nur proportional zu N statt zu N² (bei direkter Teilchen-Teilchen-Wechselwirkung) ansteigt.

Der "Bulge" wurde nicht mit einem physikalisch korrekten Modell, sondern mit einem einfachen Trick erzeugt: Die anfangs kreisförmigen Bahnen der zufallsverteilten Teilchen werden abhängig vom Radius um einen Winkel

\Delta \theta ={f \over a^{q}}

gekippt, wobei f und q positive Zahlen sind und a der Betrag des Bahnradius des jeweiligen Testteilchens in Kiloparsec (kpc) ist. Der Inklinationswinkel wird dabei als "Skalenhöhe" für eine exponentielle Abnahme der Teilchendichte mit zunehmendem Abstand (Winkel) zur galaktischen Ebene verwendet, so dass eine weiche Grenze entsteht.

In der Simulation bewegen sich die Galaxien auf parabolischen Bahnen um den gemeinsamen Schwerpunkt, wobei der kleinste Abstand 20 kpc beträgt. Die Anfangspositionen und -geschwindigkeiten der Galaxien wurden zunächst um θ um die y-Achse, dann um φ um die (neue) z-Achse (Bildebene = x-y-Achse) gedreht, um die Neigungen der Ebenen zueinander zu erreichen (mathematischer Drehsinn). Dabei wurde diese Matrix verwendet:

\left({\begin{matrix}\cos \theta \cos \phi &-\sin \phi &\sin \theta \cos \phi \\\cos \theta \sin \phi &\cos \phi &\sin \theta \sin \phi \\-\sin \theta &0&\cos \theta \end{matrix}}\right)

Anschließend wurde das gesamte System noch um die Eulerwinkel α, β, γ gedreht, um den Blickwinkel zu verändern. Dafür wurden die folgenden Drehmatrizen nacheinander ausgeführt:

\left({\begin{matrix}\cos \alpha &-\sin \alpha &0\\\sin \alpha &\cos \alpha &0\\0&0&1\end{matrix}}\right)\circ \left({\begin{matrix}\cos \beta &0&\sin \beta \\0&1&0\\-\sin \beta &0&\cos \beta \end{matrix}}\right)\circ \left({\begin{matrix}\cos \gamma &-\sin \gamma &0\\\sin \gamma &\cos \gamma &0\\0&0&1\end{matrix}}\right)

In diesen Bildern ist &alpha=-160°, β=40° und γ=-30°. Die Zeitpunkte der Bilder relativ zum Zeitpunkt des Perigals t = 0 sind -152,0, -35,6 und +120,2 Mio Jahre.

Daten zu den Bildern
Galaxie	N	a_max	f	q	φ	θ
1 (blau)	15000	15 kpc	0,005	2,8	-30°	+70°
2 (rot)	5000	9 kpc	0,01	4,0	-30°	-85°