DaDaPlayground

Kommentare anzeigen Worklog abonnieren
Gehe zu Seite Zurück  1, 2, 3

Worklogs DaDaPlayground

Numerik

Montag, 28. Dezember 2009 von darth
Hallo,

mal was ganz anderes. So zur Abwechslung.
Im Moment bin ich "nur" am lernen (sollte ich jedenfalls *seufz*) und zwar ein lustiges Fach, das sich "Numerik" nennt. Das Fach ist ziemlich toll und voll elitär (naja, nicht wirklich). Eigentlich geht es darum, dass man sich einen Dreck um algebraische Korrektheit schwert und alles rein numerisch (HAHA! sprich basierend auf Zahlen) löst. Klingt spassig? Ist es auch!

Und im Zuge meiner Lernerei habe ich auch ein paar hübsche Algorithmen geschrieben, ich weiss nicht bei allen wer genau sie erfunden hat, aber ich kann nicht behaupten sie selber gefunden zu haben.

Zuerst zu etwas, das vielleicht sogar Verwendung finden könnte in irgend einem Projekt.

------------------------------

Problemstellung:
Gegeben seien eine Anzahl an Stützstellen (x_i , y_i) für i=0,2,..,(n-1). Finde ein Polynom p(x), so dass p(x_i)=y_i für alle i von 0 bis (n-1).

Ansatz 1: Lagrange
Der böse Mathematiker Mr. Lagrange hat eine Möglichkeit gefunden das zu machen. Sie ist ziemlich aufwändig von Hand und relativ mühsam, wenn man die Funktion um eine weitere Stützstelle erweitern will. Aber dennoch relativ einfach zu verstehen, man berechnet zuerst die Lagrangepolynome indem man einfach alle Terme (x-x_i)/(x_i-x_k) (für i!=k) multipliziert. Das Polynom ist dann die Summe aller Lagrangepolynome multipliziert mit dem dazugehörigen Y-Wert der Stützstelle.
Fragt mich nicht warum es funktioniert, aber es geht.
BlitzBasic: [AUSKLAPPEN]
Function LagrangeInterpolate#(X#[100], Y#[100], N, Px#)
Local L#[100], P#

For i=0 To N-1
L[i]=1

For k=0 To N-1
If i<>k
L[i]=L[i]*(Px-X[k])/(X[i]-X[k])
EndIf
Next

P=P+L[i]*Y[i]
Next

Return P
End Function

Diese Funktion liefert den interpolierten Y-Wert eines Punktes mit gegebenem X-Wert zurück.

Ansatz 2: Newton
Newtons Ansatz ist einiges eleganter und viel einfacher um weitere Stützstellen zu erweitern (von Hand, dem Computer ist das relativ egal :>). Man berechnet sukzessiv die dividierten Differenzen und multipliziert dann mit dem Stützstellenpolynom. Es ist schwer in Worte zu fassen, und um es intelligent hinzuschreiben fehlt mir die Geduld Very Happy, wer will, kann hier nachlesen LINK 1 oder hier LINK 2.
BlitzBasic: [AUSKLAPPEN]
Function NewtonInterpolate#(X#[100], Y#[100], N, Px#)
Local A#[100], P#

For i=0 To N-1
A[i]=Y[i]
Next

For k=1 To N-1
For i=N-1 To k Step -1
A[i]=(A[i]-A[i-1])/(X[i]-X[i-k])
Next
Next

P=A[N-1]
For k=(N-2) To 0 Step -1
P=A[k]+(Px-X[k])*P
Next

Return P
End Function

Auch diese Funktion liefert den Y-Wert zu einem bestimmten X-Wert zurück.

Satz: Falls das Interpolationspolynom p(x) existiert ist es eindeutig.

Zu deutsch heisst das, dass es eigentlich völlig egal ist ob ihr nun Newton oder Lagrange bevorzugt, das ergebnis ist das gleiche.

Weiter im Text. Nullstellen suchen ist ziemlich wichtig in der Analysis. Um einen Extremalwert zu finden sucht man die Nullstellung in der Ableitung z.B. Allerdings können diese Funktionen ziemlich hässlich werden (sucht euch selber ein Beispiel für hässliche Funktionen, packt ein paar Cosinuse, Sinuse und exps rein und habt Spass).

------------------------------

Problemstellung:
Finde x, so, dass f(x)=0.

Es gibt Möglichkeiten die Nullstelle einer Funktion numerisch zu bestimmen. Auch hier hat der liebe Newton seine Finger im Spiel, Nullstellensuche mit dem "Newton-Rhapson-Verfahren" (kurz NewtonIteration, keine Ahnung warum man den armen Rhapson rauslässt :/ ).
Dazu sucht man sich einen Startpunkt, und geht dann tangential zur Kurve auf die X-Achse, und das, bis die Funktion nahe genug bei 0 ist.
Das Verfahren hat allerdings seine Probleme. Die Konvergenz ist nur lokal, aber da ist es O(h^2) und das ist afaik das idealste das man bekommen kann. Aber es kann halt eben vorkommen, dass es gar nicht konvergiert.
Das Verfahren ist auch in unendlich Dimensionen anwendbar, allerdings sollte man da eine Armijo-Schrittweitensteuerung einfügen, und sicherstellen, dass die Hessematrize positiv definit ist (indem man eine Cholesky-Zerlegung davon erzwingt zum Beispielt). Ich habe das in Matlab implementiert, aber in BB - NO WAY CHOSE!
Aber für normale f(x)=y ist es kein Problem. Hier der Code dazu:
BlitzBasic: [AUSKLAPPEN]
Function NewtonIter(X0#)
Local Tol#, Abl#, h#, X#

Tol=10^(-8)
h=10^(-4)
X=X0

While Abs(f(X))>Tol
Abl=(f(X+h)-f(X-h))/(2*h)

If Abl<Tol
Exit
EndIf

X=X-f(X)/Abl
Wend

Return X
End Function

Function f#(X#)
Return X^2-3*X+2 ;beispiel
End Function

Die Funktion f#(X#) sollte eure Funktion beinhalten, von der ihr die (also eine davon, falls es mehrere gibt) Nullstelle haben wollt. Die Ableitung braucht ihr nicht, die wird automatisch (numerisch) geschätzt. Der Rückgabewert ist die X-Koordinate der Nullstelle.

------------------------------

Und zum grossen Finale (vorläufig jedenfalls.. ich bin etwa in der Hälfte des Skripts): Eine LU-Zerlegung!
Für alle die, die nicht wissen was das ist: Man zerlegt eine Matrix A zu einer unteren Dreiecksmatrize (L) und einer oberen Dreiecksmatrize (U), so dass A=L*U. Der Sinn davon ist, dass man so ein lineares Gleichungssystem ziemlich einfach berechnen kann.

Problemstellung:
Gegeben ist ein lineares Gleichungssystem A*x=b. Finde x.

Zuerst löst man die LU-Zerlegung. Wie das geht werde ich hier nicht genau erklären, aber es ist im Prinzip eine Gauss-Elimination, genauer nachzulesen hier LINK 1 oder hier LINK 2.
Danach löst man L*y=b, das ist ziemlich einfach, da L eine untere Dreiecksmatrize ist, und man nur vorwärts substituieren muss. Und nachher löst man U*x=y, selbes Prinzip, im Kopfstand. Und wie man sieht, hat man dann mit x seine Lösung und kann glücklich sein.
Hier der Code:
BlitzBasic: [AUSKLAPPEN]
;Maximales N: 100
;
; +-----> j
; |
; | A_ij
; |
; i v
;
; Matrizen und Vektoren Indizes beginnen standardmässig bei 1
; (nicht bei 0, es sind keine Arrays)

Function GetI(i, j, N) ;2D -> 1D
Return j+N*i
End Function

Function GetIJ(idx, N, idx[2]) ;1D -> 2D
idxI=idx/N
idxJ=idx-idxI*N

idx[1]=idxI
idx[2]=idxJ
End Function

Function LU(A#[1000], N, L#[1000], U#[1000])
For i=1 To N
For j=1 To N
L[GetI(i, j, N)]=0

If i=j
L[GetI(i, j, N)]=1
EndIf
Next
Next

For k=1 To (N-1)
For i=(k+1) To N
L[GetI(i, k, N)]=A[GetI(i, k, N)]/A[GetI(k, k, N)]

For j=k To N
A[GetI(i, j, N)]=A[GetI(i, j, N)]-L[GetI(i, k, N)]*A[GetI(k, j, N)]
Next
Next
Next

For i=1 To N
For j=1 To N
U[GetI(i, j, N)]=A[GetI(i, j, N)]
Next
Next
End Function

Function ForwardSolve(A#[1000], b#[100], N, X#[100])
Local sum#

X[1]=b[1]/A[GetI(1, 1, N)]

For i=2 To N
sum=0
For j=1 To (i-1)
sum=sum+A[GetI(i, j, N)]*X[j]
Next

X[i]=(b[i]-sum)/A[GetI(i, i, N)]
Next
End Function

Function BackwardSolve(A#[1000], b#[100], N, X#[100])
Local sum#

X[N]=b[N]/A[GetI(N, N, N)]

For i=(N-1) To 1 Step -1
sum=0
For j=i+1 To N
sum=sum+A[GetI(i, j, N)]*X[j]
Next

X[i]=(b[i]-sum)/A[GetI(i, i, N)]
Next
End Function

Function Solve(A#[1000], b#[100], N, X#[100])
Local L#[1000], U#[1000], Y#[100]

LU(A, N, L, U)

ForwardSolve(L, b, N, Y)
BackwardSolve(U, Y, N, X)
End Function

Function PrintMatrix(A#[1000], N)
Local z$

For i=1 To N
For j=1 To N
z$=z$+A[GetI(i, j, N)]+" "
Next

Print z$
z$=""
Next
End Function

Function PrintVector(X#[100], N)
For i=1 To N
Print X[i]
Next
End Function

Der Matlab Code ist einiges schöner, da ich dört einfacher mit Matrizen hantieren kann und auch mehrere Werte auf einmal zurückliefern kann.
Ein paar Dinge zum meinem Code: Man kann maximal 100x100 Matrizen lösen (dank BBs interessantem Arraysystem). Ich muss ständig von 2D Matrixkoordinaten zu 1D Arraykoordinaten casten, das kostet Zeit. Matrizen/Vektor Indizes fangen bei 1 an, es sind KEINE Arrays und beginnen daher nicht bei 0! (Und vor allem wird es in Matlab so gehandhabt :> ).

Ein paar Dinge zum Algorithmus allgemein: Die LU-Zerlegung ist ziemlich effektiv für kleine Matrizen die man genau lösen will. Die Problematik dabei ist, dass man nicht irgendwelche Matrizen nehmen kann. Zuersteinmal müssen die Matrizen regulär sein, d.h invertierbar (das ist ziemlich logisch, nicht reguläre Matrizen haben nur eine triviale Lösung oder gar keine). Zudem sollten die Matrizen strikt diagonal dominant sein (kurz: SDD) und natürlich R(nxn). Es gibt eine Möglichkeit auch nicht SDD Matrizen zu lösen, dazu führt man eine Pivotsuche ein, aber so weit hab ich noch nicht gelesen Smile das kommt in den nächsten paar Seiten und ist absolut MÜHSAM!
Für grössere Matrizen sollte man sich von der LU-Zerlegung fernhalten, die Zerlegung ist O(n^3) (iirc). Ausserdem gibt es selbst für sparse (dünnbesetzte) Matrizen für L und U einen riesen fill-in, so dass der Speicherverbrauch bis zum Overflow ansteigt. Für grössere Matrizen benutzt man im Normalfall das CG-Verfahren (hier müssen die Matrizen SPD sein, aber das Verfahren konvergiert sehr gut) oder das GMRES-Verfahren (für beliebige Matrizen, aber langsamer), je nach Wunsch kann man die Matrix auch noch (mit Genauigkeitsverlust) vorkonditionieren.

------------------------------

Soviel für heute. Vielleicht kommen noch ein paar andere Algorithmen hinzu, das Schreiben dieser Einträge ist eine lustige Übung für die kommende mündlich Prüfung. Ausserdem macht es Spass mit der "Luft" zu fachsimpeln Very Happy
Ich weiss nicht wie gut die Algorithmen hier ankommen, weil ich den Hintergrund der Forenbenutzer nicht kenne. Allerdings vermute ich, dass einige Dinge dieses Eintrags den mathematischen Erfahrungshorizont einiger Leute hier weit übersteigt. Dabei muss man sich nicht schämen, ich hatte das Zeug im zweiten Semester an der Universität und hatte damals einige Mühe es zu verstehen. Aber nun muss ich es ja, weil die Prüfung sonst ziemlich peinlich wird Very Happy

Das wars von meiner Seite, ich wünsche der Bevölkerung eine schöne Nacht,
MfG,
Darth

Jayks und Hurray!

Sonntag, 20. Dezember 2009 von darth
Hallo,

ich habe nun etwas Zeit investiert neben dem ausgehendem Studium von "Justice League of America" (sehr gute Serie übrigens!) um meine beiden Stücke zu vereinheitlichen. Das ganze sieht dann etwa so aus:

user posted image

Weil ich scheinbar nicht fähig bin in 100% der Fällen links von rechts zu unterscheiden habe ich beschlossen, die ganzen Zahlen einfach anzuschreiben, ein Hoch auf Paint!

Aber genug davon, auf zu etwas genau Gleichem:
Ich habe die ganze Physik eigentlich von Anfang an darauf ausgelegt, um sie in mein Spiel zu integrieren. Das wird aber (wahrscheinlich) noch ein ziemliches Stück arbeit werden, weil ich auch da wieder zig verschiedene Datenarten verwendet hab *seufz* ... Designfail nennt man das oder in meinem Fall: überhaupt kein Design, und das rächt sich nun :> Wie auch immer!

Um zu sehen wie das ganze funktionieren könnte, habe ich kurz ein Testspielchen geschrieben, es ist nur ein Level und wird als solches vielleicht zu Amusement-Zwecken von mir mal weitergeführt, aber das steht in den Sternen.
PS: Wusstet ihr, dass man aufgrund der Präzession der Erde sein Sternzeichen zwei Monate später wählen müsste? Es ist ein Zyklus von 26000 Jahren, bis wir wieder am selben Ort stehen. Und da seit der Einführung der Sternzeichen durch die Babylonier vor etwa 4000 Jahren etwa 2/12 dieses Zyklus vergangen sind kommt es auf diese 2 Monate :> Mechanik ist doch etwas Tolles!

Wo war ich? Ach ja, bei meinem Spielchen. Hier ist ein mal Bild wie das Ding aussieht:
user posted image

Das Ziel des Spiels ist es den blöde hüpfenden Pfeil zu erreichen, damit er nicht weiter herumhypert! (Stattdessen wird es durch eine nervöse Schrift ersetzt, die euch erklärt dass ihr gewonnen habt :> Phun!)
Die Steuerung ist denkbar einfach. Mit den Pfeilen <- und -> kann man seine Kugel ... eher seine Scheibe, wie auch immer ... hin und her bewegen, solange ihr Boden unter den Füssen habt geht es ziemlich schnell, in der Luft ist die Steuerung begrenzt. Springen kann man mit der Pfeil rauf (ich hab leider kein schönes Zeichen dafür, î sieht doof aus!) Taste, solange man irgendwo aufm Boden steht.
Wenn ihr euch so verfranst habt, dass ihr nichtmehr weiterkommt, dann könnt ihr mit R das Level resetten.

Der Link dazu ist hier.

Soviel sollte vorerst reichen. Mehr dazu später, nach den Ferien - wahrscheinlich. Ich wünsche denjenigen die feiern schöne Feiertage (immer schön religionsneutral bleiben, obwohl ich bezweifle dass sich Andersgläubige über zusätzliche Ferien ärgern) und demnächst ein schönes neues Jahr.

MfG,
Darth

Altes und Neues

Freitag, 18. Dezember 2009 von darth
Hallo,

ich habe im Moment keine Vorlesungen mehr und muss eigentlich "nur" noch auf Prüfungen lernen, aber das hat Zeit bis nach Weihnachten. Bis dahin habe ich vor, wieder etwas an meiner Physik zu arbeiten.

Zuerst zu etwas älterem, nach Monty Python: "And now for something completely similar..". Ich hatte ursprünglich vor, meine Physik aus Geschwindigkeitsgründen zu C++ als dll zu konvertieren, aber da hab ich im Moment nicht die Motivation dazu (ich hab noch anders zu tun das ich in etwa gleich wenig kann :> ), daher dachte ich mir, ich machs mal in Java, dann seh ich wenigstens wie das OOP-mässig aussehen könnte. Ich bin mir nicht ganz sicher was der Stand des letzen Uploads war, aber möglich sind hier:
- 3D
- konkave und konvexe Objekte (die konkaven Objekte werden intern automatisch in konvexe zerlegt)
- das wars

Ich bin nicht ganz sattelfest mit *.jar Files, ich hoffe ich habe es richtig gemacht, das ganze benötigt JOGL, eine freie OGL Library für Java, auch bei Java und DLLs bin ich nicht ganz sicher, aber bei mir läuft es, ich hoffe das geht bei euch auch.
ANM: Ich weiss dass es auf MacOS Probleme macht (ein Kollege hat es getestet), das liegt daran, dass die JOGL Lib gewisse C Stücke enthält, die für Windows ausgelegt sind.

Bild und Link:
user posted image
Links oben: Anz. Objekte, rechts unten: (theoretische) FPS

Und nun zu etwas Neuem. Mein Wasser hat ein kleines Update erhalten, auch da weiss ich nichtmehr auf welchem Stand der letzte Upload war, von daher nochmal: Das Wasser ist nun Grid-basierend und von daher ziemlich viel schneller.
Ausserdem bin ich dabei, die Interaktion zwischen Wasser und Objekten einzufügen, das funktioniert eigentlich schon, mit dem kleinen Nachteil, dass ich bei Physik und Wasser zwei verschiedene "Datentypen" verwende, von daher interagieren die Objekte nicht miteinander.

Vorerst mal ein Bild:
user posted image
Links oben: Anz. Partikel, Mitte oben links: (theoretische FPS), Mitte oben rechts: (gebremste) FPS

Das Dreieck schwimmt auf dem Wasser, die anderen Objekte sind unbeweglich.

Hier der Link dazu.
Mit der linken Maustaste kann man neues Wasser erstellen, mit der rechten kann man Objekte packen und verschieben (hierbei spielt es keine Rolle, ob man bewegliche oder statische Objekte packt).

Ausblick auf Folgendes:
Was jetzt noch kommt ist die Vereinheitlichung -.- Wie oben erwähnt benutze ich im Wassercode einen Type "Polygone" und im Physikcode "Body", auch haben die beiden ganz andere Ansätze, die Polygone sind Punktbasiert und Bodys sind, ähnlich wie Meshes, Objektbasiert, also mit Zentrum, Rotation und den Punkten relativ dazu. Es wird ziemlich lustig werden den ganzen Gerümpelcode zu einem einzelnen Stück zu verschmelzen, ich freu mich drauf \o/
Ausserdem ist das Wasser eher "empirisch" aufgebaut, Partikel suchen sich einfach freie Stellen und schlüpfen dort hin, ich habe eigentlich vor, da gewisse Gleichungen der SPH-Methodik einzubauen, ich hab zwar eine lauffähige Version, die sich ziemlich stark an Noobodys Variante orientiert, aber die hat einige grobe Macken, und im Moment gefällt mir mein Ansatz besser..

Und zum Schluss:
Wer Lust hat ist natürlich gerne dazu eingeladen Bugs zu suchen und allgemein Belastungstests oder ähnliches durchzuführen.

Ich wünsche schon mal schöne Feiertage und verabschiede mich.
MfG,
Darth

Grund und so

Sonntag, 15. November 2009 von darth
Hallo,

ich habe auch einen Worklog \o/ endlich habe ich den Knopf gefunden, mit dem man Worklogs erstellen kann. Uhm, nun ja, um etwas ernster zu werden:

Wozu dient dieser Worklog?
Ich hab mittlerweilen keine Lust mehr, meine Dinge zwischen Projektforum/WiP/Codearchiv/Abfalleimer zu sortieren. Deshalb habe ich mich entschlossen, einen Worklog zu schreiben, indem ich einfach mal sammeln kann, was ich auf die Beine stelle. Ausserdem tut hier wohl niemand blöd, wenn man etwas unfertiges postet, es ist ja (tadaa!) ein Worklog.

Was kann ich von dem Worklog erwarten?
Im Prinzip nicht viel. Dieser Worklog dient in erster Linie mir selber. Ich brauche halt einfach irgend einen Ort wo ich meinen Fortschritt dokumentieren kann, da ich im Moment meine Zeit zwischen Uni und Freizeit teilen muss und manchmal den Faden verliere. Darum dieser Worklog...
Andere sind natürlich herzlich eingeladen dem Worklog zu folgen, sollten sie Interesse daran finden.

Man sollte nicht damit rechnen, dass alle Dateien bis in alle Ewigkeit verfügbar sein werden, ich lade das meiste auf den Forenspace hoch, und das ist auf 10MB begrenzt, ich muss also hin und wieder löschen. Aber meistens sind das dann alte Dinge, die sowieso nicht länger von Interesse sind.

Was gibt es hier zu sehn?
In erster Linie werde ich hier wohl meine Versuche einer Physikengine hineinstellen. Dies beinhaltet meinen Ansatz der Wassersimulation, sowie meine RigidBody-Physics. Gesammelt wird das ganze unter DaDaPhysics (... hihi ...), v.a weil es mal zusammen verbunden werden soll.
Wenn die Physik soweit ist wie ich sie gerne hätte (und da müssen leider noch ein paar Dinge verändert werden) werde ich die Entwicklung meines Spiels wieder aufnehmen und wohl auch hier posten, weil ich keine Lust auf 2 Worklogs habe (Grund: siehe "Wozu dient dieser Worklog?").
Vielleicht werden hier auch so zwischendurch ein paar andere Versuche veröffentlicht, das wird sich noch zeigen.

Was gibt es jetzt schon?
Wie gesagt habe ich bereits Ansätze von Wasser und Physik Simulationen, das findet man irgendwo verstreut zwischen Codearchiv und WiP, aber hier nochmal die Sammlung:

Wasser:
user posted image

Physik:
user posted image

Das meiste was hier landen wird, wird Opensource sein. (bis auf das Spiel denk ich mal... es macht keinen Spass wenn man den Code sehen kann und weiss wie die Geschichte laufen wird oder so :/ ausserdem könnte man Spielstände manipulieren). Die Physik soll eigentlich als eine Engine entwickelt werden, die man ohne grosse Umstände in BB verwenden kann. Ich arbeite an einer C++ DLL, aber ich bin kein Fan von C++, deshalb wird das wohl eine Weile dauern Razz

GG and out!

MfG,
Darth

Gehe zu Seite Zurück  1, 2, 3