A&A Title Image Startseite | Programme: blauerText | Algorithmen, Lektionen, Analysen | Unsere Vorläufer | Kontakt, Datenschutz

VdS-Journal 85 und 86: Der Dunklen Materie auf der Spur

Schnell-Anwahl:
VdS-Journal 85 – Umlaufgeschwindigkeiten der Sterne - Modellkurve der leuchtenden Materie
VdS-Journal 86 – Umlaufgeschwindigkeiten mit leuchtender und  Dunkler Materie

VdS-Journal 85: Umlaufgeschwindigkeiten in Modellkoordinaten

# rotationGx.py - der dunklen Materie auf die Spur kommen
from turtle import *
from random import *
from math import *
from plot import *



def rho(r): # dichte: Quotient 2 pi wird benötigt, damit die Gesamtmasse 1 ist
    return 1.0 / (2 * pi) * exp(-r)

# Hauptprogramm. Rechenzeit ca. 30 Sekunden (PC im Jahr 2020)
rMax = 8.0             # bis zu welchem Radius?
Schritte = 100         # wieviele Schritte in r-Richtung?
initKoor(0,rMax,0.5, 0,1,0.1) 
d = rMax / Schritte    # Größe eines Quadrates, dessen Masse summiert wird
R = -d / 2.0
for i in range(1,Schritte + 1): # alle Ausgaben in r-Richtung
    R = R + d                   # groß-R ist der Radius, für den berechnet werden soll
    aX = 0                      # die Beschleunigung ist erst mal 0
    aY = 0
    x = -2 * rMax - d / 2       # x auf das äußert links liegende Feld
    for j in range(1,4 * Schritte + 1):
        x = x + d               # x einen Schritt in x-Richtung weiter
        y = -2 * rMax - d / 2   # y auf das äußerst unten liegende Feld
        for k in range(1, 4 * Schritte + 1):
            y = y + d           # einen Schritt in y-Richtung weiter
            r = sqrt(x*x + y*y) # Radius des Summanden-Feldes
            if (r <= 2 * rMax): # ist das innerhalb eines Kreises?
    		    # ja!
                p = rho(r)                   # Dichte lesen
                dx = x - R                   # Abstand vom Zielpunkt ...
                dy = y
                dr = sqrt(dx * dx + dy * dy) # ... ist dr
                a = p / dr / dr      # Beschleunigung nach Gravitationsgesetz
                ax = a * dx / dr     # koordinatenweise aufteilen
                ay = a * dy / dr     # ay kann entfallen, weil insgesamt 0
                aX = aX + ax * d * d # Gesamtbeschleuinigun koordinatenweise summieren
                aY = aY + ay * d * d # aY könnte entfallen: es kompensiert sich
    aX = abs(aX)                # nur der Betrag interessiert (ist nach innen gerichtet)
    w = sqrt(aX / R)            # Kreisfrequenz omega.
    #print(R, w*R, w, aX)
    Plot(R, w * R, 3, "black")  # w*R: omega * R = absolute Geschwindigkeit
update()
done()

VdS-Journal 86: Leuchtende Materie und Dunkle Materie

Das Programm gestatte es, eine Halo aus Dunkler Materie hinzuzufügen. Hierfür sind der Faktor der Masse fM und der Faktor des Radius fR einzusetzen. Eine Halo mit fM=10 und fR=5 hat die 10-fache Masse der leuchtenden Materie, die der Abfall nach dem Exponentialgesetz ist 10-fach langsamer.


# rotationGxDm.py - der dunklen Materie auf die Spur kommen
from turtle import *
from random import *
from math import *
from plot import *

n3198 = [[0.229403, 0.241472],
[0.516157, 0.403917],
[0.74556, 0.482944],
[1.03231, 0.535629],
[1.26172, 0.579533],
[1.49112, 0.619046],
[1.77787, 0.632218],
[2.00728, 0.645389],
[2.23668, 0.649779],
[2.52343, 0.66295],
[2.75284, 0.676122],
[2.98224, 0.680512],
[3.4984, 0.684902],
[3.9572, 0.667341],
[4.47336, 0.667341],
[4.98951, 0.671731],
[5.39097, 0.667341],
[5.96448, 0.65856],
[6.99679, 0.65417],
[7.4556, 0.645389],
[7.97175, 0.649779],
[8.43056, 0.65417],
[8.94672, 0.65417],
[9.46287, 0.656365],
[9.86433, 0.660755],
[10.4378, 0.660755],
[10.8966, 0.65856]]

def rho(r, fM, fR): 
    # Dichte der leuchtenden Materie : Quotient 2 pi wird benötigt, damit die
    # Gesamtmasse 1 ist
    pLM = 1.0 / (2 * pi) * exp(-r)
    # Dichte der dunklen Materie: Quotient 2 pi · fr² wird benötigt, damit die
    # Gesamtmasse 1 ist
    pDM = 1.0 / (2 * pi * fR * fR) * exp(-r / fR)
    # die Dunkle Materie geht mit einem Massenfaktor ein
    return pLM + fM * pDM

# Hauptprogramm.  Rechenzeit ca.  90 Sekunden (PC im Jahr 2020)
rMax = 11.0             # bis zu welchem Radius?
Schritte = 100         # wieviele Schritte in r-Richtung?
fM = 6.5 # Gesamtmasse der Dunklen Materie, bezogen auf die Masse der leuchtenden
         # Materie
fR = 6.9 # wievielfach größer ist der charakteristische Radius der Dunklen Materie?
initKoor(0,rMax,1.0, 0,1,0.1) 
d = rMax / Schritte    # Größe eines Quadrates, dessen Masse summiert wird
R = -d / 2.0
# die Messwerte zeichnen
for i in range(27):
    Plot(n3198[i][0], n3198[i][1], 5, "red")
# die Simulation rechnen und zeichnen
for i in range(1,Schritte + 1): # alle Ausgaben in r-Richtung
    R = R + d                   # groß-R ist der Radius, für den berechnet werden soll
    aX = 0                      # die Beschleunigung ist erst mal 0
    aY = 0
    x = -2 * rMax - d / 2       # x auf das äußert links liegende Feld
    for j in range(1,4 * Schritte + 1):
        x = x + d               # x einen Schritt in x-Richtung weiter
        y = -2 * rMax - d / 2   # y auf das äußerst unten liegende Feld
        for k in range(1, 4 * Schritte + 1):
            y = y + d           # einen Schritt in y-Richtung weiter
            r = sqrt(x * x + y * y) # Radius des Summanden-Feldes
            if (r <= 2 * rMax): # ist das innerhalb eines Kreises?
    		    # ja!
                p = rho(r, fM, fR)           # Dichte lesen
                dx = x - R                   # Abstand vom Zielpunkt ...
                dy = y
                dr = sqrt(dx * dx + dy * dy) # ...  ist dr
                a = p / dr / dr      # Beschleunigung nach Gravitatiionsgesetz
                ax = a * dx / dr     # koordinatenweise aufteilen
                ay = a * dy / dr     # ay kann entfallen, weil insgesamt 0
                aX = aX + ax * d * d # Gesamtbeschleuinigun koordinatenweise summieren
                aY = aY + ay * d * d # aY könnte entfallen: es kompensiert sich
    aX = abs(aX)                # nur der Betrag interessiert (ist nach innen gerichtet)
    w = sqrt(aX / R)            # Kreisfrequenz omega.
    #print(R, w*R, w, aX)
    Plot(R, w * R, 3, "black")  # w*R: omega * R = absolute Geschwindigkeit
update()
done()