nicht angemeldet
Gunnar Bittersmann
JürgenB
Tabellenkalk
JürgenB
Camping_RIDER
Camping_RIDER
Camping_RIDER
Hallo liebe Forianer,
ich stehe gerade vor einem Problem, dass ich mir einen Lösungsalgorithmus für ein lineares Gleichungssystem schreiben muss (in einer Fortran77-basierten Sprache). Und da ich durchs Mitlesen in diesem Forum hier schon mitbekommen habe, dass hier doch einige Programmierer/Mathematiker o.ä. anwesend sind, dachte ich, ich stelle mein Problem hier kurz einmal vor - eventuell hat der eine oder andere einen (oder mehrere) Hinweise für mich.
Nun handelt es sich dabei ein 4x4-Gleichungssystem, mit zum Teil sehr "ungeschickt" skalierten Matritzen, wie (wobei A = Koeffizientenmatrix und b = Störvektor) z.B.:
A = [ 3.91126953294906e-13 0.000290530407713016 0.118719995662381 0.999861140092325;
-7.47034897265332e-13 -0.000158174186504138 -0.0169130360585405 -9.28246904902663e-06;
1.42680307002891e-12 8.61151625166872e-05 0.00240945754016845 8.61761980650444e-11;
-2.72512972030656e-12 -4.68838903437727e-05 -0.000343255085472549 -8.00038984641107e-16 ]
b = [ -0.000272415115940619; -9.83039417624391e-07; 1.42984225458907e-07; 3.29280620635692e-09 ]
Hierfür (aus Matlab: rank(A) = rank([A, b]) = 4, daher eindeutig lösbar) habe ich ein bestehendes Skript (Gauss-Elimination mit Pivoting) modifiziert bzw. in diese Fortran77-basierte Sprache portiert. Das funktioniert soweit auch ganz gut, wenngleich zunaechst die Matrix etwas skaliert werden muss, um die Zahlen dem float-Format etwas zugänglicher zu machen.
Das Problem tut sich nun auf, dass sich die Einträge des Gleichungssystem sukzessive verändern, bis irgendwann ein folgender Gestalt gelöst werden muss:
AA = [ 9.71932695068104e-21 1.95542531016539e-06 0.0319594622483047 0.999775498752868;
-7.90937295921570e-21 -4.53931579859449e-07 -0.00194341016115096 -3.96425686254142e-06;
6.43647249705794e-21 1.05375479248690e-07 0.000118176051434193 1.57188613761892e-11;
-5.23785873026922e-21 -2.44618178588262e-08 -7.18612026002067e-06 -6.23275967050862e-17 ]
bb = [ 4.42616984297712e-05; -1.36475494887893e-06; 9.48249135038941e-08; 3.66966454607991e-09 ]
Hierfür erhalte ich in Matlab rank(AA) = 3 und rank([AA, bb]) = 4, das GLS ist also überbestimmt.
Nun stellen sich mir folgende Fragen:
- Ist das GLS AA,bb tatsächlich überbestimmt, oder ergibt sich dieses Problem aus der ungünstigen Skalierung der Matrix bzw. dieser schon "sehr kleinen" Einträge? (Eine zunächste Skalierung mit 10e+22 führt aber zum selben Ergebnis; weshalb das GLS tatsächlich überbestimmt sein dürfte.)
- Welches Verfahren eignet sich für die Lösung des GLS AA,bb. Ich dachte hier an eine QR-Zerlegung mit kleinsten Fehlerquadraten zum Finden des Lösungsvektors.
- Gibt es einen fertigen Sourcecode o.ä., mit dem beide Probleme behandelt/gelöst werden können? Ich bin leider kein Programmierer/Mathematiker und weiß daher nicht genau, ob ich selbst einen Löser so schreiben kann, dass dieser effizient und sicher die beiden GLS behandeln kann (noch dazu, wenn man auf das float-Zahlenformat beschränkt ist).
Über Rückmeldungen wäre ich sehr dankbar!
Liebe Grüße jakob
Hallo,
zum Teil sehr "ungeschickt" skalierten Matritzen,
Meinst du damit den Wertebereich von e-1 bis e-21? Wo kommen denn diese Werte her, wenn man fragen darf? Mir scheint da ein Konzeptproblem vorzuliegen.
Gruß Kalk
Hallo,
zum Teil sehr "ungeschickt" skalierten Matritzen,
Meinst du damit den Wertebereich von e-1 bis e-21? Wo kommen denn diese Werte her, wenn man fragen darf? Mir scheint da ein Konzeptproblem vorzuliegen.
Gruß Kalk
Tja, dabei handelt es sich um Funktionsparameter die das zeitabhängige Materialverhalten eines Festkörpers beschreiben, welche wiederum in einer Differentialgleichung auftreten und aus einem Umfassenden Versuchsprogramm im Rahmen einer Datenregressionsanalyse gewonnen wurden ...
Gruß, jakob78
Tja, dabei handelt es sich um Funktionsparameter die das zeitabhängige Materialverhalten eines Festkörpers beschreiben …
Physiker! Das sind dann wohl auch die einzigen, die noch FORTRAN verwenden. ;-)
LLAP
Hallo Gunnar,
Physiker! Das sind dann wohl auch die einzigen, die noch FORTRAN verwenden. ;-)
und sie wissen, warum.
Gruß Jürgen
Hallo,
Physiker!
und sie wissen, warum.
Aber warum verwenden sie Zahlen mit einer scheinbaren Genauigkeit von 14 oder mehr Stellen?
Gruß Kalk
Hallo
Aber warum verwenden sie Zahlen mit einer scheinbaren Genauigkeit von 14 oder mehr Stellen?
wenn es um die Simulation von Experimenten geht, werden in den meisten Fällen 3 Stellen reichen. Präzisionsmessungen erfordenr etwa 10 oder 11 Stellen.
Da beim Lösen von Differentialgleichungen meistens iterative Verfahren eingesetzt werden und in Folge dessen sehr viele sequentielle Rechenschritte durchgeführt werden müssen, muss die Genauigkeit der einzelnen Schritte schon deutlich höher sein. Stell dir vor, du hast eine Iteration U_n+1 = F(U_n) und lässt n bis 10^10 laufen, dann bleiben von 14 Stellen u.U. nur noch 4 übrig.
Gruß Jürgen
PS Warum „scheinbaren Genauigkeit“?
PPS Die mir bekannten Programmiersprache unterstützen alle das Datenformat „Single Precision“ mit 32 Bit Speicherbelegung und 7 Stellen Genauigkeit sowie „Double Precision“ mit 64 Bit Speicherbelegung und 15 Stellen Genauigkeit.
Hallo,
PS Warum „scheinbaren Genauigkeit“?
Weil Zahlen mit sovielen geltenden Ziffern suggerieren, dieser Wert ist so. Dabei könnte man z.B. grade bei Materialbeschaffenheit froh sein, wenn die Meßgenauigkeit im Promillbereich liegt.
Gruß Kalk
Hallo
PS Warum „scheinbaren Genauigkeit“?
Weil Zahlen mit sovielen geltenden Ziffern suggerieren, dieser Wert ist so. Dabei könnte man z.B. grade bei Materialbeschaffenheit froh sein, wenn die Meßgenauigkeit im Promillbereich liegt.
ach, du beziehst dich auf die Daten von Jakob und nicht auf die Rechengenauigkeit von Fortran bzw. von Programmiersprachen im Allgemeinen. Das habe ich falsch verstanden.
Gruß Jürgen
Aloha ;)
Weil Zahlen mit sovielen geltenden Ziffern suggerieren, dieser Wert ist so. Dabei könnte man z.B. grade bei Materialbeschaffenheit froh sein, wenn die Meßgenauigkeit im Promillbereich liegt.
ACK. Aber bedenke: Die Werte, die in der Matrix vorkommen, müssen nicht unbedingt in dieser Form real messbar sein. Trotzdem ist es u.U. im entsprechenden Verfahren nötig, die zwischen Messdaten und Ergebnissen liegenden Zwischenschritte entsprechend genau durchzuführen. Es ist einfach - wie an vielen Stellen in der Wissenschaft - eine Frage der begründeten Abschätzung, welche Genauigkeit an welcher Stelle nötig ist.
Grüße,
RIDER
Hallo,
ACK. Aber bedenke: Die Werte, die in der Matrix vorkommen, müssen nicht unbedingt in dieser Form real messbar sein. Trotzdem ist es u.U. im entsprechenden Verfahren nötig, die zwischen Messdaten und Ergebnissen liegenden Zwischenschritte entsprechend genau durchzuführen. Es ist einfach - wie an vielen Stellen in der Wissenschaft - eine Frage der begründeten Abschätzung, welche Genauigkeit an welcher Stelle nötig ist.
Da hast du in der Hinsicht recht. Allerdings fehlt mir diese Begründung für die Beispieldaten, bzw. sind ja schon vom TO selbst Zweifel angemeldet worden. Da werden in der Matrix nicht Daten miteinander verrechnet, die mal eben um den Faktor 20 differieren, sondern es geht um 20 Größenordnungen.
Und das halte ich halt für ein Konzeptproblem.
Gruß Kalk
Aloha ;)
Da hast du in der Hinsicht recht. Allerdings fehlt mir diese Begründung für die Beispieldaten, bzw. sind ja schon vom TO selbst Zweifel angemeldet worden. Da werden in der Matrix nicht Daten miteinander verrechnet, die mal eben um den Faktor 20 differieren, sondern es geht um 20 Größenordnungen.
Und das halte ich halt für ein Konzeptproblem.
Ja und Nein. Auf den ersten Blick ja. Auf den zweiten Blick fällt mir allerdings auf, dass die Größenordnungen so schon in Ordnung sein könnten. Die Werte mit der kleinsten Größenordnung stehen alle in der ersten Spalte und sind unter sich betrachtet innerhalb einer Größenordnung. Zwar ist der Unterschied zu den anderen Spalten sehr groß, das kann aber egal sein - Matrixmultiplikation wird ja als Zeilen/Spalten-Skalarprodukt ausgeführt. "Miteinander verrechnet" ist also imho nicht unbedingt der passendste Ausdruck. Die Daten in einer Matrix sind viel eher als "in Bezug zueinander stehend" wahrzunehmen, als tatsächlich als Werte, die miteinander verrechnet werden. Gerade bei physikalischen Problemen ist die Skalierung von Spalten oft auch eine Frage der verwendeten Einheit, die reine Zahl bietet da nicht unbedingt Aufschluss.
Ist aber unerheblich, der TO wird schon wissen, ob seine Daten für sein Problem valide sind oder nicht, wollte das nur angemerkt wissen ;)
Grüße,
RIDER
Aloha ;)
Spontane Frage: Warum nicht mit LR-Zerlegung und Vorwärts-/Rückwärtseinsetzen?
Immerhin ist das eine 4x4-Matrix, die ist vom Rechenaufwand her eigentlich recht unproblematisch...
Wenn Interesse besteht: Ich dürfte aus meinen Numerik-Übungs-Beständen eine fertige Implementierung eines LGS-Lösers mit LR-Zerlegung für Matlab haben. Ich such das mal raus...
Grüße,
RIDER
Aloha ;)
Ups..
Spontane Frage: Warum nicht mit LR-Zerlegung und Vorwärts-/Rückwärtseinsetzen?
Weil die Matrix nicht diagonaldominant ist und die LR-Zerlegung gar nicht existiert... ja, ich gehe mich schämen und schau trotzdem meine Unterlagen mal durch, ob vllt. was brauchbares dabei ist.
Grüße,
RIDER
Aloha ;)
Zweiter Anlauf, hoffentlich hilfreicher.
Ich frage mich gerade: Warum so kompliziert und die Löserei nicht direkt durch Matlab machen lassen?
Du hast die Gleichung A * x = b
Du hättest gerne x, also x = inv(A) * b (wenn inv(A) die Linksinverse von A ist)
In Matlab-Syntax ist das dann: x = A \ b
Voila ;) Und imho dürfte das auch so ziemlich die genauste Berechnung sein, die du mithilfe von Matlab erhalten kannst. Kann mich auch täuschen, aber meiner Erfahrung nach sind die internen Implementationen schon sehr genau. Zumal imho keine numerische Auslöschung auftreten sollte, da weder bei Matrixmultiplikation noch bei Inversenberechnung großartig Differenzen gebildet werden - höchstens vielleicht bei der Berechnung der Determinante. Alles andere ist Maschinengenauigkeit.
Grüße,
RIDER
Hallo,
Ich frage mich gerade: Warum so kompliziert und die Löserei nicht direkt durch Matlab machen lassen?
Du hast die Gleichung A * x = b
Du hättest gerne x, also x = inv(A) * b (wenn inv(A) die Linksinverse von A ist)
In Matlab-Syntax ist das dann: x = A \ b
Jaa ... in Matlab ist das ja kein Problem. Soweit hab ich das auch. - Es geht nur darum, dass ich diese ganze Gleichungslösereiprozedur in eine andere Sprache (auf Fortran77 basiert) portieren muss, die keinerlei Funktionalitäten kennt, d.h. weder Gleichungslöser o.ä. und ich einen selbst schreiben muss.
Nun war ich am überlegen, was sich hierfür am besten anbietet um sowohl bestimmte als auch überbestimmte LGS zu lösen ... - und da ich leider kein Mathematiker bzw. Programmierer bin bringen einen vor allem diese ungeschickt skalierten Matritzen ins Grübeln ...
LG jakob78
Aloha ;)
Jaa ... in Matlab ist das ja kein Problem. Soweit hab ich das auch. - Es geht nur darum, dass ich diese ganze Gleichungslösereiprozedur in eine andere Sprache (auf Fortran77 basiert) portieren muss, die keinerlei Funktionalitäten kennt, d.h. weder Gleichungslöser o.ä. und ich einen selbst schreiben muss.
Ah, entschuldige. Dann hatte ich dein Problem misverstanden. Okay. Du suchst also einen möglichst einfach zu implementierenden Algorithmus, der dir ein lineares Gleichungssystem mit einer 4x4-Matrix löst. Okay. Was spricht gegen einen Einsatz der Cramerschen Regel? Das ist zwar keine sehr performante Lösung, das spielt bei 4x4-Matrizen allerdings noch keine wirkliche Rolle - in dem Fall geht es ja vor allem darum, eine Lösung zu erhalten.
Der Aufwand im Programmieren beschränkt sich darauf, eine Funktion zu schreiben, die die Determinante einer 4x4-Matrix berechnet, z.B. über den Laplaceschen Entwicklungssatz. Die Funktion zum Lösen der Gleichung tut dann nichts anderes, als (nach der cramerschen Regel) je zwei Determinanten pro Variable zu berechnen.
Dass das Verfahren unperformanter ist als Gauß-Elimination, QR-Zerlegung und Co. wirst du bei einer 4x4-Matrix noch nicht wirklich spüren und der Aufwand zur Implementierung ist bei der Cramerschen Regel wirklich deutlich (!) geringer.
Zum Thema Genauigkeit gilt in etwa das Gleiche, was ich vorher bei der Matlab-Lösung gesagt hatte.
Ich bin ziemlich sicher, dass das Problem mit relativ wenig Zeilen Programmcode abzufrühstücken ist - Ich würds ja für dich machen, nur hab ich leider keinerlei Dunst von Fortran ;)
Grüße,
RIDER
Servus
Ah, entschuldige. Dann hatte ich dein Problem misverstanden. Okay. Du suchst also einen möglichst einfach zu implementierenden Algorithmus, der dir ein lineares Gleichungssystem mit einer 4x4-Matrix löst. [...]
mmmh, fast. - Also so einen hätte ich schon; ich habe mir vor einiger Zeit eine Gauß-Elimination mit pivoting geschrieben.
Das eigentliche Problem ist jenes, dass das Gleichungssystem zunächst bestimmt (d.h. eindeutig lösbar) ist, sich im Laufe der Zeit verändert und überbestimmt wird, sodass mit Gauß-Elimination keine Lösung mehr auffindbar ist (wie auch; da für ein überbestimmtes LGS nicht zwangsläufig eine Lösung existieren muss.)
Meine Idee war z.B. ein Lösen mit den kleinsten Fehlerquadraten. Das funktioniert mit Sicherheit für die überbestimmten GLS; ich weiß aber nicht, ob man damit auch bestimmte GLS lösen kann. - Von der Grundüberlegung her müsste das auch funktionieren, allerdings habe ich hier keinerlei Literaturstellen gefunden die explizit sagen: "Mit den kleinsten Fehlerquadraten können sowohl bestimmte als auch überbestimmte GLS gelöst werden; für bestimmte GLS erhält man dann die 'exakte' Lösung; bei überbestimmten jene für die die Norm minimal ist."
LG jakob78
Die Funktion zum Lösen der Gleichung tut dann nichts anderes, als (nach der cramerschen Regel) je zwei Determinanten pro Variable zu berechnen.
Eine. Plus die Determinante der Koeffizientenmatrix.
Bei einem Gleichungssystem aus n Gleichungen mit n Variablen musst du nicht 2_n_ Determinanten berechnen, sondern n + 1.
LLAP
Bei einem Gleichungssystem aus n Gleichungen mit n Variablen musst du nicht 2_n_ Determinanten berechnen, sondern n + 1.
Markdown funktioniert – oder auch nicht.
Schön, dass es die Editierfunktion gibt. Schade nur, dass die 10 Minuten schon um sind.
LLAP
PS: Wo ist die Zitatesammlung‽
PPS: Wo ist das Interrobang in der Zeichenliste über dem Eingabefeld?
Tach!
Markdown funktioniert – oder auch nicht.
Schön, dass es die Editierfunktion gibt. Schade nur, dass die 10 Minuten schon um sind.
Die Vorschaufunktion befindet sich übrigens links neben dem Absendebutton.
dedlfix.
Aloha ;)
Bei einem Gleichungssystem aus n Gleichungen mit n Variablen musst du nicht 2_n_ Determinanten berechnen, sondern n + 1.
Das stimmt natürlich ;) Aber wenn wir schon beim Erbsen zählen sind, selbst bei Anwendung des Satz von Sarrus müssen mindestens 4 Determinanten pro 4x4-Determinante berechnet werden, also doch eher insgesamt 4n+4 Determinanten, nur mit Entwicklungssatz dann 12n + 12 Determinanten :P
Grüße,
RIDER
Sieh dir mal an wie Solver arbeiten. Vielleicht reicht es zur Ideenfindung.