LR-Zerlegung < Matlab < Mathe-Software < Mathe < Vorhilfe
|
Status: |
(Frage) reagiert/warte auf Reaktion | Datum: | 17:22 Mo 29.10.2012 | Autor: | Mija |
Aufgabe | Wir sollen die LR-Zerlegung implementieren. |
Hallo,
ich habe die LR-Zerlegung bereits implementiert. Hier ist mal mein Code:
function [L,R,p] = LRZerl(A)
% Das ist die LR-Zerlegung der Matrix A.
n=size(A,1); % Anzahl Zeilen
p=[1:n]';
for j=1:n-1
m=j;
%Pivotisierung
for i=j:n % i-te Zeile
if abs(A(p(i),j)>abs(A(p(j),j)))
m=i;
absajj=abs(A(p(i),j));
end
end
%vertausche Zeile j und Zeile m in A
pm=p(m);
p(m)=p(j);
p(j)=pm;
%Elimination
for i=j+1:n % i-te Zeile
aij=A(p(i),j)/A(p(j),j);
for k=j+1:n
A(p(i),k)=A(p(i),k)-aij*A(p(j),k);
end
A(p(i),j)=aij;
end
end
R=triu(A(p,:));
L=eye(n)+tril(A(p,:),-1);
end
Nun soll ich allerdings die p's innerhalb der A's wegmachen und stattdessen Zeilenvertauschungen noch einfügen.
Ich habe schonmal die p's weggemacht, allerdings stecke ich fest, wie ich die Zeilenvertauschungen hinbekomme, so dass auch alles wieder ordentlich funktioniert.
Hier mein Code ohne p's in den A's:
function [L,R,p] = LRZerl(A)
% Das ist die LR-Zerlegung der Matrix A.
n=size(A,1); % Anzahl Zeilen
p=[1:n]';
for j=1:n-1
m=j;
%Pivotisierung
for i=j:n % i-te Zeile
if abs(A(i,j)>max(abs(A(j,j))))
m=i;
absajj=abs(A(i,j));
end
end
%vertausche Zeile j und Zeile m in A
pm=p(m);
p(m)=p(j);
p(j)=pm;
%Elimination
for i=j+1:n % i-te Zeile
aij=A(i,j)/A(j,j);
for k=j+1:n
A(i,k)=A(i,k)-aij*A(j,k);
end
A(i,j)=aij;
end
end
R=triu(A(p,:));
L=eye(n)+tril(A(p,:),-1);
end
Ich würde mich sehr freuen, wenn mir jemand weiterhelfen könnte! :)
|
|
|
|
Status: |
(Frage) beantwortet | Datum: | 18:50 Mo 29.10.2012 | Autor: | Mija |
Ich habe es bis jetzt folgendermaßen gelöst:
function [L,R,p] = LRZerl(A)
% Das ist die LR-Zerlegung der Matrix A.
n=size(A,1); % Anzahl Zeilen
p=[1:n]';
for j=1:n-1
m=j;
%Pivotisierung
for i=j:n % i-te Zeile
if abs(A(i,j))==max(abs(A(i,j)))
m=i;
%absajj=abs(A(i,j));
end
end
for i=j:n
pm=A(m,i);
A(m,i)=A(j,i);
A(j,i)=pm;
end
%vertausche Zeile j und Zeile m in A
pm=p(m);
p(m)=p(j);
p(j)=pm;
%Elimination
for i=j+1:n % i-te Zeile
aij=A(i,j)/A(j,j);
for k=j+1:n
A(i,k)=A(i,k)-aij*A(j,k);
end
A(i,j)=aij;
end
end
R=triu(A);
L=eye(n,n)+tril(A,-1);
end
Allerdings stimmt die Reihenfolge der Einträge in der ersten Spalte bei mir nicht. Ich habe folgendes Beispiel genommen:
>> A=[2 4 1; 1 2 4; 4 1 2];
>> L=LRZerl(A)
R =
4.0000 1.0000 2.0000
0 3.5000 0
0 0 3.5000
L =
1.0000 0 0
0.2500 1.0000 0
0.5000 0.5000 1.0000
Eigentlich müsste hier
L =
1.0000 0 0
0.5000 1.0000 0
0.2500 0.5000 1.0000
rauskommen.
Wo ist mein Fehler?
|
|
|
|
|
Status: |
(Mitteilung) Reaktion unnötig | Datum: | 10:27 Di 30.10.2012 | Autor: | Mija |
Kann mir denn wirklich niemand weiterhelfen?
|
|
|
|
|
hi,
ich zerpflück einfach mal deinen Programmcode.
Die Zeile
abs(A(i,j))==max(abs(A(i,j)))
macht nichts weiteres als, einen Wert mit sich selber vergleichen. Du pivotisierst immer die letzte Zeile nach oben.
1: | for i=j:n
| 2: | pm=A(m,i);
| 3: | A(m,i)=A(j,i);
| 4: | A(j,i)=pm;
| 5: | end
|
Hier kopierst du wirklich jeden Eintrag? Man kann komplette Zeilen kopieren oder gar tauschen.
Ehrlicherweise weiß ich auch nicht genau, wo in dem Teil der Elimination der Fehler steckt. Ich habe manch einmal auch große Problem die Zugriffe auf das Array nachzuvollziehen.
Zum Glück muss man nicht mit den ganzen Indizies herumhandtieren.
Ich habe deinen Code etwas abgeändert, sodass er läuft:
http://pastebin.de/30566
gruß
wieschoo
|
|
|
|
|
Status: |
(Frage) beantwortet | Datum: | 15:55 Sa 23.02.2013 | Autor: | henrymrl |
Aufgabe | Erstellen Sie eine Prozedur zur LU-Zerlegung bei fester Mantissenlänge, ohne bzw.
mit Spaltenpivotstrategie. |
Hallo,
ich habe die LR-Zerlegung ebenfalls selbst programmiert und mir ein paar Ideen aus den oben beschriebenen Quelltexten geholt (da ich Matlab erst seid einer Woche benutze).
Ich habe ebenfalls ein Programm geschrieben, mit dem die Mantissenlänge einer Float-Zahl "beschnitten" werden kann.
Das Problem ist, diese beiden Programme zusammenzuführen.
Meine Vermutung ist, das Matlab die Operation zeilenweise ausführt, mein Rundungs-Programm aber nur einzelne Zahlen verarbeitet.
Das ist der Quelltext für die Rundungs-Routine:
function [y]=runden(x,t)
k=0;
while floor(x)==0
x=x*10;
k=k-1;
end
x=x/10;
k=k+1;
if x>0
while (x>1)
x=x/10;
k=k+1;
end
[mm] x=x*(10^t);
[/mm]
x=round(x);
[mm] y=x/(10^t);
[/mm]
elseif x<0
x=x*(-1);
while (x>1)
x=x/10;
k=k+1;
end
[mm] x=x*(10^t);
[/mm]
x=round(x);
[mm] y=(x/(10^t))*(-1);
[/mm]
end
[mm] y=y*(10^k);
[/mm]
end
und so wird die Routine aufgerufen:
U(k+1:n,k)=runden((U(k+1:n,k)/U(k,k)),t);
U(k+1:n,k+1:n)=runden((U(k+1:n,k+1:n)-runden((U(k+1:n,k)*U(k,k+1:n)),t)),t);
Für jede Hilfe bedanke ich mich bereits im Voraus.
|
|
|
|
|
Moin,
dann ändere doch deine Funktion zum Runden, indem du alles, was du bisher gemacht hast in eine Schleife
for i=1:length(x)
% bearbeite x(i)
end
steckst und x zurück lieferst.
|
|
|
|