Anonymous user
LU decomposition: Difference between revisions
→{{header|Perl 6}}
(Added Sidef) |
|||
Line 2,164:
lu_backsub(lup, transpose([1, 1, -1, -1]));
/* matrix([-204], [2100], [-4740], [2940]) */</lang>
=={{header|PARI/GP}}==
<lang parigp>matlup(M) =
{
my (P = matid(#M), L = P, U = M);
for (i = 1, #M-1, \\ pivoting
p = M[z=i,i];
for (k = i, #M, if (M[k,i] > p, p = M[z=k,i]));
if (i != z, \\ swap rows
k = U[i,]; U[i,] = U[z,]; U[z,] = k;
k = P[i,]; P[i,] = P[z,]; P[z,] = k;
);
);
for (i = 1, #M-1, \\ decompose
for (k = i+1, #M,
L[k,i] = U[k,i] / U[i,i];
for (j = i, #M, U[k,j] -= L[k,i] * U[i,j])
)
);
[L,U,P] \\ return L,U,P triple matrix
}</lang>
Output:
<pre>
gp > [L,U,P] = matlup([1,3,5;2,4,7;1,1,0]);
gp > L
[ 1 0 0]
[1/2 1 0]
[1/2 -1 1]
gp > U
[2 4 7]
[0 1 3/2]
[0 0 -2]
gp > P
[0 1 0]
[1 0 0]
[0 0 1]
gp > [L,U,P] = matlup([11,9,24,2;1,5,2,6;3,17,18,1;2,5,7,1]);
gp > L
[ 1 0 0 0]
[3/11 1 0 0]
[1/11 23/80 1 0]
[2/11 37/160 1/278 1]
gp > U
[11 9 24 2]
[ 0 160/11 126/11 5/11]
[ 0 0 -139/40 91/16]
[ 0 0 0 71/139]
gp > P
[1 0 0 0]
[0 0 1 0]
[0 1 0 0]
[0 0 0 1]
</pre>
=={{header|Perl 6}}==
| |||