====== std.decomp — Matrixzerlegungen ======
Echte SVD nach einseitigem Jacobi, Pseudoinverse, Rang-k-Näherung, polare Zerlegung und allgemeine Eigenwerte über QR-Iteration.
→ [[lyx_-_programmiersprache:units:mathematik|Mathematik-Units]] · [[lyx_-_programmiersprache:units:linalg|std.linalg]] · [[lyx_-_programmiersprache:units:sparse|std.sparse]]
import std.decomp;
Schicht 1, hängt an [[lyx_-_programmiersprache:units:linalg|std.linalg]], ''std.math'' und ''std.alloc''. 14 Funktionen — die kleinste Schnittstelle der Sammlung bei 913 Zeilen Umsetzung. Alle Beispiele mit ''lyxc 1.0.21A'' übersetzt und ausgeführt.
----
===== Verhältnis zu ''std.linalg'' =====
''std.decomp'' ist die **Fortsetzung** von ''std.linalg'', nicht dessen Ersatz. Die Trennung verläuft entlang „Lösen" gegen „Zerlegen":
^ Aufgabe ^ Unit ^
| Quadratisches System lösen, Determinante, Inverse | [[lyx_-_programmiersprache:units:linalg|std.linalg]] |
| Symmetrische Eigenwerte (Jacobi) | ''std.linalg'' |
| **SVD**, Pseudoinverse, Rang-k-Näherung | ''std.decomp'' |
| **Allgemeine** (nicht symmetrische) Eigenwerte | ''std.decomp'' |
Beides zusammen wären über 2000 Zeilen in einer Unit. Getrennt bleibt die Aufteilung verständlich, und wer nur Gleichungssysteme löst, zieht die schwerere Unit nicht mit.
----
===== SVD und was daraus folgt =====
unit main;
import std.io;
import std.string;
import std.alloc;
import std.linalg;
import std.decomp;
fn ZF(t: pchar, v: f64): void { Print(t); PrintF64(v); }
fn ZI(t: pchar, v: int64): void { Print(t); PrintLn(IntToStr(v)); }
fn main(): int64 {
var m: int64 := 3; var n: int64 := 2;
var a: int64 := alloc(m * n * 8);
LinSet(a, n, 0,0, 1.0); LinSet(a, n, 0,1, 2.0);
LinSet(a, n, 1,0, 3.0); LinSet(a, n, 1,1, 4.0);
LinSet(a, n, 2,0, 5.0); LinSet(a, n, 2,1, 6.0);
var s: int64 := alloc(n * 8);
DecompSvdValues(a, m, n, s);
ZF("Singulaerwert 0 : ", peekf64(s));
ZF("Singulaerwert 1 : ", peekf64(s + 8));
ZF("Norm2 : ", DecompNorm2(a, m, n));
ZF("Cond2 : ", DecompCond2(a, m, n));
ZI("Rank : ", DecompRank(a, m, n, DecompEpsilon()));
var b: int64 := alloc(m * 8);
pokef64(b, 1.0); pokef64(b + 8, 2.0); pokef64(b + 16, 3.0);
var x: int64 := alloc(n * 8);
var ok: bool := DecompLstSq(a, m, n, b, DecompEpsilon(), x);
Print("LstSq ok="); Print(IntToStr(ok as int64)); PrintLn("");
ZF("x0 : ", peekf64(x));
ZF("x1 : ", peekf64(x + 8));
var c: int64 := alloc(4 * 8);
LinSet(c, 2, 0,0, 4.0); LinSet(c, 2, 0,1, 1.0);
LinSet(c, 2, 1,0, 2.0); LinSet(c, 2, 1,1, 3.0);
ZF("Spektralradius : ", DecompSpectralRadius(c, 2));
var evr: int64 := alloc(2 * 8);
var evi: int64 := alloc(2 * 8);
DecompEigenGeneral(c, 2, evr, evi);
ZF("Eigenwert 0 Re : ", peekf64(evr));
ZF("Eigenwert 1 Re : ", peekf64(evr + 8));
return 0;
}
Singulaerwert 0 : 9.525518
Singulaerwert 1 : 0.514300
Norm2 : 9.525518
Cond2 : 18.521305
Rank : 2
LstSq ok=1
x0 : 0.000000
x1 : 0.499999
Spektralradius : 5.000000
Eigenwert 0 Re : 5.000000
Eigenwert 1 Re : 2.000000
Sämtliche Werte stimmen mit der Referenzrechnung überein.
**Drei Größen fallen mit den Singulärwerten ab, ohne dass man sie einzeln rechnen müsste:**
* ''DecompNorm2'' ist der **größte** Singulärwert (9,525518).
* ''DecompCond2'' ist das **Verhältnis** von größtem zu kleinstem (9,525518 / 0,514300 = 18,52).
* ''DecompRank'' ist die **Anzahl** der Singulärwerte oberhalb der Toleranz.
----
===== Rang ist eine Frage der Toleranz =====
**''DecompRank'' braucht ein ''tol'', und die Wahl ist eine Entscheidung.** In exakter Arithmetik ist der Rang eindeutig; in ''f64'' ist ein Singulärwert von ''1e−16'' nicht von null zu unterscheiden. ''DecompEpsilon()'' liefert ''1e−15'' als Vorgabe.
Bei gemessenen Daten ist die richtige Schwelle **größer**: dort trennt sie Signal von Rauschen, und ein Singulärwert unterhalb der Messgenauigkeit ist keine echte Dimension. Wer die Vorgabe nimmt, bekommt bei verrauschten Daten regelmäßig vollen Rang — mathematisch korrekt, sachlich nutzlos.
----
===== Pseudoinverse und Ausgleichsrechnung =====
^ Funktion ^ Zweck ^
| ''DecompPseudoInverse(a, m, n, tol, out)'' | Moore-Penrose-Inverse |
| ''DecompLstSq(a, m, n, b, tol, x)'' | überbestimmtes System, kleinste Quadrate |
| ''DecompLowRank(a, m, n, k, out)'' | beste Näherung vom Rang ''k'' |
Die Pseudoinverse gibt es auch für **nicht quadratische** und **singuläre** Matrizen, wo ''LinInverse'' scheitert. Für ein überbestimmtes System liefert ''DecompLstSq'' die Lösung mit dem kleinsten Fehlerquadrat — oben ''x = (0, 0.5)'' für drei Gleichungen mit zwei Unbekannten.
> **''out'' darf nicht auf ''a'' zeigen.** ''DecompPseudoInverse'' und ''DecompPolar'' prüfen das und geben ''false'' zurück; verlassen sollte man sich darauf nicht, sondern getrennte Puffer verwenden. ''DecompPseudoInverse'' setzt außerdem ''m >= n'' voraus.
''DecompLowRank'' ist die Grundlage von Kompression und Entrauschen über die SVD: die ''k'' größten Singulärwerte tragen den Hauptteil der Information, der Rest ist oft Rauschen.
----
===== Eigenwerte allgemeiner Matrizen =====
^ Funktion ^ Zweck ^
| ''DecompHessenberg(a, n, h, q)'' | Vorbereitung für die QR-Iteration |
| ''DecompEigenGeneral(a, n, evRe, evIm)'' | Eigenwerte, **komplex** |
| ''DecompEigenVector(a, n, lambda, vec)'' | Eigenvektor zu einem bekannten Eigenwert |
| ''DecompSpectralRadius(a, n)'' | größter Betrag aller Eigenwerte |
| ''DecompPolar(a, n, q, p)'' | polare Zerlegung ''A = Q·P'' |
> **Die Eigenwerte kommen in zwei Feldern — Real- und Imaginärteil.** Eine nicht symmetrische reelle Matrix kann komplexe Eigenwerte haben, und zwar immer paarweise konjugiert. Wer nur ''evRe'' liest, bekommt bei einer Drehmatrix zwei gleiche Zahlen und übersieht, dass die eigentliche Information im Imaginärteil steckt.
>
> Für **symmetrische** Matrizen ist ''LinEigenSym'' aus [[lyx_-_programmiersprache:units:linalg|std.linalg]] die bessere Wahl: dort sind alle Eigenwerte reell, das Jacobi-Verfahren ist stabiler und liefert die Eigenvektoren gleich mit.
''DecompSpectralRadius'' beantwortet die Frage, ob eine Iteration ''x ← A·x'' konvergiert: unter 1 ja, darüber nicht.
----
Letzte Aktualisierung: 2026-08-16 · alle Beispiele mit ''lyxc 1.0.21A'' übersetzt und ausgeführt; Singulärwerte, Kondition, Rang, Ausgleichslösung und Eigenwerte gegen eine unabhängige Referenzrechnung geprüft.