Schema gesucht: Sehr teure Gleichungen parallel lösen, ohne viel Rechenzeit zu verschwenden

SeppJ

Sehr viel dazugekommen. eine lange Antwort:

rapso schrieb:

deine beschreibung ist auch ein wenig verwirrend. da du ueber stetigkeit/differenzierbarkeit sprichst, suchst du naehrungsverfahren die besser sind als z.B. newton rhapson?

Die Stetigkeit und Differenzierbarkeit ist nur Nebeninformation, falls es etwas nützt. Ich kann leider nicht sehr viel über die Funktion aussagen.

da du von parallelisierung der funktion sprichst, suchst du eine einfach loesung um sie mehrmals parralel fuer verschiedene aufgabenstellungen abzuarbeiten?

Das kann ich schon. Ich kann viele Funktionswerte gleichzeitig bestimmen, aber die Berechnung des jeweiligen Einzelwertes ist nicht weiter verbesserbar.

zudem sagst du, dass du das verfahren fuer einen wert nicht parallelisieren kannst, dann schreibst du dass du keine 10k schritte machen willst (und zudem ist die funktion stetig), somit koenntest du schon mehrere punkte gleichzeitig fuer ein problem berechnen und dadurch die naehrung verbessern? (z.b. newton-rhapson mit einer spline statt nur der funktions-tangente, wenn du z.b. immer 4 samples berechnen wuerdest).

Ja, das ist möglich und etwas, woran ich noch nicht gedacht hatte. (Erstaunlicherweise muss ich sagen. Jetzt wo du es sagst ist es dumm, dass ich daran nicht gedacht habe)

suchst du komplett neue loesungsansaetze oder optimierungen fuer was du schon hast?

Jain. Dieses Problem ist neu für mich. Ich habe noch nichts funktionierendes, da ich bei allen meinen bisherigen Ansätzen zu dem Schluss gekommen bin, dass sie zu lange dauern würden.

TGGC schrieb:

Hast du mal darueber nachgedacht E(x) auf einer GPU zu implementieren?

Kommt drauf an: Zur Berechnung eines Wertes würde es nichts bringen, es würde sogar wesentlich langsamer, da der einzelne GPU-Kern nicht mit einer CPU mithalten kann. Es könnte was bringen, um parallel mehrere Werte zu berechnen. Aber der Portierungsaufwand wäre enorm und ich kann bereits dutzende Einzelwerte parallel berechnen (aber nicht hunderte, wie bei einer GPU). Und die vorhandene Hardware ist schon auf ein klassisches CPU-lastiges Programm ausgelegt. Aus diesen Gründen habe ich die Idee verworfen.

knivil schrieb:

Wie waere es mit einer Taylorentwicklung ...

Taylorentwicklung ist etwas für theoretische Mathematik, nicht für angewandte Numerik. Oder möchtest du mir etwas konkretes mit dem Stichwort sagen, dass ich gerade nicht verstehe?

Bitte ein Bit schrieb:

Als Nebenbedingung soll dabei gelten $F(Exi,Exi+1)=C$ wobei C eine vorgegebene Konstante ist.

Diese Bedingung ist irgentwie ein Knackpunkt.

Wie wäre es wenn du E(x) durch eine einfache Funktion E'(x) approximierst und danach erst die Funktion F(x, y) mit E'(x) bestimmst ?

Du musst du wohl erst einige E(xi)'s berechen und mit denen eine passende Approximation für E(x) suchen. Dabei würde ich versuchen den Grundtyp von E'(x) festzustellen. Nicht das E(x) vom Typ her eine Schwingung ist und du als Ansatz für E'(x) ein Polynom nimmst. Ist vermutlich eine Optimierungsaufgabe über Approximationen, abhängig von dem simulierten System.

Hmm, das ist eine Idee. Das Problem ist, dass ich derzeit nicht viel über die Funktion aussagen kann, erst recht nicht über einen möglichen Grundtyp. Und es ist auch nicht nur eine Funktion, sondern ich werde das Problem oftmals mit vielen ähnlichen Funktionen lösen müssen (darum frag ich ja auch nach einem möglichst allgemeinen Algorithmus, sonst könnte ich es auch einmal mit viel Geduld von Hand machen). Aber ich könnte mal über das Wochenende für ein paar der Funktionen mal ein paar Werte über den ganzen Bereich verteilt angucken. Vielleicht gibt es da ja schöne Muster.

knivil

Oder möchtest du mir etwas konkretes mit dem Stichwort sagen

Nimm als Approximation die Taylorentwicklung.

otze

mal ein anderer Ansatz, vielleicht hilft der dir:

da du ja nicht nur einen Wert, sondern eine Kette (x_0,...,x_n) brauchst, kannst du eigentlich nach all diesen Punkten parallel suchen

das heißt du wirfst einen Optimierungsalgorithmus auf der Funktion

$\sum\limits_{i=1}^n (F(E(x_{i-1}),E(x_i))-C)^2$
an. An der Stelle könntest du dann klassische multidimensionale Optimierungsverfahren wählen. Eventuell brauchst du noch einen Regularisierungsterm der dir versichert, dass die x_i monoton steigend sind oder musst das Problem passend umformulieren (zum Beispiel als $x\_i=x\_{i-1}+y\_i, y\_i > 0$ ).

Wenn du nach vielen x_i suchst, könntest du dann einfach einen effizienten Suchalgorithmus wie die CMA anwerfen. die CMA hätte den Vorteil, dass du alle Evaluationen eines Optimierungsschritts parallelisieren kannst. Andererseits brauchst du nur ca $\sqrt{n}$ parallele Evaluationen.

SeppJ

@otze: Das klingt schon einmal cool, daran hatte ich noch nicht gedacht. Wird eine Weile dauern bis ich mir die Stichworte angelesen habe. Aber wenn deine Beschreibung stimmt, dann geht dein Vorschlag mein Hauptproblem an, was natürlich recht günstig ist.

cooky451

SeppJ schrieb:

aber das was du wissen möchtest:

Mich würde tatsächlich auch die Funktion interessieren.

SeppJ

Achtung: Offtopic!

cooky451 schrieb:

SeppJ schrieb:

aber das was du wissen möchtest:

Mich würde tatsächlich auch die Funktion interessieren.

Naja, das ist dann so etwas wie:
E(x) = (<H(x)>, $\sigma_H(x)$ )
Ich kann mit Latex leider keine spitzen Klammern schreiben. Das soll heißen: E(x) ist das Paar aus Erwartungswert von H(x) und Varianz von H(x), wobei der Erwartungswert von H(x) definiert ist als:
<H> = $-\frac{\partial}{\partial \beta} \ln Z(\beta)$ bei vorgegebener Konstante $\beta$ (Einheit: Inverse Energie). Das Z ist wiederum
Z(\beta) = \frac{1}{\text{Vorfaktor}} \int e^{-\beta E(\vec r\_1, \vec r\_2, \dots, \vec r\_N, V, x)}\text{d}^3 \vec r\_1 \text{d}^3 \vec r\_2 \dots \text{d}^3 \vec r\_N \text{d} V
Das E ist jetzt ein anderes E als oben. Ich will mich aber an die gängige Konvention halten. Daher heißt das hier jetzt auch E, denn es ist eine Energie. Das E ist:
E(\vec r\_1, \vec r\_2, \dots, \vec r\_N, V, x) = \sum\_{i,j} \text{Potential1_{i,j}}(\vec r\_i, \vec r\_j) + \sum_{i,j} \text{Potential2_{i,j}}(\vec r\_i, \vec r\_j) + \sum_{i,j,k}\text{Potential3_{i,j,k}}(\vec r\_i, \vec r\_j, \vec{i,j,k}) + E_{\text{V}}(V) + E_{\text{x}}(x, \vec r\_1, \vec r\_2, \dots, \vec{r_N})
Das Z hat nun aber 3N+1 Freiheitsgrade, das Integral kann man niemals ausrechnen. Daher macht man einen Trick: Die meisten Zahlen in dem Exponenten sind sehr groß, tragen also kaum zum Integral bei. wenn man nun herausfinden könnte, welche E die kleinsten sind, dann kann man mit diesen das Z ziemlich gut abschätzen. Außerdem ist hier mein <H> (bis auf Vorfaktor, der dann passend gewählt ist) identisch mit dem Erwartungswert der E (das zweite E), jeweils gewichtet mit dem $e^{-\beta E}$ . wenn ich nun also Konfigurationen der $\vec r_i$ und $V$ erzeugen kann, die mir eine Reihe von E's (das zweite E) liefern, die schon nach diesem Faktor gewichtet sind, dann kann ich einfach den Mittelwert und Varianz über diese E's (das zweite E) nehmen um mein erstes E zu erhalten. Das kann man aber recht einfach erreichen, wenn man eine gegebene Konfiguration (\vec r\_1, \dots ,\vec r\_N, V) leicht verändert, jeweils das E (das zweite) berechnet (nennen wir die Werte $E_1$ und $E_2$ ) und die neue Konfiguration mit einer Wahrscheinlichkeit von $P=\text{max}(1, e^{-\beta(E\_2 -E\_1)})$ akzeptiert. Wenn man so eine Kette von E_i erzeugt, dann kann man beweisen, dass deren Verteilung (bei unendlich vielen Schritten) genau den gewünschten Kriterien entspricht. Aber mit ziemlich vielen Schritten erhält man auch schon ganz gute Werte :p .

Und das war klassische statistische Physik im kanonischen Ensemble für Einsteiger mit einem kleinen Exkurs in Monte Carlo Simulation. Vermutlich hat niemand, der es nicht schon kennt, auch nur ein Wort verstanden . Und die, die es schon kennen, kreiden mir jetzt die ganzen Vereinfachungen und Ungenauigkeiten an, die ich gemacht habe, um einigermaßen verständlich zu bleiben.
**
Ich will es daher nochmal in Worte fassen:** Die Funktion E, um die es in diesem Thread geht, ist der Erwartungswert und die Schwankung einer Messung der Energie eines klassischen Systems von sehr vielen Teilchen mit gewissen Wechselwirkungen. Das kann im Allgemeinen niemand exakt berechnen, da man hundertausende Freiheitsgerade hat, aber die Computerphysik kennt sehr gute Näherungslösungen, die aber auch eine ganze Weile zum Rechnen benötigen.

otze

Sowas hatte ich mir schon gedacht. Ja, da lässt sich leider nicht viel dran drehen.

Vielleicht hilft dir das:
http://www.machinelearning.org/archive/icml2009/papers/229.pdf

es wird zwar eine andere Problemkategorie betrachtet, aber der Algorithmus ist speziell für verrauschte Evaluationsergebnisse(wie du sie hast). Eventuell kannst du damit sogar die Simulationszeiten der Funktion E runter drehen, da sich der Algorithmus dann selbst stabilisiert.

Bitte ein Bit

knivil schrieb:

Nimm als Approximation die Taylorentwicklung.

Eine Taylorreihe als Approximation ist an sich eine gute Sache da man schön den Restfehler abschätzen kann.

Aber wir befinden uns hier im Bereich der Simulation wo die Berechnung eines Wertes f(x) einen Tag benötigt. Und ich fürchte egal welche Funktion das auch sein mag (siehe SeppJ letzten Beitrag), wird das Berechnen der Ableitungen nicht trivialer Natur (im akademischen Sinne) sein.

Hmm, ohne jetzt das Problem wirklich bis ins Detail verstanden zu haben... Waere es moeglich sich zu einem Satz von gut verteilten Parametern einen Satz Ergebnisse zu berechnen und dann zu versuchen eine Loesung in der Mannigfaltigkeit, in der die Ergebnisse liegen, zu approximieren? Quasi so eine Art Reduced-Basis-Ansatz, wie man den z.B. aus dem Finite-Elemente-Kontext kennt.

SeppJ

Nochmal Rückmeldung: Vielen Dank für die Anregungen, mir ist da tatsächlich eine Idee gekommen, die mir sehr gefällt. Das kann man tatsächlich als ein Optimierungsproblem mit wechselnder Dimension und ein paar einfachen Nebenbedingungen ansehen, wodurch schon einmal die Parallelität der Rechnung gegeben wäre. Aber CMA-ES ist da vermutlich eine zu dicke Kanone, die ganzen Approximierer in diesem Thread haben mich ein bisschen genauer über die Funktion nachdenken lassen und ich denke sie ist gutartig genug für eine Interpolation. Mein Algorithmus sähe dann so aus, mit einer Stützpunktverteilung zu starten, mit Nebenbedingung, dass der linkester Stützpunkt mein Ausgangspunkt ist und alle Punkte geordnet sind. Dann wird so gut wie möglich interpoliert und die Stützpunkte nach dieser Interpolation auseinander oder zusammen geschoben, wobei bei Bedarf rechts Punkte in oder aus dem Definitionsintervall wandern, so dass gerade ein Wert jenseits meiner angestrebten rechten Grenze liegt. Wenn die Funktionswerte aller Stützpunktpaare eine gewisse Mindestgüte erreicht hat, bin ich fertig. Verworfene Stützpunkte verbessern immer noch die Interpolation und waren somit nicht völlig vergebens. Ich schätze, das Schema wird nach wenigen (5?) Durchgängen eine ausreichende Genauigkeit erreicht haben, wodurch auch die Gesamtzeit hinreichend klein wäre.

Wird natürlich etwas dauern, bis ich bestätigen kann, ob das funktioniert oder nicht.