Anonymous user
LU decomposition: Difference between revisions
Updated D code
m (Wrong sizeof()) |
(Updated D code) |
||
Line 585:
=={{header|D}}==
{{trans|Common Lisp}}
<lang d>import std.stdio, std.algorithm, std.typecons
std.array, std.conv;
bool isRectangular(T)(in T[][] M) pure nothrow {
//return !canFind!((r){ return r.length != M[0].length; })(M);
foreach (row; M)
if (row.length != M[0].length)
return false;
return true;
}
string prettyPrint(T)(in T[][] M) {
return "[" ~ array(map!text(M)).join(",\n ") ~ "]";
}
T[][] matrixMul(T)(in T[][] A, in T[][] B) pure nothrow
in {
assert(A.length && B.length &&
isRectangular(A) && isRectangular(B) &&
A[0].length == B.length);
} body {
auto result = new T[][](A.length, B[0].length);
auto aux = new T[B.length];
foreach (j; 0 .. B[0].length) {
foreach (k, row; B)
aux[k] = row[j];
foreach (i, ai; A)
result[i][j] = dotProduct(ai, aux);
}
return result;
}
/// Creates a nxn identity matrix.
in {
assert(n >= 0);
Line 603 ⟶ 632:
/// Creates the pivoting matrix for m.
in {
// is square
assert(isRectangular(m) && m[0].length == m.length);
} body {
immutable int n = m.length;
auto P = identityMatrix!T(n);
foreach (
T max = m[
int row =
foreach (
if (m[
max = m[
row =
}
if (
swap(P[
}
Line 627 ⟶ 656:
/// Decomposes a square matrix A by PA=LU and returns L, U and P.
lu(T)(in T[][] A) pure nothrow
in {
// is square
assert(isRectangular(A) && A[0].length == A.length);
} body {
immutable int n = A.length;
auto L = new T[][](n, n);
auto U = new T[][](n, n);
Line 639 ⟶ 669:
U[i][0 .. i] = 0;
}
const auto P = pivotize!T(A);
auto A2 = matrixMul!T(P, A);
Line 658 ⟶ 688:
}
return
}
void main() {
enum double[][] a = [[1, 3, 5],
[2, 4, 7],
[1, 1, 0]];
foreach (part; lu(a).tupleof)
writeln(prettyPrint(part), "\n");
writeln();
enum double[][] b = [[11, 9, 24, 2
[
[3, 17, 18, 1],
[2, 5, 7, 1]];
foreach (part; lu(b).tupleof)
writeln(prettyPrint(part), "\n");
| |||