Programm für Abschlussarbeit

Wenn du dir nicht sicher bist, in welchem der anderen Foren du die Frage stellen sollst, dann bist du hier im Forum für allgemeine Fragen sicher richtig.
Antworten
marteand
User
Beiträge: 7
Registriert: Montag 8. Juli 2013, 18:29

Hi @ all,

ich schreibe grade eine Abschlussarbeit in der ich ein Programm schreiben muss. Titel der Arbeit: "Effiziente Berechnung der Schallausbreitung.."

Das Programm besteht bis jetzt nur für konkrete Werte und soll noch parametrisiert werden. Außerdem muss es noch für 3d Probleme erweitert werden, bis jetzt nur 2d. Die Ausgabe beim plotten könnte ebenfalls noch verbessert werden. Außerdem muss ich noch das ganze Programm in eine Schleife packen mit passender Abbruchbedingung.

Wäre über jede hilfe wirklich sehr dankbar.

Ich habe noch nie vorher mit Python gearbeitet das meiste aus dem Programm hat mir mein Prof. geschrieben, aber muss jetzt auch mal selbst was hinkriegen.
Falls jemand das Programm in irgend einer Weise verbessern kann würde dass schon helfen.
Danke im vorraus.

Programmtext:

Code: Alles auswählen

from numpy import *                     #importiert alle Dateien von numerischem Python
import matplotlib.pyplot as plt         #importiert Dateien zum ploten
from scipy.integrate import odeint      
from scipy.optimize import *            #importiert Funktionen von Scipy zur 
from mayavi import mlab                 #zur 3D-Darstellung
from time import time                   #zur Ermittlung der Programmlaufzeit

c=343                                   #globale Konstante c

n = 15                                   # Exponent Geschwindigkeitsfunktion;globale Konstante Dimensionslos

v0 = 50.                                # Geschwindigkeit Stroemung in m/s



def v(xx):                              # definition des Geschwindigkeitsprofils
    x,y,z = xx                          #
    v = array((v0*(1-tanh(y**n)), 0., 0.))  #Geschwindigkeitsprofil
    dv = array( ((0,0,0),(-n*v0*y**(n-1)*(1-tanh(y**n)**2),0,0),(0,0,0)) )      #Ableitung des Geschwindigkeitsprofils
    return v,dv                         #

def f(y,t):                             
    x=y[0:3]                            
    s=y[3:6]
    vv,dv = v(x)
    sa = sqrt((s*s).sum())              # Betrag von s-Vektor
    f0 = c*s/sa + vv
    f1 = dot(s,dv)
    return hstack((f0,f1))

t  = linspace(0, 0.008, 2000)

x0 = array((0,0,0.))
n0 = array((1,   0.1,  0.1))/sqrt(2)
s0 = n0 / (c+dot(v(x0)[0],n0))
y0 = hstack((x0,s0))
#ti = time()
soln = odeint(f, y0, t)
ort = soln[:,0:3]

def abstd(n):
        n0 = n/sqrt((n*n).sum())
        s0 = n0 / (c+dot(v(x0)[0],n0))
        y0 = hstack((x0,s0))
        soln = odeint(f, y0, t)
        ort = soln[:,0:3]
        abstand = sqrt((abs(ziel-ort)**2).sum(1))
        indmin = abstand.argmin()
        tmin = t[indmin]
        a = ort[indmin]
        u = ort[indmin+1] - a
        d = cross(ziel-a,u)
        d = sqrt((d*d).sum())/sqrt((u*u).sum())
        return abstand[indmin],tmin#,d



def abstd2(n):
        return abstd(n)[0]
#ziel = array(1,1,0)
ziel  = array((1,1,0.01))
ti = time()
res = fmin(abstd2, n0,maxfun=300,full_output=0)
print time()-ti
print res
print abstd(res)
 
print sqrt(((x0-ziel)**2).sum())/c



plt.plot(t,ort[:,0],label='x')
plt.plot(t,ort[:,1],label='y')
plt.plot(t,ort[:,2],label='z')
plt.legend()
plt.show()
pl=mlab.plot3d(ort[:,0],ort[:,1],ort[:,2],color=(1,0,0))
mlab.axes(pl)
mlab.show()


for y in arange(-1,1,0.1):
 print v( array((0,y,0)) )
y = arange(-2,2,0.01)
vx = 1-tanh(y**10)
vxs = -10*y**9*(1-tanh(y**10)**2) 
plt.plot (y,vxs)
plt.plot (y,vx)
plt.show()


for y in arange(-1,1,0.1):
 print v( array((0,y,0)) )
y = arange(-2,2,0.01)
vx = 1-tanh((y**10)
vxs = -10*y**9*(1-tanh(y**10)**2) 
plt.plot (y,vxs)
plt.show()
Zuletzt geändert von Anonymous am Montag 8. Juli 2013, 18:47, insgesamt 1-mal geändert.
Grund: Quelltext in Python-Code-Tags gesetzt.
Benutzeravatar
darktrym
User
Beiträge: 784
Registriert: Freitag 24. April 2009, 09:26

Das muss Hass sein, wieso stimmen nicht mal die Einrückungen?
„gcc finds bugs in Linux, NetBSD finds bugs in gcc.“[Michael Dexter, Systems 2008]
Bitbucket, Github
marteand
User
Beiträge: 7
Registriert: Montag 8. Juli 2013, 18:29

Wieso Hass??

Prog ist ja noch nicht fertig ... ich weiß die Form ist nicht soo schön.

Und Hilfe oder eine idee?
BlackJack

@marteand: Erster Verbesserungsvorschlag: Klammern so setzen das der Compiler nicht mit einem Syntax-Fehler abbricht. Dann läuft das vielleicht besser. ;-)

Dann könnte man an der Form ein wenig arbeiten. Konsistente Einrückung, ein wenig mehr Leerzeichen, damit es lesbarer wird. Also allgemein den Style Guide for Python Code berücksichtigen.

Der Mix zwischen Berechnungen und Definitionen von Funktionen ist auch etwas unübersichtlich. Normalerweise definiert man erst alle Funktionen und das „Hauptprogramm” kommt auch in eine Funktion, damit man nicht versehentlich (oder absichtlich) in den Funktionen auf etwas zurückgreift, was nicht als Argument übergeben wurde aber keine Konstante ist.

Bei Sternchenimporten besteht die Gefahr von Namenskollisionen. Ausserdem weiss man nicht mehr so einfach aus welchem Modul ein Name stammt ohne dass man das auswendig weiss, oder nachschlägt. Solange bis man es gefunden hat.

Kommentare sollten dem Leser einen Mehrwert zum Quelltext geben. Bei ``#globale Konstante c`` zu ``c=343`` ist der Leser kein Stück schlauer. Dass das eine globale Konstante ist und das sie `c` heisst, sieht man ganz prima an der Zuweisung und das sie auf Modulebene steht. Viel interessanter ist doch was diese Konstante im Kontext des Programms *bedeutet*. Faustregel: In Kommentaren nicht schreiben *was* der Code macht, das kann man am Code selbst schon ablesen, sondern *warum* er das macht was er tut. Bei Konstanten/Variablen wäre ein Name der nicht nur aus einem Buchstaben besteht, sondern dem Leser vermittelt was der Wert bedeutet, sogar noch besser als ein Kommentar. Das gilt auch für Funktionsnamen.

Kommentare werden üblicherweise auch vor den Code geschreiben, den sie kommentieren. Am Zeilenende hat man das Problem, dass man oft nicht mehr so einfach die 80 Zeichen-Grenze einhalten kann. Ausserdem hast Du die alle ausgerichtet, was nervig wird wenn man die Zeilen verändert. Dann muss man immer die Ausrichtung wieder reparieren. Funktionen kann man mit einem DocString dokumentieren statt mit Kommentaren.
marteand
User
Beiträge: 7
Registriert: Montag 8. Juli 2013, 18:29

Ich habe jetzt mal versucht das ganze übersichtilicher und in besserer Form zu gestalten. Wie ich das mit den Kommentaren vorne richtig realisieren soll ist mir nicht klar. Sry dass das letzte einen syntax fehler hatte. Hab da wohl irgendeine Version genommen. Hier nochmal eine bessere Version:

from numpy import* #importiert alle Dateien vom Modul numpy; ermoeglicht Vektor- und Matrizenrechnung mit arrays

from scipy.integrate import odeint #importiert aus dem Modul Scipy die Funktion odeint zum loesen der DGL's
from scipy.optimize import fmin #importiert Funktionen von Scipy zur Minimierung des Abstands zwischen Mikrofon und Strahl

import matplotlib.pyplot as plt #importiert Dateien zum ploten
from mayavi import mlab #zur 3D-Darstellung der Loesung

from time import time #zur Ermittlung der Programmlaufzeit

c=343 #globale Konstante c Schallgeschwindigkeit [m/s]; die Schallgeschwindigkeit wird wie die Lufttemperatur und Luftfeuchte als konstant betrachtet
v0 = 50. #globale Konstante Stroemungsgeschwindigkeit [m/s]

n = 10 # Exponent der Geschwindigkeitsfunktion;globale Konstante Dimensionslos;es muss mindestens n = 2 sein bei groesser werden dem n wird das geschwindigkeitsprofil rechteckfoermiger
ko= 1 # Koeffizient der Geschwindigkeitsfunktion; globale Konstante, Dimensionslos;veraendert den Durchmesser des Geschwindigkeitsprofils
zx = 2
zy = 2
zz = 2
#n0 = array((nx, ny, nz))/sqrt(nx**2+ny**2+nz**2)
#x0 = array((sx, sy, sz))
ziel = array((zx, zy, zz))

n0 = array((1, 1, 1))/sqrt(3)
x0 = array((0, 0, 0.))
#ziel = array((1, 1, 1))
#s0 = n0 / (c+dot(v(x0)[0],n0))
#y0 = hstack((x0,s0))

t = linspace(0, 0.08, 2000)

def v(xx):
'Funktion zur Beschreibung Geschwindigkeitsprofils'
x,y,z = xx
v = array((v0*(1-tanh(ko*((y**2+z**2)**n))), 0., 0.)) #Geschwindigkeitsprofil
dv = array( ((0,0,0),((-n*v0*y**(n-1)*(1-tanh(ko*((y**2+z**2)**n)))**2),0,0),((-n*v0*z**(n-1)*(1-tanh(ko*((y**2+z**2)**n)))**2),0,0)) ) #Ableitung des Geschwindigkeitsprofils
return v,dv #

def f(y, t):
'Funktion zur Beschreibung der DGLs'
x = y[0:3]
s = y[3:6]
vv, dv = v(x)
sa = sqrt((s*s).sum()) # Betrag von s-Vektor
f0 = c*s/sa + vv
f1 = dot(s, dv)
return hstack((f0, f1))

def abstd(n):
'Funktion zur Minimierung des Abstands'
n0 = n/sqrt((n*n).sum())
s0 = n0 / (c+dot(v(x0)[0], n0))
y0 = hstack((x0, s0))

soln = odeint(f, y0, t)
ort = soln[:,0:3]
abstand = sqrt((abs(ziel-ort)**2).sum(1))
indmin = abstand.argmin()
tmin = t[indmin]
a = ort[indmin]
u = ort[indmin+1] - a
d = cross(ziel-a,u)
d = sqrt((d*d).sum())/sqrt((u*u).sum())
return abstand[indmin], tmin, d

def abstd2(n):
return abstd(n)[0]

ti = time()
#while d>0.001: ;versuch Hauptprogramm in schleife zu packen'

def prog():
'Funktion für Hauptprogramm'

ti = time()

s0 = n0 / (c+dot(v(x0)[0],n0))
y0 = hstack((x0,s0))
soln = odeint(f, y0, t)
ort = soln[:,0:3]
return prog

res = fmin(abstd2, n0,maxfun=300,full_output=0)

print time()-ti
print res
print abstd(res)
print sqrt(((x0-ziel)**2).sum())/c

s0 = n0 / (c+dot(v(x0)[0],n0))
y0 = hstack((x0,s0))
soln = odeint(f, y0, t)
ort = soln[:,0:3]
plt.plot(t,ort[:,0],label='x')
plt.plot(t,ort[:,1],label='y')
plt.plot(t,ort[:,2],label='z')
plt.legend()
plt.show()
pl=mlab.plot3d(ort[:,0],ort[:,1],ort[:,2],color=(1,0,0))
mlab.axes(pl)
mlab.show()


#for y in arange(-1,1,0.8):
#print v( array((0,y,0)) )
#y = arange(-2,2,0.01)
#vx = 1-tanh(ko*((y**2+z**2)**n
#vx = array((1-tanh(ko*(y**n)), 0., 1-tanh(ko*(y**n)))
#vxy = (-n*v0*y**(n-1)*(1-tanh(ko*((y**2+z**2)**n)))**2)
#plt.plot (y,vxy)
#plt.plot (y,vx)
#plt.show()
#pl=mlab.plot3d(vx[:,0],vx[:,1],vx[:,2],color=(0,1,0))
#mlab.axes(pl)
#mlab.show()
Sirius3
User
Beiträge: 17753
Registriert: Sonntag 21. Oktober 2012, 17:20

@marteand: erst einmal, für Quelltexte gibt es oberhalb des Editfeldes den Knopf »python«. Dann kann man das auch lesen, da ohne Einrückungen Pythoncode völlig sinnlos ist.
Welchen von BlackJacks Vorschlägen hast Du jetzt versucht umzusetzen?
marteand
User
Beiträge: 7
Registriert: Montag 8. Juli 2013, 18:29

Code: Alles auswählen

from numpy import*                     #importiert alle Dateien vom Modul numpy; ermoeglicht Vektor- und Matrizenrechnung mit arrays

from scipy.integrate import odeint     #importiert aus dem Modul Scipy die Funktion odeint zum loesen der DGL's 
from scipy.optimize import fmin        #importiert Funktionen von Scipy zur Minimierung des Abstands zwischen Mikrofon und Strahl

import matplotlib.pyplot as plt        #importiert Dateien zum ploten
from mayavi import mlab                #zur 3D-Darstellung der Loesung

from time import time                  #zur Ermittlung der Programmlaufzeit

c=343                                  #globale Konstante c Schallgeschwindigkeit [m/s]; die Schallgeschwindigkeit wird wie die Lufttemperatur und Luftfeuchte als konstant betrachtet
v0 = 50.                               #globale Konstante Stroemungsgeschwindigkeit [m/s]

n = 10                                 # Exponent der Geschwindigkeitsfunktion;globale Konstante Dimensionslos;es muss mindestens n = 2 sein bei groesser werden dem n wird das geschwindigkeitsprofil rechteckfoermiger
ko= 1                                  # Koeffizient der Geschwindigkeitsfunktion; globale Konstante, Dimensionslos;veraendert den Durchmesser des Geschwindigkeitsprofils
zx = 2
zy = 2
zz = 2
#n0 = array((nx, ny, nz))/sqrt(nx**2+ny**2+nz**2)
#x0 = array((sx, sy, sz))
ziel  = array((zx, zy, zz))             

n0 = array((1, 1, 1))/sqrt(3)
x0 = array((0, 0, 0.))
#ziel  = array((1, 1, 1))
#s0 = n0 / (c+dot(v(x0)[0],n0))
#y0 = hstack((x0,s0))

t  = linspace(0, 0.08, 2000)

def v(xx): 
    'Funktion zur Beschreibung Geschwindigkeitsprofils'
    x,y,z = xx                          
    v = array((v0*(1-tanh(ko*((y**2+z**2)**n))), 0., 0.))  #Geschwindigkeitsprofil
    dv = array( ((0,0,0),((-n*v0*y**(n-1)*(1-tanh(ko*((y**2+z**2)**n)))**2),0,0),((-n*v0*z**(n-1)*(1-tanh(ko*((y**2+z**2)**n)))**2),0,0)) )      #Ableitung des Geschwindigkeitsprofils
    return v,dv                         #

def f(y, t):
    'Funktion zur Beschreibung der DGLs'
    x = y[0:3]                            
    s = y[3:6]
    vv, dv = v(x)
    sa = sqrt((s*s).sum())              # Betrag von s-Vektor
    f0 = c*s/sa + vv
    f1 = dot(s, dv)
    return hstack((f0, f1))

def abstd(n):
    'Funktion zur Minimierung des Abstands'
    n0 = n/sqrt((n*n).sum())
    s0 = n0 / (c+dot(v(x0)[0], n0))
    y0 = hstack((x0, s0))
        
    soln = odeint(f, y0, t)
    ort = soln[:,0:3]
    abstand = sqrt((abs(ziel-ort)**2).sum(1))
    indmin = abstand.argmin()
    tmin = t[indmin]
    a = ort[indmin]
    u = ort[indmin+1] - a
    d = cross(ziel-a,u)
    d = sqrt((d*d).sum())/sqrt((u*u).sum())
    return abstand[indmin], tmin, d
    
def abstd2(n):
    return abstd(n)[0]

ti = time()
#while d>0.001: ;versuch Hauptprogramm in schleife zu packen'
    
def prog():
    'Funktion fuer Hauptprogramm'
    
    ti = time()
    
    s0 = n0 / (c+dot(v(x0)[0],n0))
    y0 = hstack((x0,s0))
    soln = odeint(f, y0, t)
    ort = soln[:,0:3]
    return prog
    
res = fmin(abstd2, n0,maxfun=300,full_output=0)
    
print time()-ti
print res
print abstd(res)
print sqrt(((x0-ziel)**2).sum())/c

s0 = n0 / (c+dot(v(x0)[0],n0))
y0 = hstack((x0,s0))
soln = odeint(f, y0, t)
ort = soln[:,0:3]
plt.plot(t,ort[:,0],label='x')
plt.plot(t,ort[:,1],label='y')
plt.plot(t,ort[:,2],label='z')
plt.legend()
plt.show()
pl=mlab.plot3d(ort[:,0],ort[:,1],ort[:,2],color=(1,0,0))
mlab.axes(pl)
mlab.show()


#for y in arange(-1,1,0.8):
 #print v( array((0,y,0)) )
#y = arange(-2,2,0.01)
#vx = 1-tanh(ko*((y**2+z**2)**n
#vx = array((1-tanh(ko*(y**n)), 0., 1-tanh(ko*(y**n)))
#vxy = (-n*v0*y**(n-1)*(1-tanh(ko*((y**2+z**2)**n)))**2) 
#plt.plot (y,vxy)
#plt.plot (y,vx)
#plt.show()
#pl=mlab.plot3d(vx[:,0],vx[:,1],vx[:,2],color=(0,1,0))
#mlab.axes(pl)
#mlab.show()
Zuletzt geändert von marteand am Dienstag 9. Juli 2013, 22:00, insgesamt 1-mal geändert.
marteand
User
Beiträge: 7
Registriert: Montag 8. Juli 2013, 18:29

Ich hab versucht mehr Leerzeichen einzubauen sowie die Funktionen mit Docstrings zu erläutern. Ausserdem habe ich versucht das ganze in eine Hauptfunktion zu schreiben wobei mir nicht klar ist ob das so richtig und sinnvoll war. Das es eine Hauptfunktion geben muss ist mir bewusst hab da mal versucht eine einzubauen. Diese Hauptfunktion sollte auch noch in eine Schleife wobei ich für diese noch keine sinnvolle Abbruchbedingung habe. Außerdem habe ich versucht mehr globale Variablen zu deklarieren. Sorry wenn vielleicht manche Begriffe nicht ganz richtig sind. Ich arbeite das erste mal mit Python und meine Informatik bzw. Programmier Kenntnisse sind nicht so dolle...

Das ganze sollte wenns geht auch 3d sein, habe das auch versucht umzusetzen allerdings ist es mir nicht gelungen die Geschwindigkeitsfunktion in 3d zu plotten. Den Strahlenverlauf der Schallausbreitung konnte ich in 3d darstellen.
Benutzeravatar
snafu
User
Beiträge: 6740
Registriert: Donnerstag 21. Februar 2008, 17:31
Wohnort: Gelsenkirchen

@marteand: Du wirst dir zur Lösung dieser Aufgabe wohl (Online-)Literatur schnappen und dich gründlich in Python, sowie NumPy und SciPy einarbeiten müssen (den mathematischen Background bringst du ja hoffentlich mit). Hattet ihr im ganzen Studium denn überhaupt nichts mit diesen Themen zu tun? Und nur aus Interesse: Bis wann muss der Kram denn fertig werden?
marteand
User
Beiträge: 7
Registriert: Montag 8. Juli 2013, 18:29

Das ganze hätte natürlich schon gestern fertig werden sollen.
Ich hatte einen Info Kurs im studium aber das was objektorientierte Programmierung ansonsten hatte man auch 1-2 mal matlab gebraucht aber war nur nebenbei in irgendwelchen kursen,mehr war da nicht. Ansonsten hatten wir vor allem den umgang mit cad und fem Programmen... Aber nix selber proggen...

Mit der online Literatur habe ich schon versucht. Und manchmal findet mann auch gute Erklärungen von anderen Unis oder eben python.docs...
ansonsten habe ich auch viel mit der Hilfe in Python selbst versucht.
Die Mathematik verstehe ich schon so halbwegs. Aber dem programm das beizubringen bzw. verstehen wie die das meinen .... ist für mich wie chinesisch.
Das bisher entstandene war auch mehr try and error und viel hilfe vom Prof.
Antworten