====== std.linalg — dichte lineare Algebra ======
Matrizen in ''f64'': LU, QR, Cholesky, Gleichungssysteme lösen, Determinante, Inverse, Rang, Ausgleichsrechnung und symmetrische Eigenwerte nach Jacobi.
→ [[lyx_-_programmiersprache:units:mathematik|Mathematik-Units]] · [[lyx_-_programmiersprache:units:decomp|std.decomp]] · [[lyx_-_programmiersprache:units:sparse|std.sparse]] · [[lyx_-_programmiersprache:units:vector|std.vector]]
import std.linalg;
Schicht 0, hängt an ''std.math'' und ''std.alloc''. 35 Funktionen. Alle Beispiele mit ''lyxc 1.0.21A'' übersetzt und ausgeführt und gegen eine unabhängige Referenzrechnung geprüft.
----
===== Darstellung =====
Eine ''m×n''-Matrix ist ein Speicherblock aus ''m·n'' ''f64''-Werten, **zeilenweise** abgelegt. Der Zugriff läuft über ''LinGet''/''LinSet'', die die Spaltenzahl als Argument bekommen:
var n: int64 := 3;
var a: int64 := alloc(n * n * 8);
LinSet(a, n, 0, 1, -2.0); // Zeile 0, Spalte 1
var v: f64 := LinGet(a, n, 0, 1);
> **Die Spaltenzahl wird nicht mitgeführt.** ''LinGet(a, cols, i, j)'' rechnet ''a + (i·cols + j)·8'' und prüft nichts. Ein falsches ''cols'' liest stillschweigend an der falschen Stelle — kein Fehler, nur falsche Zahlen. Bei nicht-quadratischen Matrizen ist das die häufigste Fehlerquelle; die Maße gehören neben die Matrix.
----
===== Ein vollständiges Beispiel =====
unit main;
import std.io;
import std.string;
import std.alloc;
import std.linalg;
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 n: int64 := 3;
var a: int64 := alloc(n * n * 8);
LinSet(a, n, 0,0, 4.0); LinSet(a, n, 0,1, -2.0); LinSet(a, n, 0,2, 1.0);
LinSet(a, n, 1,0, -2.0); LinSet(a, n, 1,1, 4.0); LinSet(a, n, 1,2, -2.0);
LinSet(a, n, 2,0, 1.0); LinSet(a, n, 2,1, -2.0); LinSet(a, n, 2,2, 4.0);
ZF("Det : ", LinDet(a, n));
ZF("Trace : ", LinTrace(a, n));
ZI("Rank : ", LinRank(a, n, n));
ZI("IsSymmetric : ", LinIsSymmetric(a, n, LinEpsilon()) as int64);
ZI("IsPosDef : ", LinIsPosDef(a, n) as int64);
ZF("NormFro : ", LinNormFro(a, n, n));
ZF("Cond2 : ", LinCond2(a, n, n));
var b: int64 := alloc(n * 8);
pokef64(b, 1.0); pokef64(b + 8, 2.0); pokef64(b + 16, 3.0);
var x: int64 := alloc(n * 8);
var ok: bool := LinSolve(a, n, b, x);
Print("Solve ok="); Print(IntToStr(ok as int64)); PrintLn("");
ZF("x0 : ", peekf64(x));
ZF("x1 : ", peekf64(x + 8));
ZF("x2 : ", peekf64(x + 16));
var ev: int64 := alloc(n * 8);
var evec: int64 := alloc(n * n * 8);
var ok2: bool := LinEigenSym(a, n, ev, evec);
LinEigenSortDesc(ev, evec, n);
Print("Eigen ok="); Print(IntToStr(ok2 as int64)); PrintLn("");
ZF("lambda0 : ", peekf64(ev));
ZF("lambda1 : ", peekf64(ev + 8));
ZF("lambda2 : ", peekf64(ev + 16));
return 0;
}
Det : 36.000000
Trace : 12.000000
Rank : 3
IsSymmetric : 1
IsPosDef : 1
NormFro : 8.124038
Cond2 : 4.529210
Solve ok=1
x0 : 0.666666
x1 : 1.500000
x2 : 1.333333
Eigen ok=1
lambda0 : 7.372281
lambda1 : 2.999999
lambda2 : 1.627718
Sämtliche Werte stimmen mit der Referenzrechnung überein. ''lambda1'' erscheint als ''2.999999'' statt ''3'': das Jacobi-Verfahren ist **iterativ** und endet an einer Schwelle, nicht exakt. Für Eigenwerte gilt daher immer ein Vergleich mit Toleranz (→ [[lyx_-_programmiersprache:units:feq|std.feq]]), nie ''==''.
----
===== Gleichungssysteme =====
^ Weg ^ Voraussetzung ^ Wann ^
| ''LinSolve(a, n, b, x)'' | quadratisch, regulär | der Normalfall, ein System |
| ''LinLuDecompose'' + ''LinLuSolve'' | dasselbe | **mehrere** rechte Seiten zur selben Matrix |
| ''LinCholesky'' + ''LinCholeskySolve'' | symmetrisch positiv definit | rund doppelt so schnell |
| ''LinLstSq(a, m, n, b, x)'' | überbestimmt (''m > n'') | Ausgleichsrechnung |
**Die Zerlegung einmal, das Lösen mehrfach.** ''LinSolve'' zerlegt bei jedem Aufruf neu; wer dasselbe ''a'' mit zehn verschiedenen ''b'' löst, ruft einmal ''LinLuDecompose'' und danach zehnmal ''LinLuSolve''. ''LinLuDecompose'' schreibt **in die Matrix hinein** und braucht zusätzlich ein Permutationsfeld (''n'' Einträge) sowie eine Adresse für das Vorzeichen — wer ''a'' danach noch braucht, kopiert vorher mit ''LinCopy''.
''LinIsPosDef'' beantwortet vorab, ob der Cholesky-Weg überhaupt gangbar ist.
----
===== Zustand einer Matrix =====
^ Funktion ^ Liefert ^
| ''LinDet(a, n)'' | Determinante über die LU-Zerlegung |
| ''LinRank(a, m, n)'' | Rang |
| ''LinRref(a, m, n)'' | reduzierte Zeilenstufenform, **verändert ''a''** |
| ''LinCond2(a, m, n)'' | Konditionszahl |
| ''LinNormFro'' · ''LinNormOne'' · ''LinNormInf'' | Normen |
| ''LinTrace'' · ''LinIsSymmetric'' · ''LinEquals'' | Spur, Symmetrie, Vergleich mit ''eps'' |
> **''LinCond2'' ist die Frage „darf ich dem Ergebnis trauen?"** Die Konditionszahl sagt, um welchen Faktor sich ein Fehler in den Eingabedaten im Ergebnis verstärkt. Oben ist sie 4,5 — unkritisch. Bei 1e12 bleiben von 16 Stellen noch vier; bei 1e16 ist das Ergebnis wertlos, obwohl ''LinSolve'' brav ''true'' zurückgibt. **Eine erfolgreiche Rückgabe ist keine Aussage über die Genauigkeit.**
''LinEquals'' nimmt ein ''eps'', weil ein Vergleich zweier gerechneter Matrizen mit ''=='' praktisch immer scheitert. ''LinEpsilon()'' liefert den Vorgabewert 1e−12.
----
===== Zerlegungen und Eigenwerte =====
^ Funktion ^ Ergebnis ^
| ''LinLuDecompose'' | LU mit Zeilentausch |
| ''LinQr(a, m, n, q, r)'' | orthogonales ''Q'', obere Dreiecksmatrix ''R'' |
| ''LinCholesky(a, n, out)'' | ''L'' mit ''L·Lᵀ = A'' |
| ''LinEigenSym(a, n, evals, evecs)'' | Eigenwerte **und** -vektoren, nur symmetrisch |
| ''LinEigenSortDesc'' | sortiert beide gemeinsam absteigend |
| ''LinSingularValues(a, m, n, out)'' | Singulärwerte |
''LinEigenSym'' setzt eine **symmetrische** Matrix voraus — das Jacobi-Verfahren beruht darauf. Für allgemeine Matrizen und für die vollständige SVD ist [[lyx_-_programmiersprache:units:decomp|std.decomp]] zuständig.
''LinEigenSortDesc'' sortiert Werte und Vektoren **gemeinsam**; wer nur die Werte sortiert, ordnet die Vektoren falsch zu.
----
===== Rechnen =====
''LinAdd'', ''LinSub'', ''LinScale'', ''LinMul'', ''LinMulVec'', ''LinTranspose'' — dazu ''LinFill'', ''LinZero'', ''LinIdentity'' und ''LinCopy'' zum Anlegen.
''LinMul(a, b, m, k, n, out)'' multipliziert eine ''m×k''- mit einer ''k×n''-Matrix. Die drei Maße stehen einzeln in der Signatur, weil sich das aus den Zeigern nicht ableiten lässt. **''out'' darf nicht auf ''a'' oder ''b'' zeigen.**
----
===== Abgrenzung =====
^ Aufgabe ^ Unit ^
| Dichte Matrizen, Lösen und Zerlegen | ''std.linalg'' |
| SVD, Pseudoinverse, allgemeine Eigenwerte | [[lyx_-_programmiersprache:units:decomp|std.decomp]] |
| Dünnbesetzte **Vektoren** | [[lyx_-_programmiersprache:units:sparse|std.sparse]] |
| Kleine feste Vektoren, Festkomma | [[lyx_-_programmiersprache:units:vector|std.vector]] |
Dünnbesetzte **Matrizen** und iterative Löser (CG, GMRES) gibt es derzeit nicht. ''std.matrix'' aus dem Bestand existiert, hat aber keine eigene Wiki-Seite; für Rechnungen in ''f64'' ist ''std.linalg'' der Weg.
----
Letzte Aktualisierung: 2026-08-16 · alle Beispiele mit ''lyxc 1.0.21A'' übersetzt und ausgeführt; Determinante, Norm, Kondition, Lösung und Eigenwerte gegen eine unabhängige Referenzrechnung geprüft.