Gauss Algorithmus mit Pivotisierung |
16.09.2010, 16:02 | donmemento | Auf diesen Beitrag antworten » |
Gauss Algorithmus mit Pivotisierung ich möchte eine Dreieckszerlegung durchführen, meine Matrix M ist 16x16, relativ dünn besetzt. Hab mir den Algorithmus programmiert aber die Eigenwerte liefern den falschen Wert für die Determinante (berechne ich dann als Produkt der Hauptdiagonalelemente). Die Werte habe ich mit Maple überprüft. Nachdem ich das mittlerweile drei mal neu programmiert habe , seid hir meine letzte Hoffnung Hier ist mal mein Code, sieht jemand den Fehler ? ! Elimination der Elemente unter der Hauptdiagonalen (Dreieckszerlegung) do k = 1,(n - 1) ! Pivotsuche zet = 0 DO i = k,n IF (abs(M(i,k)) > zet) THEN zet = abs(M(i,k)) s = i END IF END DO p(k) = s ! Vertauschung IF (k < s) THEN DO j = 1,n zet = M(k,j) M(k,j) = M(s,j) M(s,j) = zet END DO END IF ! Berechnung der lik DO i = (k + 1),n M(i,k) = M(i,k)/M(k,k) END DO ! Spaltenweise Berechnung von aufdatiertem A DO j = (k + 1),n DO i = (k + 1),n M(i,j) = M(i,j) - M(i,k)*M(k,j) END DO END DO END DO !Berechnen der Determinante det = 1 DO j = 1,n det = det*M(j,j) END DO Schon mal vielen Dank für eure Hilfe |
||
17.09.2010, 11:48 | Formelmonster | Auf diesen Beitrag antworten » |
Wenn ich das richtig sehe hast Du vergessen Vorzeichen der Determinante zu beachten: |
||
17.09.2010, 13:26 | donmemento | Auf diesen Beitrag antworten » |
Das ist auf jeden Fall korrekt, allerdings ändert das ja nur das Vorzeichen der Determinante. Wenn ich den Eigenwert suche, für den Det=0 ist, macht das allerdings nichts aus... |
||
17.09.2010, 14:04 | tigerbine | Auf diesen Beitrag antworten » |
Kannst du deine MAtrix mal matlab konform einstellen? |
||
17.09.2010, 16:07 | donmemento | Auf diesen Beitrag antworten » |
kurze Info zu Matrix, es handelt sich um ein Eigenwertproblem (Stabknickung), die Werte EI1-EI4 (Biegesteifigkeiten), F1, F2 (Federsteifigkeiten), F4 (Drehfedersteifigkeit), l1-l4 (Längen) sind Eingabewerte. lambda ist der Eigenwert, für den die Det=0 werden soll... [0,0,0,0,1,0,0,0,1,0,0,0,0,0,0,0; 0,-(EI2*cos(lambda*l2)*lambda^3)/(EI3),lambda^3,0,0,(EI2*sin(lambda*l2)*lambda^3)/(EI3),0,0,0,0,0,0,0,0,0,0; (EI1*sin(lambda*l1)*lambda^2)/(EI2),0,0,0,(EI1*cos(lambda*l1)*lambda^2)/(EI2),-lambda^2,0,0,0,0,0,0,0,0,0,0; 0,(EI2*sin(lambda*l2)*lambda^2)/(EI3),0,0,0,(EI2*cos(lambda*l2)*lambda^2)/(EI3),-lambda^2,0,0,0,0,0,0,0,0,0; 0,0,(EI3*sin(lambda*l3)*lambda^2)/(EI4),0,0,0,(EI3*cos(lambda*l3)*lambda^2)/(EI4),-lambda^2,0,0,0,0,0,0,0,0; 0,0,0,-sin(lambda*l4)*lambda^2,0,0,0,-cos(lambda*l4)*lambda^2,0,0,0,0,0,0,0,0; (F4*lambda)/(EI1),0,0,0,-lambda^2,0,0,0,0,0,0,0,(F4)/(EI1),0,0,0; -(EI1*cos(lambda*l1)*lambda^3)/(EI2),lambda^3,0,0,(EI1*sin(lambda*l1)*lambda^3)/(EI2),-(F1)/(EI2),0,0,0,-(F1)/EI2),0,0,0,0,0,0; 0,0,-(EI3*cos(lambda*l3)*lambda^3)/(EI4),lambda^3,0,0,(EI3*sin(lambda*l3)*lambda^3)/(EI4),-(F3)/(EI4),0,0,0,-(F3)/(EI4),0,0,0,0; 0,0,0,sin(lambda*l4),0,0,0,cos(lambda*l4),0,0,0,1,0,0,0,l4; -cos(lambda*l1)*lambda,lambda,0,0,sin(lambda l1)*lambda,0,0,0,0,0,0,0,-1,1,0,0; 0,-cos(lambda*l2)*lambda,lambda,0,0,sin(lambda*l2)*lambda,0,0,0,0,0,0,0,-1,1,0; 0,0,-cos(lambda*l3)*lambda,lambda,0,0,sin(lambda*l3)*lambda,0,0,0,0,0,0,0,-1,1; -sin(lambda*l1),0,0,0,-cos(lambda*l1),1,0,0,-1,1,0,0,-l1,0,0,0; 0,-sin(lambda*l2),0,0,0,-cos(lambda*l2),1,0,0,-1,1,0,0,-l2,0,0; 0,0,-sin(lambda*l3),0,0,0,-cos(lambda*l3),1,0,0,-1,1,0,0,-l3,0] |
||
17.09.2010, 16:26 | tigerbine | Auf diesen Beitrag antworten » |
Ohne Konkrete Zahlen nützt mir das in Matlab natürlich nichts... Wollte das zum Vergleich ausrechnen. |
||
Anzeige | ||
|
||
17.09.2010, 17:36 | Formelmonster | Auf diesen Beitrag antworten » |
Nur mal so eine prinzipielle Frage: Musst Du diese Herangehensweise benutzen um die Eigenwerte zu bestimmen? Da Du im wesentlichen Nullstellen eines Polynoms zehnten Grades suchst läuft das ganze auf ein iteratives Verfahren raus. Dieses müsstest Du mindestens zehn mal mit verschiedenen – gut geratenen – Startwerten starten um alle Nullstellen zu finden. In jedem einzelnen Iterationsschritt müsstest Du die Determinante für gegebenes Lambda neu berechnen. Das ist ziemlich teuer! Wenn Du einfach nur einen beliebigen Algorithmus zur Eigenwertberechnung brauchst solltest Du Dich mal mit dem „QR-Verfahren“ auseinandersetzen. |
||
17.09.2010, 18:19 | donmemento | Auf diesen Beitrag antworten » |
@tigerbine: werte sind frei wählbar, ich hatte ei1-ei4 = 2,1 *10^6, die steifigkeiten alle 0, dann ist der EW = pi/(l1+l2+l3+l4), sind die einzelnen Eulerfälle (man kann die drehfedersteifigkeit -> unendl. gehen lassen für den nächsten euler... @formelmonster: ich hab das so programmiert, dass ich mit lambda=0 anfange und die det berechne, anschließend lambda=lambda+1, wenn sich das Vorzeichen der Det ändert weiß ich, dass ich an der Nullstelle vorbei bin, dann gehe ich einen schritt zurück mache mit der nächst kleineren zehnerpotenz weiter. So kann ich mich mit beliebiger genauigkeit dem kleinsten EW nähern (ist auch der einzige den ich brauche) werd mir das QR-Verfahren auf jeden Fall mal ansehen |
||
22.09.2010, 12:55 | donmemento | Auf diesen Beitrag antworten » |
erstmal danke für eure Mithilfe. Bin mittlerweile fest überzeugt, dass der Fehler nicht im Algorithmus, sondern in meiner Matrix ist. |
||
22.09.2010, 13:00 | tigerbine | Auf diesen Beitrag antworten » |
Eine LR-Zerlegung hatte ich hier mal programmiert. Daher nach matlab gefragt. Da du aber Parameter hast, konnte ich das nicht testen. EW Programme habe ich aktuell keine. [WS] - Programmsammlung Numerik Für Gegencheck fehlt mir die Zeit. Da sieht es bei Formelmonster vielleicht besser aus. Gruß und viel Erfolg. |
|
Verwandte Themen
Die Beliebtesten » |
|
Die Größten » |
|
Die Neuesten » |
|