Auswahl Aktion im Forum Suche Kontakt Für Mitglieder Mathematisch für Anfänger Wer ist Online | |
| Autor |
  C++: LAPACK: Unterbestimmtes LGS lösen |
|
Longfinger
Aktiv 
Dabei seit: 01.06.2005
Mitteilungen: 529
Aus: Wuppertal
 |      Themenstart: 2012-06-26 19:01
|
Hallo,
wie der Betreff schon sagt, möchte ich ein unterbestimmtes LGS Ax=b lösen, also die Koeffizientenmatrix A beseitzt einen Eigenwert 0 und ist damit singulär. Da ich an einem Grundgerüst eines alten Codes bereits arbeite wäre LAPACK optimal. Gibt es da eine Funktion, die mir das gewünschte liefert? Ich dachte zunächst ich kann es mit dgels_ machen, aber entweder ich mach da was falsch, oder das geht damit nicht.
Das Problem könnte evtl. folgendes sein:
* DGELS solves overdetermined or underdetermined real linear systems
* involving an M-by-N matrix A, or its transpose, using a QR or LQ
* factorization of A. It is assumed that A has full rank.
*
* The following options are provided:
*
* 2. If TRANS = 'N' and m < n: find the minimum norm solution of
* an underdetermined system A * X = B.
Bei mir ist A eine MxN Matrix, aber M ist gleich N, man sieht also nicht direkt, dass A singulär ist.
EDIT: Jetzt seh ichs erst. A soll vollen Rang haben. Das hat es ja dann bei mir nicht. Jetzt frage ich mich, ob ich leicht die linear abhängigen Zeilen bestimmenn kann und diese manuell vor Aufruf der Funktion rauslöschen kann.
Je länger ich darüber nachdenke, destmo mehr verwirrt mich das ganez.
Ganz konkret:
Dann ist dieses LGS Ax=b doch gar nicht lösbar, da sich Zeile 1 und 3 doch widersprechen, oder?
Viele Grüße
Longfinger
[ Nachricht wurde editiert von Longfinger am 26.06.2012 19:28:20 ]
[ Nachricht wurde editiert von Longfinger am 27.06.2012 11:32:46 ]
|
 Profil
 Quote
 Link |
viertel
Senior 
Dabei seit: 04.03.2003
Mitteilungen: 21549
Aus: Hessen
 | Beitrag No.1, eingetragen 2012-06-26 20:18
|
Profil
Quote
Link |
davidhigh
Senior
Dabei seit: 10.03.2007
Mitteilungen: 2482
Aus: Kiel
 | Beitrag No.2, eingetragen 2012-06-27 00:19
|
Hallo,
du könntest natürlich erstmal eine Singular-Value-Decomposition machen, und dann eine Matrix mit vollem Rang konstruieren, auf die du dann dgels loslässt. Das ist aber natürlich alles andere als angenehm...
Gruß, David
-----------------
Eine Hauptursache der Armut in den Wissenschaften ist meist eingebildeter Reichtum
(Bertolt Brecht, Leben des Gallilei)
|
Profil
Quote
Link |
Longfinger
Aktiv
Dabei seit: 01.06.2005
Mitteilungen: 529
Aus: Wuppertal
| Beitrag No.3, vom Themenstarter, eingetragen 2012-06-27 11:32
|
Hallo
@ Viertel: Ich habe mich verschrieben, die 2x2 Matrix sollte nochmal invertiert werden.
@ alle: Es liegt anscheinend im Moment kein Implementierungsfehler vor, sondern vermutlich ein Verständisfehler im Algorithmus. Denn er tut genau das, was er soll nur, dass dieses LGS nicht lösbar ist, sollte nicht sein. Ich melde mich nochmal wenn ich mehr weiß.
Viele Dank.
[ Nachricht wurde editiert von Longfinger am 27.06.2012 11:33:07 ]
|
Profil
Quote
Link |
viertel
Senior
Dabei seit: 04.03.2003
Mitteilungen: 21549
Aus: Hessen
| Beitrag No.4, eingetragen 2012-06-27 18:14
|
2012-06-26 20:18 - viertel in Beitrag No. 1 schreibt:
Und was meinst du mit „Zeilen widersprechen sich“? Ok, mit den richtigen Zahlen seh ich es jetzt auch:
  
0.25*(2, 4; 2, 1; 1, 2)*(5, 3; 3, 10)^(-1)*(2, 2, 1; 4, 1, 2)=(18/41, 15/82, 9/41; 15/82, 33/164, 15/164; 9/41, 15/164, 9/82)
Und damit ist
(18/41, 15/82, 9/41; 15/82, 33/164, 15/164; 9/41, 15/164, 9/82)*(x_1;x_2;x_3)=(0;0;-3)
nicht lösbar, da die erste Zeile der Marix das Doppelte der dritten ist, was mit der rechten Seite nicht zusammen paßt.
|
Profil
Quote
Link |
Longfinger
Aktiv
Dabei seit: 01.06.2005
Mitteilungen: 529
Aus: Wuppertal
| Beitrag No.5, vom Themenstarter, eingetragen 2012-06-27 20:03
|
Danke Viertel, diese Situation darf in meinem Algorithmus nicht auftauchen, glaube ich. Ich weiß nicht wieso das so passiert. Es handelt sich um den Active Set Algorithmus für Quadratische Programme. Ich habe ein Beispiel durchgerechnet, bei dem der programmierte Algorithmus den gleichen 'Fehler' liefert wie mit Hand durchgerechnet. Wäre jemand bereit meine Lösung mal auf Fehler zu durchsuchen? Ich habe alles gefühlte 100 mal durchgerechnet und weiß langsam nicht mehr weiter.
|
Profil
Quote
Link |
viertel
Senior
Dabei seit: 04.03.2003
Mitteilungen: 21549
Aus: Hessen
| Beitrag No.6, eingetragen 2012-06-28 01:22
|
2012-06-27 20:03 - Longfinger in Beitrag No. 5 schreibt:
diese Situation darf in meinem Algorithmus nicht auftauchen, glaube ich. Da niemand weiß, wie du zu den Matrizen kommst, kann auch schlecht jemand was dazu sagen.
|
Profil
Quote
Link |
Longfinger
Aktiv
Dabei seit: 01.06.2005
Mitteilungen: 529
Aus: Wuppertal
| Beitrag No.7, vom Themenstarter, eingetragen 2012-06-28 11:17
|
Ich würde dann einen neuen Thread eröffnen, da das nichts mehr mit Programmieren zu tun hat, sondern eher mit Optimierung. Dann würde ich natürlich auch den Algorithmus wiedergeben. Die Frage war nur, ob jemand bereit wäre so etwas zu kontrollieren, da die Rechnungen zum Teil aufwendig sind (auf den ersten Blick).
|
Profil
Quote
Link |
lehn
Junior 
Dabei seit: 09.02.2008
Mitteilungen: 12
Aus: Ulm
 | Beitrag No.8, eingetragen 2012-10-05 00:57
|
(2012-06-26 19:01 - Longfinger
* DGELS solves overdetermined or underdetermined real linear systems
* involving an M-by-N matrix A, or its transpose, using a QR or LQ
* factorization of A. It is assumed that A has full rank.
*
* The following options are provided:
*
* 2. If TRANS = 'N' and m < n: find the minimum norm solution of
* an underdetermined system A * X = B.
Bei mir ist A eine MxN Matrix, aber M ist gleich N, man sieht also nicht direkt, dass A singulär ist.
EDIT: Jetzt seh ichs erst. A soll vollen Rang haben.
Wenn die Matrix nicht vollen Rang hat, dann kannst du DGELSY statt DGELS benutzen. Die Routine benutzt einen Schätzer für den Rang und benutzt dann eine QR-Zerlegung mit Privatisierung.
Und natürlich habe ich das auch in FLENS eingebaut ;-)
Aber selbst wenn du FLENS nicht benutzen willst kann dieses Beispiel ganz nützlich sein. Da wird der Algorithmus skizziert.
|
Profil
 www
Quote
Link |
|
Wie unterstütze ich den Matheplaneten mit Geld?
[Neues Thema]