| 20.12.2007, 18:25 |
ProAl |
Auf diesen Beitrag antworten » |
Cholesky Verfahren - Programieren Mit Gauß'sche Normalen Gleichung
hallo leute ich soll ein programm schreiben der mir ein ausgleichsproblem
A*x=b
nach gauß'sch Normalen gleichen mit cholesky erfahren ausgibt.
wobei A ist ein mxn matrix ist, m>=n
gauß'sche Normalen gleichung ist definiert durch
B=A^T*A
und c=A^T*b
nun ich habe das program bis ausgabe von obere/untere dreiecksmatrix sowie ausgabe von B=A^T*A und c=A^T*b geschrieben.
ich muss nun durch vorwärt/rückwärts-seinsetzen x werte bestimmen aber leider komme ich nicht weiter bzw. nicht weiß wie ich vorwärt/rückwärts-seinsetzen programieren soll bzw. kann.
hat eventuel jemand ein idee wie ich dies lösen kann.
mein quell code sieht so aus.
| code: |
1:
2:
3:
4:
5:
6:
7:
8:
9:
10:
11:
12:
13:
14:
15:
16:
17:
18:
19:
20:
21:
22:
23:
24:
25:
26:
27:
28:
29:
30:
31:
32:
33:
34:
35:
36:
37:
38:
39:
40:
|
function gausCholesky(A,b);
B=A'*A; %Gauß'sche Normalengleichung um
%nxn zu Erzegen
C=A'*b;
[n,n] = size(B);
L=0; % Initialization von L
L(1,1)=sqrt(B(1,1));
for k=1:n-1 %Schleife
for i=1:k
sum=0;
for j=1:i-1 %Schleife Für Berechnung Der Summe
sum=sum+L(i,j)*L(k+1,j);
end
L(k+1,i)=1/L(i,i)*(B(k+1,i) - sum); %Nicht Diagonale Elemente von L
end
sum=0; %Intialisierung Der Variable Summe "0"
for j=1:k %Schleife Für Berechnung Der Summe
sum=sum+L(k+1,j)^2;
end
L(k+1,k+1)=sqrt(B(k+1,k+1) - sum); %Diagonale Elemente von L
end
%--------------------------------------------------------------
% Ab Hier Ausgabe Daten
%--------------------------------------------------------------
disp('Ausgabe Der Eigentliche Matrix A')
A
disp('Ausgabe Der Transponierter Eigentliche Matrix A^T')
AT=A'
disp('Ausgabe B=A^T*A Matrix(nxn)')
B
disp('Ausgabe von C Vektor')
C
disp('Ausgabe von Obere Dreiecks Matrix L (nxn)')
L
disp('Ausgabe von Untere Dreiecks Matrix L^T (nxn)')
LT=L' |
|
|
| 20.12.2007, 18:45 |
ProAl |
Auf diesen Beitrag antworten » |
kann etwa so gehen
| code: |
1:
2:
3:
4:
5:
6:
7:
8:
|
x(1)=C(1)/B(1,1);
for i=2:n
x(i)=-B(i,1)*B(1);
for k=2:i-1
x(i)=B(i)-L(i,k)*y(k);
x(i)=(c(i)+B(i))/L(i,i);
end;
end; |
|
theoretisch sollte gehen aber geht trotzdem nicht. |
| 20.12.2007, 18:49 |
tigerbine |
Auf diesen Beitrag antworten » |
Warum ist einmal c und einmal C. Matlab unterscheidet das.
Wie lautet die Fehlermeldung? |
| 20.12.2007, 18:52 |
ProAl |
Auf diesen Beitrag antworten » |
| Zitat: |
Original von tigerbine
Warum ist einmal c und einmal C. Matlab unterscheidet das.
Wie lautet die Fehlermeldung? |
ach danke das habe ich übersehen, fehler meldung gib nicht aus.
und wenn ich c richtig schreibe dann gibt mir folgender fehler meldung
| code: |
1:
2:
|
Error in ==> gausCholesky at 27
x(1)=c(1)/B(1,1); |
|
ok nun habe ich alle c richtig als klein c geschrieben und zwar es waren zwei stück, aber geht trotzdem nicht und immer noch gleich fehler meldung wie oben. |
| 20.12.2007, 19:50 |
ProAl |
Auf diesen Beitrag antworten » |
so ich bin ein stück weiter, aber immer noch nicht am ziel und zwar durch diese for schleife wird mir x werte ausgegeben aber dummer weise als zeile und nicht als spalte.
| code: |
1:
2:
3:
4:
5:
6:
7:
8:
|
for i=1:n
temp=0;
for j=1:i-1
temp = temp + x(j)*L(i,j);
end
x(i)=(c(i)-temp)/L(i,i);
end
x |
|
was ist denn jetz falsch? |
| 20.12.2007, 19:55 |
ProAl |
Auf diesen Beitrag antworten » |
wenn ich bei der ausgabe x transformiere dann passt es damit? aber ist es so richtig? |
| Anzeige | |
|
|
| 20.12.2007, 19:59 |
tigerbine |
Auf diesen Beitrag antworten » |
| Zitat: |
| x werte ausgegeben aber dummer weise als zeile und nicht als spalte. |
Matlab erzeugt doch selbständig die Dimensionen von x, oder? Schreib doch statt x(i) => x(i,1) |
| 20.12.2007, 20:00 |
ProAl |
Auf diesen Beitrag antworten » |
somit sieht mein ganze quell code so aus
| code: |
1:
2:
3:
4:
5:
6:
7:
8:
9:
10:
11:
12:
13:
14:
15:
16:
17:
18:
19:
20:
21:
22:
23:
24:
25:
26:
27:
28:
29:
30:
31:
32:
33:
34:
35:
36:
37:
38:
39:
40:
41:
42:
43:
44:
45:
46:
47:
48:
49:
50:
51:
|
function gausCholesky(A,b);
B=A'*A; %Gauß'sche Normalengleichung um
%nxn zu Erzegen
c=A'*b;
[n,n] = size(B);
L=0; % Initialization von L
L(1,1)=sqrt(B(1,1));
for k=1:n-1 %Schleife
for i=1:k
sum=0;
for j=1:i-1 %Schleife Für Berechnung Der Summe
sum=sum+L(i,j)*L(k+1,j);
end
L(k+1,i)=1/L(i,i)*(B(k+1,i) - sum); %Nicht Diagonale Elemente von L
end
sum=0; %Intialisierung Der Variable Summe "0"
for j=1:k %Schleife Für Berechnung Der Summe
sum=sum+L(k+1,j)^2;
end
L(k+1,k+1)=sqrt(B(k+1,k+1) - sum); %Diagonale Elemente von L
end
for i=1:n
temp=0;
for j=1:i-1
temp = temp + x(j)*L(i,j);
end
x(i,1)=(c(i)-temp)/L(i,i);
end
x
%--------------------------------------------------------------
% Ab Hier Ausgabe Daten
%--------------------------------------------------------------
disp('Ausgabe Der Eigentliche Matrix A')
A
disp('Ausgabe Der Transponierter Eigentliche Matrix A^T')
AT=A'
disp('Ausgabe B=A^T*A Matrix(nxn)')
B
disp('Ausgabe von C Vektor')
c
disp('Ausgabe von Obere Dreiecks Matrix L (nxn)')
L
disp('Ausgabe von Untere Dreiecks Matrix L^T (nxn)')
LT=L'
|
|
also wenn ich x werte vektor bei der ausgabe transformiere pass es und wenn man probe auch macht klappt auch? aber die frage ist ob es so wirklich richtig ist? |
| 20.12.2007, 20:01 |
ProAl |
Auf diesen Beitrag antworten » |
| Zitat: |
Original von tigerbine
| Zitat: |
| x werte ausgegeben aber dummer weise als zeile und nicht als spalte. |
Matlab erzeugt doch selbständig die Dimensionen von x, oder? Schreib doch statt x(i) => x(i,1) |
omg danke vielmals, es stimmt es sollte x(i,1) sein. |
| 20.12.2007, 20:19 |
ProAl |
Auf diesen Beitrag antworten » |
RE: Cholesky Verfahren - Programieren Mit Gauß'sche Normalen Gleichung
ich habe doch noch einer frage und zwar wenn ich probe innerhalb der funktion definiere etwa so
| code: |
1:
2:
3:
4:
5:
6:
7:
|
if
b-(L*x) < 1e-8
disp('Ja, Stimmt!')
else
disp('Nein, Stimmt Nicht:')
end |
|
bekomme ich folgender fehler
| code: |
1:
2:
3:
4:
|
??? Error: File: gausCholesky.m Line: 55 Column: 6
Expression or statement is incomplete or incorrect.
|
|
was könnte ursache dafür sein? |
| 28.12.2007, 01:18 |
ProAl |
Auf diesen Beitrag antworten » |
hier die endgültige programm code, es lauft wie es sein sollte, leider probe habe ich per hand nicht durchgeführt, daher kann ich nichts sagen ob es wirklich so wie es sein sollte.
immer hin funktioniert es.
| code: |
1:
2:
3:
4:
5:
6:
7:
8:
9:
10:
11:
12:
13:
14:
15:
16:
17:
18:
19:
20:
21:
22:
23:
24:
25:
26:
27:
28:
29:
30:
31:
32:
33:
34:
35:
36:
37:
38:
39:
40:
41:
42:
43:
44:
45:
46:
47:
48:
49:
50:
51:
52:
53:
54:
55:
56:
57:
58:
59:
|
function gausCholesky(A,b);
%Benutzungs Anweisung Der Programm
% - Eingabe Einer Matrix A sovie vektor b
% dann in der command windows "gausCholesky(A,b)" eigeben
B=A'*A; %Gauß'sche Normalengleichung um
%nxn Matrix zu Erzegen
c=A'*b; %Vektor
[n,n] = size(B); %Dimension von Matrix
L=0; %Initialization von L
L(1,1)=sqrt(B(1,1));
for k=1:n-1 %Schleife
for i=1:k
sum=0;
for j=1:i-1 %Schleife Für Berechnung Der Summe
sum=sum+L(i,j)*L(k+1,j);
end
L(k+1,i)=1/L(i,i)*(B(k+1,i) - sum); %Nicht Diagonale Elemente von L
end
sum=0; %Intialisierung Der Variable Summe "0"
for j=1:k %Schleife Für Berechnung Der Summe
sum=sum+L(k+1,j)^2;
end
L(k+1,k+1)=sqrt(B(k+1,k+1) - sum); %Diagonale Elemente von L
end
%--------------------------------------------------------------
% Ab Hier Schleife zum Berechnen von X-Werte
%--------------------------------------------------------------
for i=1:n
temp=0;
for j=1:i-1
temp = temp + x(j)*L(i,j);
end
x(i,1)=(c(i)-temp)/L(i,i);
end
%--------------------------------------------------------------
% Ab Hier Ausgabe Daten und Probe
%--------------------------------------------------------------
disp('Ausgabe Der Eigentliche Matrix A')
A
disp('Ausgabe Der Transponierter Eigentliche Matrix A^T')
AT=A'
disp('Ausgabe B=A^T*A Matrix(nxn)')
B
disp('Ausgabe von C Vektor')
c
disp('Ausgabe von Obere Dreiecks Matrix L (nxn)')
L
disp('Ausgabe von Untere Dreiecks Matrix L^T (nxn)')
LT=L'
disp('Ausgabe von x - Werte')
x
|
|
|
