Seite 1 von 2

Verfasst: Samstag 3. November 2007, 19:24
von gurke111
Nehmen wir mal diese beiden Funktionen aus dem Skript:

Code: Alles auswählen

def NumerovUpdate(E, PSI):
  n = len(PSI)
  newPsi = (2*(1-5.0/12*H*H*k(n-1,E))*PSI[n-1]-1*(1+H*H/12*k(n-2,E))*PSI[n-2])/(1+H*H/12*k(n,E))
  return newPsi

  
def NumerovInt(drawWin, sgn, E):
  if sgn == 1:
    PSI = [0, 1]
  else:
    PSI = [-H, H]
  for n in xrange(int(N)):
    newPsi = NumerovUpdate(E, PSI)
    PSI.append(newPsi)
    #if abs(newPsi) < PSICUT:
      #if n%2000 == 0:
        #drawWin.Point((n*BREITE/N),newPsi,2)
    if abs(newPsi) > PSIMAX:
      logit ("... abgebrochen (PSIMAX ueberschritten)")
      PSI[:] = [] # leere die Liste
      if newPsi > 0: return 1
      else: return -1
Jetzt komme ich ohne `global` aus. Aber ist es jetzt nicht so, dass in beiden Funktionen je eine Kopie der LANGEN Liste `PSI`existiert? Oder ist das in Python nicht merh so zu trennen... also call by reference, call by value. Call by reference würde hier nur Sinn machen. Macht Python das in dem Fall automatisch?

edit: @ blackjack: oh, jetzt erst deinen letzten post gelesen.. ich gucks mir mal an, moment :)

Verfasst: Samstag 3. November 2007, 19:31
von BlackJack
Es wird keine Kopie gemacht. Es werden *immer* Referenzen auf Objekte übergeben.

Die Unterscheidung "by reference" und "by value" macht nur Sinn wenn es wirklich beide Arten in einer Sprache gibt und wenn es "Referenz" als Datentyp gibt. In einem Java-Buch von Sun steht zum Beispiel Java macht immer "Call by Value" wobei die Werte bei nicht-primitiven Typen eben Referenzen auf Objekte sind.

Verfasst: Samstag 3. November 2007, 19:41
von gurke111
Danke fuer Deine Muehe!

Wegen der Funktionen x, V und k: Naja.. in der physikalischen Problemstellung ist V(x) das Potential. Wenn die Funktion getrennt definiert ist, ist das Skript schneller an ein anderes Problem anzupassen. Genauso ist `k` ein feststehender Ausdruck der "Numerov"-Iteration für die Schrödingergleichung. Bringt Deine Aufteilung einen Geschwindigkeitsvorteil? Wenn nicht, dann nimm es mir nicht übel, wenn ich es so stehen lasse. Trotzdem nochmal danke :)
Es wird keine Kopie gemacht. Es werden *immer* Referenzen auf Objekte übergeben.
Ah, gut zu wissen. Liegt wohl in der Philisophie von Python, ja?
Edit: Aber was ist denn nun, wenn ich eine Variable an eine Funktion übergebe (welches ja, wie ich nun weiß, per Referenz geschieht) und diese dann in der neuen Funktion ändere? Dann sollte ja in der übergebenden Funktion diese Variable nicht mitgeändert werden. Wird in dem Moment dann doch eine Kopie (die dann ja keine echte Kopie mehr ist) angelegt?
Dann habe ich die GUI "ausgebaut" und zwei Testläufe gemacht. Ohne Psyco: 1 Minute und 6 Sekunden, mit Psyco nur 15 Sekunden Laufzeit.
Ich möchte ja eine graphische Ausgabe haben. Auch während der Iteration. Ist einfach anschaulich :) Aber Du hast völlig recht. Die GUI bremst das Skript sehr aus. Ohne einen Punkt zu zeichnen (durch Auskommentieren von `window.Point`) braucht die Iteration bei mir 73s (also Deine Größenordnung ohne Psyco). Wenn ich jeden 2000. Punkt zeichne reicht das für die Anschauung noch aus und die Zeit geht hoch auf 87s. Ich denke, dass das für mich ein Kompromiss ist.

Dieses Psyco schau ich mir bei Gelegenheit gerne mal an. Hört sich ja magisch an. Hab eben mal bisschen gelesen. 10x - 100x Speedup bei rechenintesiven Algorithmen versprechen die, ja? Wow. Problem: Dieses Skript muss auf dem Uni-System laufen. Die haben, denke ich, nur Standardbibliothek. Dann kann ich es gleich vergessen oder?

Verfasst: Samstag 3. November 2007, 20:01
von Leonidas
gurke111 hat geschrieben:Edit: Aber was ist denn nun, wenn ich eine Variable an eine Funktion übergebe (welches ja, wie ich nun weiß, per Referenz geschieht) und diese dann in der neuen Funktion ändere? Dann sollte ja in der übergebenden Funktion diese Variable nicht mitgeändert werden. Wird in dem Moment dann doch eine Kopie (die dann ja keine echte Kopie mehr ist) angelegt?
Es werden nie Kopien implizit angelegt. Du kannst Kopien mit dem Modul `copy` anfertigen.

Es sieht so aus: du hast zwei Typen von Objekten: mutable (änderbare) und immutable (nicht-änderbare). Wenn du einer Funktion ein immutables Objekt übergibst, dann *kann* sie das Objekt nicht verändern. Wenn es ein mutables Objekt ist, kann sie dieses verändern. Beispiel:

Code: Alles auswählen

vier = 3
drei = [1, 2]

def fix_content(drei, vier):
    vier = 4
    drei.append(3)
    drei = [4, 5, 6]

fix_content(drei, vier)
print vier
print drei

Verfasst: Samstag 3. November 2007, 20:10
von gurke111
Okay, Ausgabe:
3
[1, 2, 3]

ich sehe:
`vier` wird nicht verändert, `drei` schon.

`drei` wurde aber auf zwei verschiedene Weisen versucht zu ändern:

Code: Alles auswählen

    drei.append(3)
    drei = [4, 5, 6]
Die erste der beiden Zeilen hat offenbar Wirkung gehabt "nach außen", die zweite nicht. Dann hängt es also auch noch von der Art der Veränderung ab, ob ein Objekt mutable oder immutable ist?

Verfasst: Samstag 3. November 2007, 20:17
von Leonidas
gurke111 hat geschrieben:`drei` wurde aber auf zwei verschiedene Weisen versucht zu ändern:

Code: Alles auswählen

    drei.append(3)
    drei = [4, 5, 6]
Die erste der beiden Zeilen hat offenbar Wirkung gehabt "nach außen", die zweite nicht. Dann hängt es also auch noch von der Art der Veränderung ab, ob ein Objekt mutable oder immutable ist?
Deswegen habe ich auch die zwei Arten aufgeführt. Bei der ersten wird das Objekt verändert (wenn du `id(drei)` davor und danach aufrufst, wirst du sehen dass das das gleiche Objekt ist -> gleiche ID). Bei der zweiten wird an drei ein *neues* Objekt gebunden, d.h. das Objekt was vorher an drei gebunden war, wird gar nicht geändert. Nun ist also im Funktions-Scope drei an ein neues Objekt gebunden. Außerhalb der Funktion hat aber `drei` noch das alte Objekt.

Heißt also: `=` ändert keine Objekte. Nie. Die Objektmethoden ändern Objekte, aber auch nicht immer (bei immutablen Objekten ändern die Objekte sich nach der Erstellung nie, siehe `vier`).

Verfasst: Sonntag 4. November 2007, 00:50
von BlackJack
@gurke111: Meine Aufteilung bzw. Zusammenfassung auf die beiden Funktionen bringt einen kleinen Geschwindigkeitsvorteil, wahrscheinlich hauptsächlich weil Funktionsaufrufe gespart werden. Deine Originalfunktionen ohne GUI haben hier 1 Minute 35 Sekunden ohne Psyco und 17 Sekunden mit Psyco gebraucht.

Die Neuaufteilung war schnell gemacht, ich habe einfach alle Formeln eingesetzt und dann war das neue `f()` ein offensichtlicher Teil, den man als eigene Funktion dort herausziehen konnte. So sieht das "Gesamtkunstwerk" aus:

Code: Alles auswählen

               /           2    /     |     3 |     \   \             /   2    /     |     5 |     \       \
              |       5 * H  * |  E - | n - - | * H  |   |           |   H  * |  E - | n - - | * H  |       |
              |                 \     |     2 |     /    |           |         \     |     2 |     /        |
  2 * PSI   * |  1 - ----------------------------------  | - PSI   * |  ------------------------------ + 1  |
         -1    \                     6                  /       -2    \               6                    /
 -------------------------------------------------------------------------------------------------------------
                                        2    /     |     1 |     \
                                       H  * |  E - | n - - | * H  |
                                             \     |     2 |     /
                                      ------------------------------ + 1
                                                    6
Du kannst `psyco` so benutzen:

Code: Alles auswählen

try:
    import psyco
    psyco.full()
except ImportError:
    pass
Wenn das Modul nicht importiert werden kann, läuft das Programm halt ohne "Turbo". :-)

Verfasst: Sonntag 4. November 2007, 15:22
von gurke111
BlackJack hat geschrieben: Die Einrückung ist jetzt wenigstens konsequent "falsch". ;-)

Da wollte ich nochmal nachfragen.. Sollte ich Deiner Meinung nach lieber 4-spaces Einrückung haben? Ich habe mir jetzt gerade vor kurzem die 2-spaces angewöhnt, weil es bei arg verschachtelten Schleifen/Abfragen sonst schnell richtig richtig abgeht... naja :)

Womit hast Du die Formel so ascii-Bild-mäßig formatiert? Das ist ne coole Sache :>

Verfasst: Sonntag 4. November 2007, 15:32
von Trundle
gurke111 hat geschrieben:
BlackJack hat geschrieben: Die Einrückung ist jetzt wenigstens konsequent "falsch". ;-)

Da wollte ich nochmal nachfragen.. Sollte ich Deiner Meinung nach lieber 4-spaces Einrückung haben?
PEP 8 sagt 4 Spaces.
gurke111 hat geschrieben: Womit hast Du die Formel so ascii-Bild-mäßig formatiert? Das ist ne coole Sache :>
aamath, pretty printer, gibt bestimmt noch massig andere..

Verfasst: Sonntag 4. November 2007, 16:34
von Leonidas
Trundle hat geschrieben:
gurke111 hat geschrieben: Womit hast Du die Formel so ascii-Bild-mäßig formatiert? Das ist ne coole Sache :>
aamath, pretty printer, gibt bestimmt noch massig andere..
Vielleicht hat BlackJack ja auch Maxima benutzt.

Verfasst: Sonntag 4. November 2007, 16:37
von BlackJack
@gurke111: Wenn die Einrückungen zu tief werden, ist das für gewöhnlich ein Zeichen dafür, dass man die Funktion in mehrere Funktionen aufteilen sollte.

Die Formel ist mit Jave gesetzt. Ich hätte auch die Ausgabe von Maxima nehmen können, aber da muss man bei wxMaxima erst auf Textdarstellung umschalten und die gefällt mir nicht so gut wie die von Jave.

Edit: Hab nur mal so aus Neugier noch eine Pyrex-Variante ausprobiert:

Code: Alles auswählen

cdef extern from "math.h":
    double fabs(double x)


cdef double f(double a, double n, double E, double H):
    return (H * H * (E - fabs(n - a) * H)) / 6.0;


def NumerovUpdate(double E, object PSI, double H):
    cdef int n
    cdef double newPsi, p1, p2
    n = len(PSI)
    p1 = PSI[-1]
    p2 = PSI[-2]
    newPsi = ((2 * p1 * (1 - 5 * f(1.5, n, E, H))
               - p2 * (f(2.5, n, E, H) + 1))
              / (f(0.5, n, E, H) + 1))
    PSI.append(newPsi)
    return newPsi
Die liegt mit 19 Sekunden Laufzeit zwischen Psyco und reinem Python. Interessant, dass Psyco hier besser ist.

Verfasst: Sonntag 4. November 2007, 18:46
von birkenfeld
Leonidas hat geschrieben: Vielleicht hat BlackJack ja auch Maxima benutzt.
Oder Maple, oder Mathematica, ... so ziemlich jedes größere CAS hat einen Textkonsolenmodus.

Verfasst: Montag 5. November 2007, 13:50
von BlackJack
Mir ist aufgefallen, das die Liste PSI eigentlich gar nicht benutzt wird, bzw. immer nur die letzten beiden Elemente. Dazu braucht man dann aber keine Liste, sondern nur zwei Namen für die beiden Elemente:

Code: Alles auswählen

def f(a, n, E, H):
    return (H * H * (E - abs(n - a) * H)) / 6.0;


def NumerovUpdate(E, n, psi_1, psi_2, H):
  newPsi = ((2 * psi_1 * (1 - 5 * f(1.5, n, E, H))
             - psi_2 * (f(2.5, n, E, H) + 1))
            / (f(0.5, n, E, H) + 1))
  return newPsi


def NumerovInt(sgn, E):
    if sgn == 1:
        psi_2, psi_1 = 0, 1
    else:
        psi_2, psi_1 = -H, H
    for n in xrange(2, int(N) + 2):
        psi_2, psi_1 = psi_1, NumerovUpdate(E, n, psi_1, psi_2, H)
#         if abs(psi_1) < PSICUT and n % 200 == 0:
#             window.Point((n * BREITE / N), psi_1, 2)
        if abs(psi_1) > PSIMAX:
            logit ("... abgebrochen (PSIMAX ueberschritten)")
            return psi_1 > 0
`psi_1` ist das jeweils letzte berechnete psi und `psi_2` das Vorletzte.

Mit den beiden ersten Funktionen in Pyrex ist die Laufzeit auf 13 Sekunden runter, mit Pyrex + Psyco sogar auf 6 Sekunden! Hier der Pyrex-Code:

Code: Alles auswählen

cdef extern from "math.h":
    double fabs(double x)


cdef double f(double a, int n, double E, double H):
    return (H * H * (E - fabs(n - a) * H)) / 6.0;


def NumerovUpdate(double E, int n, double psi_1, double psi_2, double H):
    newPsi = ((2 * psi_1 * (1 - 5 * f(1.5, n, E, H))
               - psi_2 * (f(2.5, n, E, H) + 1))
              / (f(0.5, n, E, H) + 1))
    return newPsi

Verfasst: Montag 5. November 2007, 23:21
von BlackJack
Hm, mir fiel eben auf, dass `n` nie kleiner als 2 ist, also ist nur eines der drei `abs()` nötig. Das führt dann letztendlich zu folgender "Gesamtfunktion":

Code: Alles auswählen

  //           3        2     \        /             3         2     \     \
- |\|2 n - 5| H  - 2 E H  - 12/ psi  + \(20 n - 30) H  - 20 E H  + 24/ psi |
  \                                2                                      1/
----------------------------------------------------------------------------
                                    3        2
                         (2 n - 1) H  - 2 E H  - 12
Oder in Code:

Code: Alles auswählen

cdef extern from "math.h":
    double fabs(double x)

def NumerovUpdate(double E, int n, double psi_1, double psi_2, double H):
    return (-((fabs(2*n-5)*H*H*H-2*E*H*H-12)*psi_2
              +((20*n-30)*H*H*H-20*E*H*H+24)*psi_1)
            / ((2*n-1)*H*H*H-2*E*H*H-12))
Nochmal fast eine halbe Sekunde gespart. So, nun höre ich aber auf. :-D

Verfasst: Dienstag 6. November 2007, 22:59
von gurke111
Du bist ein Freak, dankeschön :-)

Habe die Sache aber schon am Sonntag "abgegeben" :)
Also kann ich Deine extravaganten Vorschläge nicht mehr verarbeiten. Werd sie aber als Anregung definitiv im Gedächtnis behalten. Und dem Forum nach meinem Debut hier bestimmt auch treu bleiben. Super Leute!