Unterschiede

Hier werden die Unterschiede zwischen zwei Versionen gezeigt.

--- ss17:viele_dinge_fliegen_im_weltall_durcheinander:testprogramme [2017/08/23 11:55]
mtoppermann
+++ ss17:viele_dinge_fliegen_im_weltall_durcheinander:testprogramme [2017/08/28 15:41] (aktuell)
mtoppermann [Vergleich von Planetenbewegungen mit verschiedenen Methoden]
@@ Zeile 1: / Zeile 1: @@
 ====== Testprogramme ======
-===== Vergleich =====
+Diese Testprogramme helfen uns dabei gewisse Funktionen auf ihre Richtigkeit zu prüfen, wie auch zum Vergleichen welche Funktion bessere Ergebnisse liefert.
+===== Vergleich von Planetenbewegungen mit verschiedenen Methoden=====
+In diesem Skript wird die Genauigkeit des Euler-Verfahrens mit dem Leapfrog-Verfahren am Bespiel der Erdbewegung
+um die Sonne verglichen. Die Bewegung der dunkelblauen Erde wir durch das Leapfrog-Verfahren berechnet, das der hellblauen durch das Euler-Verfahren. Nach einiger Zeit sieht man bereits deutliche Unterschiede in der Bahnbewegung: Die dunkelblaue Erde verlässt die richtige Bahn langsamer, da ihre die Berechnungsfehler kleiner sind.
+<code python>
+from functions import *
+import pygame
+class Sun(object):
+    def __init__(self):
+        """initialize all necessary properties of the sun
+        """
+        self.x = 0
+        self.y = 0
+        self.radius = 20
+        self.mass = 1.989e30
+        self.colour = (255, 255, 0)
+class Earth(object):
+    def __init__(self, color):
+        """initialize all necessary properties of the earth
+        """
+        self.x = 1.49e11
+        self.y = 0.
+        self.radius = 10
+        self.mass = 5.97e24
+        self.velocity = np.array([0., 2.978e4])
+        self.force = np.array([0, 0])
+        self.colour = (color)
+    def move(self):
+        """leapfrog method for differential-equation-solving
+        """
+        self.x += 0.5 * self.velocity[0] * dt
+        self.y +=  0.5 * self.velocity[1] * dt
+        self.velocity += self.force/self.mass * dt
+        self.x += 0.5 * self.velocity[0] * dt
+        self.y +=  0.5 * self.velocity[1] * dt
+    def move2(self):
+        """euler-method for differential-equation-solving
+        """
+        self.x += self.velocity[0] * dt
+        self.y += self.velocity[1] * dt
+        self.velocity += self.force/self.mass * dt
+def draw(a, points):
+    """Draws an object on the screen
+    """
+    pygame.draw.circle(screen, a.colour, (int(a.x * 10e-10) + width/2,
+    int(a.y * 10e-10) + height/2), a.radius, 0)
+    pygame.draw.aalines(screen, a.colour, False, points, 2)
+def draw2(a):
+    """Draws an object in the middle of the screen
+    """
+    pygame.draw.circle(screen, a.colour, (width/2, height/2), a.radius, 0)
+dt = 36000
+(width, height) = (1000, 600)
+screen = pygame.display.set_mode((width, height))
+pygame.display.set_caption(" ")
+bgcolour = (10, 10, 10)
+earth = Earth((0, 0, 255))
+kerbin = Earth((50, 150, 255))
+sun = Sun()
+running = True
+while running:
+    for event in pygame.event.get():
+        if event.type == pygame.QUIT:
+            running = False
+    # display everything
+    screen.fill(bgcolour)
+    draw(erde, linepoints1)
+    draw(kerbin, linepoints2)
+    draw2(sonne)
+    pygame.display.update()
+    # calculate everything
+    erde.force = get_force(earth, sun)
+    kerbin.force = get_force(kerbin, sun)
+    erde.move()
+    kerbin.move2()
+ </code>
+=====Euler-Verfahren vs Verlet(Leapfrog)-Verfahren=====
+<code python>
+import matplotlib.pyplot as plt
+import numpy as np
+#  x .. = f(x)
+def f(x):
+    return -x
+def leapfrog(startwertx, startwertv, intervallende, anzahlschritte):
+    h = (intervallende - startwertx)/anzahlschritte
+    tpoints = np.linspace(0, intervallende, anzahlschritte)
+    xpoints = []
+    vpoints = []
+    x = startwertx
+    v = startwertv
+    for t in tpoints:
+        xpoints.append(x)
+        vpoints.append(v)
+        x += 0.5 * v * h
+        v += f(x) * h
+        x += 0.5 * v * h
+        #x += h * f(x) + f(x) * h**2/2
+    return tpoints, xpoints
+def euler(startwertx, startwertv, intervallende, anzahlschritte):
+    h = (intervallende - startwertx)/anzahlschritte
+    tpoints = np.linspace(0, intervallende, anzahlschritte)
+    xpoints = []
+    vpoints = []
+    x = startwertx
+    v = startwertv
+    for t in tpoints:
+        xpoints.append(x)
+        vpoints.append(v)
+        v_neu = v + h * f(x)
+        x += h * v
+        v = v_neu
+        #x += h * f(x) + f(x) * h**2/2
+    return tpoints, xpoints
+t1, x1 = leapfrog(0.,1., 100, 500)
+t2, x2 = euler(0.,1., 100, 500)
+pointsin = np.sin(t1)
+abweichung_1 = []
+abweichung_2 = []
+for i in range(len(t1)):
+    abweichung_1.append(pointsin[i]-x1[i])
+    abweichung_2.append(pointsin[i]-x2[i])
+fig = plt.figure(figsize = (18, 8.5))
+ax = fig.add_subplot(111)
+ax.plot(t1, x1, label="Leapfrog")
+ax.plot(t2, x2, label="Euler")
+ax.plot(t1, np.sin(t1), label="np.Sinus")
+ax.plot(t1, abweichung_1, label="Abweichung Leapfrog")
+ax.plot(t1, abweichung_2, label="Abweichung Euler")
+ax.plot(t1, 0*t1, color="black")
+plt.legend(loc=best)
+plt.show()
+</code>
+=====Asteoridengürtel nur mit Sonnengravitation=====
+<code python>
+import numpy as np
+import pygame
+from functions import get_force
+class Asteorid(object):
+	def __init__(self,phi,radius):
+		self.x = np.cos(phi) * radius
+		self.y = np.sin(phi) * radius
+		self.radius = np.random.randint(1,3)
+		self.mass = 6.687e15
+		self.velocity = np.sqrt(G * 1.989e30/radius)
+		* np.array([-np.sin(phi), np.cos(phi)])
+		self.force = np.array([0.,0.])
+		self.colour = (100,100,100)
+	def move(self):
+		self.x += 0.5 * self.velocity[0] * dt
+		self.y +=  0.5 * self.velocity[1] * dt
+		self.velocity += self.force/self.mass * dt
+		self.x += 0.5 * self.velocity[0] * dt
+		self.y +=  0.5 * self.velocity[1] * dt
+class Sun(object):
+    """Our Sun in the point of origin
+    """
+    def __init__(self):
+        self.x = 0
+        self.y = 0
+        self.radius = 20
+        self.mass = 1.989e30
+        self.colour = (255,255,0)
+def draw(a):
+    """Draws an object and its trail on the screen
+    """
+    pygame.draw.circle(screen, a.colour, (int(a.x * 0.5e-9) + width/2,
+    int(a.y * 0.5e-9) + height/2), a.radius, 0)
+def draw2(a):
+    """Draws an object in the middle of the screen
+    """
+    pygame.draw.circle(screen, a.colour, (width/2, height/2), a.radius, 0)
+G = 6.673e-11
+AE = 1.496e11
+n = 2000
+dt = 36000
+liste = []
+(width, height) = (1000,600)
+screen = pygame.display.set_mode((width, height))
+pygame.display.set_caption(" ")
+background_colour = (10,10,10)
+for i in xrange(n):
+	phi = 2 * np.random.random() * np.pi
+	radius = np.random.uniform(2.,3.4) * AE
+	obj = Asteorid(phi, radius)
+	liste.append(obj)
+sonne = Sun()
+running = True
+while running:
+	for event in pygame.event.get():
+		if event.type == pygame.QUIT:
+			running = False
+    # display everything
+	screen.fill(background_colour)
+	draw2(sonne)
+	for i in liste:
+		draw(i)
+		i.force = get_force(i, sonne)
+		i.move()
+	pygame.display.update()
+</code>
+=====Neue Berechnung des Drehimpulses nach einem Zusammenstoß=====
+<code python>
+# -*- coding: utf8 -*-
+from functions import *
+import numpy as np
+n = 2
+liste_teilchen = init_particle_list(n)
+part_1,part_2=liste_teilchen
+part_1.velocity = np.array([1.,0.])
+part_2.velocity = np.array([0.,1.])
+def get_torque(part_1, part_2):
+	"""
+	Drehimpulsberechnung für neu entstandenes Teilchen nach Kollision
+	"""
+	newpoint_x = (part_1.mass * part_1.x + part_2.mass * part_2.x)
+	/ (part_1.mass + part_2.mass)
+	newpoint_y = (part_1.mass * part_1.y + part_2.mass * part_2.y)
+	/ (part_1.mass + part_2.mass)
+	newpoint = ([newpoint_x, newpoint_y]) #Schwerpunkt zwischen beiden Teilchen
+	print part_1.x,part_1.y
+	print part_2.x,part_2.y
+	print newpoint
+	r_1 = ([part_1.x - newpoint_x, part_1.y - newpoint_y]) #Richtungsvektoren
+	r_2 = ([part_2.x - newpoint_x, part_2.y - newpoint_y])
+	print r_1
+	print r_2
+	torq_1 = part_1.mass * np.cross(r_1, part_1.velocity) #Drehimpulse der einzelnen Teilchen zum Schwerpunkt
+	torq_2 = part_2.mass * np.cross(r_2, part_2.velocity)
+	torq = torq_1 + torq_2 #Gesammtdrehimpuls (+ Eigendrehimpulse falls vorhanden)
+	print torq
+	return newpoint,torq
+for i,part_1 in enumerate(liste_teilchen):
+	for j,part_2 in enumerate(liste_teilchen):
+		if i < j:
+			get_torque(part_1, part_2)
+</code>

Mathesis Wiki

Benutzer-Werkzeuge

Webseiten-Werkzeuge

Unterschiede

Seiten-Werkzeuge