LU decomposition: Difference between revisions
Content added Content deleted
SqrtNegInf (talk | contribs) m (→{{header|Raku}}: oops, fix variable name) |
Thundergnat (talk | contribs) m (syntax highlighting fixup automation) |
||
Line 194: | Line 194: | ||
{{trans|Python}} |
{{trans|Python}} |
||
< |
<syntaxhighlight lang="11l">F pprint(m) |
||
L(row) m |
L(row) m |
||
print(row) |
print(row) |
||
Line 243: | Line 243: | ||
L(part) lu(b) |
L(part) lu(b) |
||
pprint(part) |
pprint(part) |
||
print()</ |
print()</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 279: | Line 279: | ||
{{works with|Ada 2005}} |
{{works with|Ada 2005}} |
||
decomposition.ads: |
decomposition.ads: |
||
< |
<syntaxhighlight lang="ada">with Ada.Numerics.Generic_Real_Arrays; |
||
generic |
generic |
||
with package Matrix is new Ada.Numerics.Generic_Real_Arrays (<>); |
with package Matrix is new Ada.Numerics.Generic_Real_Arrays (<>); |
||
Line 287: | Line 287: | ||
procedure Decompose (A : Matrix.Real_Matrix; P, L, U : out Matrix.Real_Matrix); |
procedure Decompose (A : Matrix.Real_Matrix; P, L, U : out Matrix.Real_Matrix); |
||
end Decomposition;</ |
end Decomposition;</syntaxhighlight> |
||
decomposition.adb: |
decomposition.adb: |
||
< |
<syntaxhighlight lang="ada">package body Decomposition is |
||
procedure Swap_Rows (M : in out Matrix.Real_Matrix; From, To : Natural) is |
procedure Swap_Rows (M : in out Matrix.Real_Matrix; From, To : Natural) is |
||
Line 366: | Line 366: | ||
end Decompose; |
end Decompose; |
||
end Decomposition;</ |
end Decomposition;</syntaxhighlight> |
||
Example usage: |
Example usage: |
||
< |
<syntaxhighlight lang="ada">with Ada.Numerics.Real_Arrays; |
||
with Ada.Text_IO; |
with Ada.Text_IO; |
||
with Decomposition; |
with Decomposition; |
||
Line 421: | Line 421: | ||
Ada.Text_IO.Put_Line ("U:"); Print (U_2); |
Ada.Text_IO.Put_Line ("U:"); Print (U_2); |
||
Ada.Text_IO.Put_Line ("P:"); Print (P_2); |
Ada.Text_IO.Put_Line ("P:"); Print (P_2); |
||
end Decompose_Example;</ |
end Decompose_Example;</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 465: | Line 465: | ||
=={{header|AutoHotkey}}== |
=={{header|AutoHotkey}}== |
||
< |
<syntaxhighlight lang="autohotkey">;-------------------------- |
||
LU_decomposition(A){ |
LU_decomposition(A){ |
||
P := Pivot(A) |
P := Pivot(A) |
||
Line 536: | Line 536: | ||
return "[" Trim(output, "`n,") "]" |
return "[" Trim(output, "`n,") "]" |
||
} |
} |
||
;--------------------------</ |
;--------------------------</syntaxhighlight> |
||
Examples:< |
Examples:<syntaxhighlight lang="autohotkey">A1 := [[1, 3, 5] |
||
, [2, 4, 7] |
, [2, 4, 7] |
||
, [1, 1, 0]] |
, [1, 1, 0]] |
||
Line 556: | Line 556: | ||
} |
} |
||
MsgBox, 262144, , % result |
MsgBox, 262144, , % result |
||
return</ |
return</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>A:= |
<pre>A:= |
||
Line 604: | Line 604: | ||
=={{header|BBC BASIC}}== |
=={{header|BBC BASIC}}== |
||
< |
<syntaxhighlight lang="bbcbasic"> DIM A1(2,2) |
||
A1() = 1, 3, 5, 2, 4, 7, 1, 1, 0 |
A1() = 1, 3, 5, 2, 4, 7, 1, 1, 0 |
||
PROCLUdecomposition(A1(), L1(), U1(), P1()) |
PROCLUdecomposition(A1(), L1(), U1(), P1()) |
||
Line 663: | Line 663: | ||
a$ = LEFT$(LEFT$(a$)) + CHR$(13) + CHR$(10) |
a$ = LEFT$(LEFT$(a$)) + CHR$(13) + CHR$(10) |
||
NEXT i% |
NEXT i% |
||
= a$</ |
= a$</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 701: | Line 701: | ||
=={{header|C}}== |
=={{header|C}}== |
||
Compiled with <code>gcc -std=gnu99 -Wall -lm -pedantic</code>. Demonstrating how to do LU decomposition, and how (not) to use macros. < |
Compiled with <code>gcc -std=gnu99 -Wall -lm -pedantic</code>. Demonstrating how to do LU decomposition, and how (not) to use macros. <syntaxhighlight lang="c">#include <stdio.h> |
||
#include <stdlib.h> |
#include <stdlib.h> |
||
#include <math.h> |
#include <math.h> |
||
Line 829: | Line 829: | ||
return 0; |
return 0; |
||
}</ |
}</syntaxhighlight> |
||
=={{header|C++}}== |
=={{header|C++}}== |
||
< |
<syntaxhighlight lang="cpp">#include <cassert> |
||
#include <cmath> |
#include <cmath> |
||
#include <iomanip> |
#include <iomanip> |
||
Line 1,015: | Line 1,015: | ||
return 0; |
return 0; |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 1,098: | Line 1,098: | ||
Uses the routine (mmul A B) from [[Matrix multiplication]]. |
Uses the routine (mmul A B) from [[Matrix multiplication]]. |
||
< |
<syntaxhighlight lang="lisp">;; Creates a nxn identity matrix. |
||
(defun eye (n) |
(defun eye (n) |
||
(let ((I (make-array `(,n ,n) :initial-element 0))) |
(let ((I (make-array `(,n ,n) :initial-element 0))) |
||
Line 1,156: | Line 1,156: | ||
;; Return L, U and P. |
;; Return L, U and P. |
||
(values L U P)))</ |
(values L U P)))</syntaxhighlight> |
||
Example 1: |
Example 1: |
||
< |
<syntaxhighlight lang="lisp">(setf g (make-array '(3 3) :initial-contents '((1 3 5) (2 4 7)(1 1 0)))) |
||
#2A((1 3 5) (2 4 7) (1 1 0)) |
#2A((1 3 5) (2 4 7) (1 1 0)) |
||
Line 1,166: | Line 1,166: | ||
#2A((1 0 0) (1/2 1 0) (1/2 -1 1)) |
#2A((1 0 0) (1/2 1 0) (1/2 -1 1)) |
||
#2A((2 4 7) (0 1 3/2) (0 0 -2)) |
#2A((2 4 7) (0 1 3/2) (0 0 -2)) |
||
#2A((0 1 0) (1 0 0) (0 0 1))</ |
#2A((0 1 0) (1 0 0) (0 0 1))</syntaxhighlight> |
||
Example 2: |
Example 2: |
||
< |
<syntaxhighlight lang="lisp">(setf h (make-array '(4 4) :initial-contents '((11 9 24 2)(1 5 2 6)(3 17 18 1)(2 5 7 1)))) |
||
#2A((11 9 24 2) (1 5 2 6) (3 17 18 1) (2 5 7 1)) |
#2A((11 9 24 2) (1 5 2 6) (3 17 18 1) (2 5 7 1)) |
||
Line 1,176: | Line 1,176: | ||
#2A((1 0 0 0) (3/11 1 0 0) (1/11 23/80 1 0) (2/11 37/160 1/278 1)) |
#2A((1 0 0 0) (3/11 1 0 0) (1/11 23/80 1 0) (2/11 37/160 1/278 1)) |
||
#2A((11 9 24 2) (0 160/11 126/11 5/11) (0 0 -139/40 91/16) (0 0 0 71/139)) |
#2A((11 9 24 2) (0 160/11 126/11 5/11) (0 0 -139/40 91/16) (0 0 0 71/139)) |
||
#2A((1 0 0 0) (0 0 1 0) (0 1 0 0) (0 0 0 1))</ |
#2A((1 0 0 0) (0 0 1 0) (0 1 0 0) (0 0 0 1))</syntaxhighlight> |
||
=={{header|D}}== |
=={{header|D}}== |
||
{{trans|Common Lisp}} |
{{trans|Common Lisp}} |
||
< |
<syntaxhighlight lang="d">import std.stdio, std.algorithm, std.typecons, std.numeric, |
||
std.array, std.conv, std.string, std.range; |
std.array, std.conv, std.string, std.range; |
||
Line 1,284: | Line 1,284: | ||
foreach (immutable m; [a, b]) |
foreach (immutable m; [a, b]) |
||
writefln(f, lu(m).tupleof); |
writefln(f, lu(m).tupleof); |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>[[1.0, 0.0, 0.0], |
<pre>[[1.0, 0.0, 0.0], |
||
Line 1,315: | Line 1,315: | ||
=={{header|EchoLisp}}== |
=={{header|EchoLisp}}== |
||
< |
<syntaxhighlight lang="scheme"> |
||
(lib 'matrix) ;; the matrix library provides LU-decomposition |
(lib 'matrix) ;; the matrix library provides LU-decomposition |
||
(decimals 5) |
(decimals 5) |
||
Line 1,354: | Line 1,354: | ||
0 0 -3.475 5.6875 |
0 0 -3.475 5.6875 |
||
0 0 0 0.51079 |
0 0 0 0.51079 |
||
</syntaxhighlight> |
|||
</lang> |
|||
=={{header|Fortran}}== |
=={{header|Fortran}}== |
||
< |
<syntaxhighlight lang="fortran">program lu1 |
||
implicit none |
implicit none |
||
call check( reshape([real(8)::1,2,1,3,4,1,5,7,0 ],[3,3]) ) |
call check( reshape([real(8)::1,2,1,3,4,1,5,7,0 ],[3,3]) ) |
||
Line 1,425: | Line 1,425: | ||
end subroutine |
end subroutine |
||
end program</ |
end program</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 1,473: | Line 1,473: | ||
===2D representation=== |
===2D representation=== |
||
{{trans|Common Lisp}} |
{{trans|Common Lisp}} |
||
< |
<syntaxhighlight lang="go">package main |
||
import "fmt" |
import "fmt" |
||
Line 1,582: | Line 1,582: | ||
u.print("u") |
u.print("u") |
||
p.print("p") |
p.print("p") |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 1,624: | Line 1,624: | ||
</pre> |
</pre> |
||
===Flat representation=== |
===Flat representation=== |
||
< |
<syntaxhighlight lang="go">package main |
||
import "fmt" |
import "fmt" |
||
Line 1,746: | Line 1,746: | ||
u.print("u") |
u.print("u") |
||
p.print("p") |
p.print("p") |
||
}</ |
}</syntaxhighlight> |
||
Output is same as from 2D solution. |
Output is same as from 2D solution. |
||
===Library gonum/mat=== |
===Library gonum/mat=== |
||
< |
<syntaxhighlight lang="go">package main |
||
import ( |
import ( |
||
Line 1,782: | Line 1,782: | ||
fmt.Printf("u: %.5f\n\n", mat.Formatted(u, mat.Prefix(" "))) |
fmt.Printf("u: %.5f\n\n", mat.Formatted(u, mat.Prefix(" "))) |
||
fmt.Println("p:", lu.Pivot(nil)) |
fmt.Println("p:", lu.Pivot(nil)) |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Pivot format is a little different here. (But library solutions don't really meet task requirements anyway.) |
Pivot format is a little different here. (But library solutions don't really meet task requirements anyway.) |
||
Line 1,819: | Line 1,819: | ||
===Library go.matrix=== |
===Library go.matrix=== |
||
< |
<syntaxhighlight lang="go">package main |
||
import ( |
import ( |
||
Line 1,845: | Line 1,845: | ||
fmt.Printf("u:\n%v\n", u) |
fmt.Printf("u:\n%v\n", u) |
||
fmt.Printf("p:\n%v\n", p) |
fmt.Printf("p:\n%v\n", p) |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 1,889: | Line 1,889: | ||
=={{header|Haskell}}== |
=={{header|Haskell}}== |
||
''Without elem-at-index modifications; doesn't find maximum but any non-zero element'' |
''Without elem-at-index modifications; doesn't find maximum but any non-zero element'' |
||
<syntaxhighlight lang="haskell"> |
|||
<lang Haskell> |
|||
import Data.List |
import Data.List |
||
import Data.Maybe |
import Data.Maybe |
||
Line 1,959: | Line 1,959: | ||
main = putStrLn "Task1: \n" >> solveTask mY1 >> |
main = putStrLn "Task1: \n" >> solveTask mY1 >> |
||
putStrLn "Task2: \n" >> solveTask mY2 |
putStrLn "Task2: \n" >> solveTask mY2 |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 2,024: | Line 2,024: | ||
===With Numeric.LinearAlgebra=== |
===With Numeric.LinearAlgebra=== |
||
< |
<syntaxhighlight lang="haskell">import Numeric.LinearAlgebra |
||
a1, a2 :: Matrix R |
a1, a2 :: Matrix R |
||
Line 2,040: | Line 2,040: | ||
main = do |
main = do |
||
print $ lu a1 |
print $ lu a1 |
||
print $ lu a2</ |
print $ lu a2</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>((3><3) |
<pre>((3><3) |
||
Line 2,072: | Line 2,072: | ||
'''Solution:''' |
'''Solution:''' |
||
<syntaxhighlight lang="idris"> |
|||
<lang Idris> |
|||
module Main |
module Main |
||
Line 2,228: | Line 2,228: | ||
putStrLn "Solution 2:" |
putStrLn "Solution 2:" |
||
printEx ex2 |
printEx ex2 |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
Line 2,252: | Line 2,252: | ||
'''Solution:''' |
'''Solution:''' |
||
< |
<syntaxhighlight lang="j">mp=: +/ .* |
||
LU=: 3 : 0 |
LU=: 3 : 0 |
||
Line 2,274: | Line 2,274: | ||
permtomat=: 1 {.~"0 -@>:@:/: |
permtomat=: 1 {.~"0 -@>:@:/: |
||
LUdecompose=: (permtomat&.>@{. , }.)@:LU</ |
LUdecompose=: (permtomat&.>@{. , }.)@:LU</syntaxhighlight> |
||
'''Example use:''' |
'''Example use:''' |
||
< |
<syntaxhighlight lang="j"> A=:3 3$1 3 5 2 4 7 1 1 0 |
||
LUdecompose A |
LUdecompose A |
||
┌─────┬─────┬───────┐ |
┌─────┬─────┬───────┐ |
||
Line 2,301: | Line 2,301: | ||
1 5 2 6 |
1 5 2 6 |
||
3 17 18 1 |
3 17 18 1 |
||
2 5 7 1</ |
2 5 7 1</syntaxhighlight> |
||
=={{header|Java}}== |
=={{header|Java}}== |
||
Translation of [[#Common_Lisp|Common Lisp]] via [[#D|D]] |
Translation of [[#Common_Lisp|Common Lisp]] via [[#D|D]] |
||
{{works with|Java|8}} |
{{works with|Java|8}} |
||
< |
<syntaxhighlight lang="java">import static java.util.Arrays.stream; |
||
import java.util.Locale; |
import java.util.Locale; |
||
import static java.util.stream.IntStream.range; |
import static java.util.stream.IntStream.range; |
||
Line 2,402: | Line 2,402: | ||
print(m); |
print(m); |
||
} |
} |
||
}</ |
}</syntaxhighlight> |
||
<pre> 1.0 0.0 0.0 |
<pre> 1.0 0.0 0.0 |
||
0.5 1.0 0.0 |
0.5 1.0 0.0 |
||
Line 2,433: | Line 2,433: | ||
=={{header|Javascript}}== |
=={{header|Javascript}}== |
||
{{works with|ES5 ES6}} |
{{works with|ES5 ES6}} |
||
< |
<syntaxhighlight lang="javascript"> |
||
const mult=(a, b)=>{ |
const mult=(a, b)=>{ |
||
let res = new Array(a.length); |
let res = new Array(a.length); |
||
Line 2,512: | Line 2,512: | ||
return [...lu(A),P]; |
return [...lu(A),P]; |
||
} |
} |
||
</syntaxhighlight> |
|||
</lang> |
|||
=={{header|jq}}== |
=={{header|jq}}== |
||
Line 2,521: | Line 2,521: | ||
'''Infrastructure''' |
'''Infrastructure''' |
||
< |
<syntaxhighlight lang="jq"># Create an m x n matrix |
||
def matrix(m; n; init): |
def matrix(m; n; init): |
||
if m == 0 then [] |
if m == 0 then [] |
||
Line 2,565: | Line 2,565: | ||
| reduce range (0;$length) as $i |
| reduce range (0;$length) as $i |
||
(""; . + reduce range(0;$length) as $j |
(""; . + reduce range(0;$length) as $j |
||
(""; "\(.) \($in[$i][$j] | right )" ) + "\n" ) ;</ |
(""; "\(.) \($in[$i][$j] | right )" ) + "\n" ) ;</syntaxhighlight> |
||
'''LU decomposition''' |
'''LU decomposition''' |
||
< |
<syntaxhighlight lang="jq"># Create the pivot matrix for the input matrix. |
||
# Use "range(0;$n) as $i" to handle ill-conditioned cases. |
# Use "range(0;$n) as $i" to handle ill-conditioned cases. |
||
def pivotize: |
def pivotize: |
||
Line 2,609: | Line 2,609: | ||
| . + [ $P ] |
| . + [ $P ] |
||
; |
; |
||
</syntaxhighlight> |
|||
</lang> |
|||
'''Example 1''': |
'''Example 1''': |
||
< |
<syntaxhighlight lang="jq">def a: [[1, 3, 5], [2, 4, 7], [1, 1, 0]]; |
||
a | lup[] | neatly(4) |
a | lup[] | neatly(4) |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{Out}} |
{{Out}} |
||
< |
<syntaxhighlight lang="sh"> $ /usr/local/bin/jq -M -n -r -f LU.jq |
||
1 0 0 |
1 0 0 |
||
0.5 1 0 |
0.5 1 0 |
||
Line 2,627: | Line 2,627: | ||
1 0 0 |
1 0 0 |
||
0 0 1 |
0 0 1 |
||
</syntaxhighlight> |
|||
</lang> |
|||
'''Example 2''': |
'''Example 2''': |
||
< |
<syntaxhighlight lang="jq">def b: [[11,9,24,2],[1,5,2,6],[3,17,18,1],[2,5,7,1]]; |
||
b | lup[] | neatly(21)</ |
b | lup[] | neatly(21)</syntaxhighlight> |
||
{{Out}} |
{{Out}} |
||
< |
<syntaxhighlight lang="sh">$ /usr/local/bin/jq -M -n -r -f LU.jq |
||
1 0 0 0 |
1 0 0 0 |
||
0.2727272727272727 1 0 0 |
0.2727272727272727 1 0 0 |
||
Line 2,646: | Line 2,646: | ||
0 0 1 0 |
0 0 1 0 |
||
0 1 0 0 |
0 1 0 0 |
||
0 0 0 1</ |
0 0 0 1</syntaxhighlight> |
||
'''Example 3''': |
'''Example 3''': |
||
<syntaxhighlight lang="jq"> |
|||
<lang jq> |
|||
# A|lup|verify(A) should be true |
# A|lup|verify(A) should be true |
||
def verify(A): |
def verify(A): |
||
Line 2,661: | Line 2,661: | ||
[0, 0, 1, -1]]; |
[0, 0, 1, -1]]; |
||
A|lup|verify(A)</ |
A|lup|verify(A)</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
true |
true |
||
Line 2,683: | Line 2,683: | ||
=={{header|Kotlin}}== |
=={{header|Kotlin}}== |
||
< |
<syntaxhighlight lang="scala">// version 1.1.4-3 |
||
typealias Vector = DoubleArray |
typealias Vector = DoubleArray |
||
Line 2,784: | Line 2,784: | ||
printMatrix("U:", u2, "%8.5f") |
printMatrix("U:", u2, "%8.5f") |
||
printMatrix("P:", p2, "%1.0f") |
printMatrix("P:", p2, "%1.0f") |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 2,846: | Line 2,846: | ||
=={{header|Lobster}}== |
=={{header|Lobster}}== |
||
< |
<syntaxhighlight lang="lobster">import std |
||
// derived from JAMA v1.03 |
// derived from JAMA v1.03 |
||
Line 2,950: | Line 2,950: | ||
print_L A |
print_L A |
||
print_U A |
print_U A |
||
print_P piv</ |
print_P piv</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 2,992: | Line 2,992: | ||
=={{header|Maple}}== |
=={{header|Maple}}== |
||
<syntaxhighlight lang="maple"> |
|||
<lang Maple> |
|||
A:=<<1.0|3.0|5.0>,<2.0|4.0|7.0>,<1.0|1.0|0.0>>: |
A:=<<1.0|3.0|5.0>,<2.0|4.0|7.0>,<1.0|1.0|0.0>>: |
||
LinearAlgebra:-LUDecomposition(A); |
LinearAlgebra:-LUDecomposition(A); |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 3,005: | Line 3,005: | ||
[0 0 1] [0.500000000000000 -1. 1.0] [0. 0. -2.] |
[0 0 1] [0.500000000000000 -1. 1.0] [0. 0. -2.] |
||
</pre> |
</pre> |
||
<syntaxhighlight lang="maple"> |
|||
<lang Maple> |
|||
A:=<<11.0|9.0|24.0|2.0>,<1.0|5.0|2.0|6.0>, |
A:=<<11.0|9.0|24.0|2.0>,<1.0|5.0|2.0|6.0>, |
||
<3.0|17.0|18.0|1.0>,<2.0|5.0|7.0|1.0>>: |
<3.0|17.0|18.0|1.0>,<2.0|5.0|7.0|1.0>>: |
||
Line 3,012: | Line 3,012: | ||
LUDecomposition(A); |
LUDecomposition(A); |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 3,042: | Line 3,042: | ||
=={{header|Mathematica}}/{{header|Wolfram Language}}== |
=={{header|Mathematica}}/{{header|Wolfram Language}}== |
||
< |
<syntaxhighlight lang="mathematica">(*Ex1*)a = {{1, 3, 5}, {2, 4, 7}, {1, 1, 0}}; |
||
{lu, p, c} = LUDecomposition[a]; |
{lu, p, c} = LUDecomposition[a]; |
||
l = LowerTriangularize[lu, -1] + IdentityMatrix[Length[p]]; |
l = LowerTriangularize[lu, -1] + IdentityMatrix[Length[p]]; |
||
Line 3,055: | Line 3,055: | ||
P = Part[IdentityMatrix[Length[p]], p] ; |
P = Part[IdentityMatrix[Length[p]], p] ; |
||
MatrixForm /@ {P.a , P, l, u, l.u} |
MatrixForm /@ {P.a , P, l, u, l.u} |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
[[File:LUex1.png]] |
[[File:LUex1.png]] |
||
Line 3,064: | Line 3,064: | ||
LU decomposition is part of language |
LU decomposition is part of language |
||
< |
<syntaxhighlight lang="matlab"> A = [ |
||
1 3 5 |
1 3 5 |
||
2 4 7 |
2 4 7 |
||
1 1 0]; |
1 1 0]; |
||
[L,U,P] = lu(A)</ |
[L,U,P] = lu(A)</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 3,091: | Line 3,091: | ||
</pre> |
</pre> |
||
2nd example: |
2nd example: |
||
< |
<syntaxhighlight lang="matlab"> A = [ |
||
11 9 24 2 |
11 9 24 2 |
||
1 5 2 6 |
1 5 2 6 |
||
Line 3,097: | Line 3,097: | ||
2 5 7 1 ]; |
2 5 7 1 ]; |
||
[L,U,P] = lu(A)</ |
[L,U,P] = lu(A)</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 3,123: | Line 3,123: | ||
===Creating a MATLAB function=== |
===Creating a MATLAB function=== |
||
<syntaxhighlight lang="matlab"> |
|||
<lang Matlab> |
|||
function [ P, L, U ] = LUdecomposition(A) |
function [ P, L, U ] = LUdecomposition(A) |
||
Line 3,181: | Line 3,181: | ||
end |
end |
||
</syntaxhighlight> |
|||
</lang> |
|||
=={{header|Maxima}}== |
=={{header|Maxima}}== |
||
< |
<syntaxhighlight lang="maxima">/* LU decomposition is built-in */ |
||
a: hilbert_matrix(4)$ |
a: hilbert_matrix(4)$ |
||
Line 3,220: | Line 3,220: | ||
lu_backsub(lup, transpose([1, 1, -1, -1])); |
lu_backsub(lup, transpose([1, 1, -1, -1])); |
||
/* matrix([-204], [2100], [-4740], [2940]) */</ |
/* matrix([-204], [2100], [-4740], [2940]) */</syntaxhighlight> |
||
=={{header|Nim}}== |
=={{header|Nim}}== |
||
Line 3,228: | Line 3,228: | ||
For display, we use the third party module "strfmt" which allows to specify dynamically the format. |
For display, we use the third party module "strfmt" which allows to specify dynamically the format. |
||
< |
<syntaxhighlight lang="nim">import macros, strutils |
||
import strfmt |
import strfmt |
||
Line 3,314: | Line 3,314: | ||
l2.print("L:", "8.5f") |
l2.print("L:", "8.5f") |
||
u2.print("U:", "8.5f") |
u2.print("U:", "8.5f") |
||
p2.print("P:", "1.0f")</ |
p2.print("P:", "1.0f")</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 3,367: | Line 3,367: | ||
=={{header|PARI/GP}}== |
=={{header|PARI/GP}}== |
||
< |
<syntaxhighlight lang="parigp">matlup(M) = |
||
{ |
{ |
||
my (L = matid(#M), U = M, P = L); |
my (L = matid(#M), U = M, P = L); |
||
Line 3,389: | Line 3,389: | ||
[L,U,P] \\ return L,U,P triple matrix |
[L,U,P] \\ return L,U,P triple matrix |
||
}</ |
}</syntaxhighlight> |
||
Output: |
Output: |
||
Line 3,454: | Line 3,454: | ||
=={{header|Perl}}== |
=={{header|Perl}}== |
||
{{trans|Raku}} |
{{trans|Raku}} |
||
< |
<syntaxhighlight lang="perl">use List::Util qw(sum); |
||
for $test ( |
for $test ( |
||
Line 3,538: | Line 3,538: | ||
print "$line\n"; |
print "$line\n"; |
||
} |
} |
||
</syntaxhighlight> |
|||
</lang> |
|||
{{out}} |
{{out}} |
||
<pre style="height:35ex">A matrix |
<pre style="height:35ex">A matrix |
||
Line 3,597: | Line 3,597: | ||
=={{header|Phix}}== |
=={{header|Phix}}== |
||
{{trans|Kotlin}} |
{{trans|Kotlin}} |
||
<!--< |
<!--<syntaxhighlight lang="phix">(phixonline)--> |
||
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span> |
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span> |
||
<span style="color: #008080;">function</span> <span style="color: #000000;">matrix_mul</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">,</span> <span style="color: #004080;">sequence</span> <span style="color: #000000;">b</span><span style="color: #0000FF;">)</span> |
<span style="color: #008080;">function</span> <span style="color: #000000;">matrix_mul</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">,</span> <span style="color: #004080;">sequence</span> <span style="color: #000000;">b</span><span style="color: #0000FF;">)</span> |
||
Line 3,674: | Line 3,674: | ||
<span style="color: #7060A8;">pp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">lu</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]),{</span><span style="color: #004600;">pp_Nest</span><span style="color: #0000FF;">,</span><span style="color: #000000;">2</span><span style="color: #0000FF;">,</span><span style="color: #004600;">pp_Pause</span><span style="color: #0000FF;">,</span><span style="color: #000000;">0</span><span style="color: #0000FF;">})</span> |
<span style="color: #7060A8;">pp</span><span style="color: #0000FF;">(</span><span style="color: #000000;">lu</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]),{</span><span style="color: #004600;">pp_Nest</span><span style="color: #0000FF;">,</span><span style="color: #000000;">2</span><span style="color: #0000FF;">,</span><span style="color: #004600;">pp_Pause</span><span style="color: #0000FF;">,</span><span style="color: #000000;">0</span><span style="color: #0000FF;">})</span> |
||
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span> |
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span> |
||
<!--</ |
<!--</syntaxhighlight>--> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 3,710: | Line 3,710: | ||
=={{header|PL/I}}== |
=={{header|PL/I}}== |
||
< |
<syntaxhighlight lang="pl/i">(subscriptrange, fofl, size): /* 2 Nov. 2013 */ |
||
LU_Decomposition: procedure options (main); |
LU_Decomposition: procedure options (main); |
||
declare a1(3,3) float (18) initial ( 1, 3, 5, |
declare a1(3,3) float (18) initial ( 1, 3, 5, |
||
Line 3,797: | Line 3,797: | ||
end LU_Decomposition; |
end LU_Decomposition; |
||
</syntaxhighlight> |
|||
</lang> |
|||
Derived from Fortran version above. |
Derived from Fortran version above. |
||
Results: |
Results: |
||
Line 3,841: | Line 3,841: | ||
=={{header|Python}}== |
=={{header|Python}}== |
||
{{trans|D}} |
{{trans|D}} |
||
< |
<syntaxhighlight lang="python">from pprint import pprint |
||
def matrixMul(A, B): |
def matrixMul(A, B): |
||
Line 3,882: | Line 3,882: | ||
for part in lu(b): |
for part in lu(b): |
||
pprint(part) |
pprint(part) |
||
print</ |
print</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>[[1.0, 0.0, 0.0], |
<pre>[[1.0, 0.0, 0.0], |
||
Line 3,913: | Line 3,913: | ||
=={{header|R}}== |
=={{header|R}}== |
||
< |
<syntaxhighlight lang="r">library(Matrix) |
||
A <- matrix(c(1, 3, 5, 2, 4, 7, 1, 1, 0), 3, 3, byrow=T) |
A <- matrix(c(1, 3, 5, 2, 4, 7, 1, 1, 0), 3, 3, byrow=T) |
||
dim(A) <- c(3, 3) |
dim(A) <- c(3, 3) |
||
expand(lu(A))</ |
expand(lu(A))</syntaxhighlight> |
||
{{Out}} |
{{Out}} |
||
Line 3,943: | Line 3,943: | ||
=={{header|Racket}}== |
=={{header|Racket}}== |
||
< |
<syntaxhighlight lang="racket"> |
||
#lang racket |
#lang racket |
||
(require math) |
(require math) |
||
Line 3,959: | Line 3,959: | ||
; #[0 -2 -3] |
; #[0 -2 -3] |
||
; #[0 0 -2]]) |
; #[0 0 -2]]) |
||
</syntaxhighlight> |
|||
</lang> |
|||
=={{header|Raku}}== |
=={{header|Raku}}== |
||
Line 3,965: | Line 3,965: | ||
Translation of Ruby. |
Translation of Ruby. |
||
<lang |
<syntaxhighlight lang="raku" line>for ( [1, 3, 5], # Test Matrices |
||
[2, 4, 7], |
[2, 4, 7], |
||
[1, 1, 0] |
[1, 1, 0] |
||
Line 4,034: | Line 4,034: | ||
say "\n$message"; |
say "\n$message"; |
||
$_».&rat-int.fmt("%7s").say for @array; |
$_».&rat-int.fmt("%7s").say for @array; |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre>A Matrix |
<pre>A Matrix |
||
Line 4,093: | Line 4,093: | ||
=={{header|REXX}}== |
=={{header|REXX}}== |
||
< |
<syntaxhighlight lang="rexx">/*REXX program creates a matrix from console input, performs/shows LU decomposition.*/ |
||
#= 0; P.= 0; PA.= 0; L.= 0; U.= 0 /*initialize some variables to zero. */ |
#= 0; P.= 0; PA.= 0; L.= 0; U.= 0 /*initialize some variables to zero. */ |
||
parse arg x /*obtain matrix elements from the C.L. */ |
parse arg x /*obtain matrix elements from the C.L. */ |
||
Line 4,154: | Line 4,154: | ||
end /*c*/ |
end /*c*/ |
||
say _ |
say _ |
||
end /*r*/; return</ |
end /*r*/; return</syntaxhighlight> |
||
{{out|output|text= when using the input of: <tt> 1 3 5 2 4 7 1 1 0 </tt>}} |
{{out|output|text= when using the input of: <tt> 1 3 5 2 4 7 1 1 0 </tt>}} |
||
<pre> |
<pre> |
||
Line 4,216: | Line 4,216: | ||
=={{header|Ruby}}== |
=={{header|Ruby}}== |
||
< |
<syntaxhighlight lang="ruby">require 'matrix' |
||
class Matrix |
class Matrix |
||
Line 4,280: | Line 4,280: | ||
l.pretty_print(" %8.5f", "L") |
l.pretty_print(" %8.5f", "L") |
||
u.pretty_print(" %8.5f", "U") |
u.pretty_print(" %8.5f", "U") |
||
p.pretty_print(" %d", "P")</ |
p.pretty_print(" %d", "P")</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 4,326: | Line 4,326: | ||
Matrix has a <code>lup_decomposition</code> built-in method. |
Matrix has a <code>lup_decomposition</code> built-in method. |
||
< |
<syntaxhighlight lang="ruby">l, u, p = a.lup_decomposition |
||
l.pretty_print(" %8.5f", "L") |
l.pretty_print(" %8.5f", "L") |
||
u.pretty_print(" %8.5f", "U") |
u.pretty_print(" %8.5f", "U") |
||
p.pretty_print(" %d", "P")</ |
p.pretty_print(" %d", "P")</syntaxhighlight> |
||
Output is the same. |
Output is the same. |
||
=={{header|Rust}}== |
=={{header|Rust}}== |
||
{{libheader| ndarray}} |
{{libheader| ndarray}} |
||
<syntaxhighlight lang="rust"> |
|||
<lang Rust> |
|||
#![allow(non_snake_case)] |
#![allow(non_snake_case)] |
||
use ndarray::{Array, Axis, Array2, arr2, Zip, NdFloat, s}; |
use ndarray::{Array, Axis, Array2, arr2, Zip, NdFloat, s}; |
||
Line 4,421: | Line 4,421: | ||
(L, U, P) |
(L, U, P) |
||
} |
} |
||
</syntaxhighlight> |
|||
</lang> |
|||
Line 4,469: | Line 4,469: | ||
=={{header|Sidef}}== |
=={{header|Sidef}}== |
||
{{trans|Raku}} |
{{trans|Raku}} |
||
< |
<syntaxhighlight lang="ruby">func is_square(m) { m.all { .len == m.len } } |
||
func matrix_zero(n, m=n) { m.of { n.of(0) } } |
func matrix_zero(n, m=n) { m.of { n.of(0) } } |
||
func matrix_ident(n) { n.of {|i| n.of {|j| i==j ? 1 : 0 } } } |
func matrix_ident(n) { n.of {|i| n.of {|j| i==j ? 1 : 0 } } } |
||
Line 4,542: | Line 4,542: | ||
say_it(a, b) |
say_it(a, b) |
||
} |
} |
||
}</ |
}</syntaxhighlight> |
||
<pre style="height: 40ex; overflow: scroll"> |
<pre style="height: 40ex; overflow: scroll"> |
||
A Matrix |
A Matrix |
||
Line 4,604: | Line 4,604: | ||
See [http://www.stata.com/help.cgi?mf_lud LU decomposition] in Stata help. |
See [http://www.stata.com/help.cgi?mf_lud LU decomposition] in Stata help. |
||
< |
<syntaxhighlight lang="stata">mata |
||
: lud(a=(1,3,5\2,4,7\1,1,0),l=.,u=.,p=.) |
: lud(a=(1,3,5\2,4,7\1,1,0),l=.,u=.,p=.) |
||
Line 4,637: | Line 4,637: | ||
2 | 1 | |
2 | 1 | |
||
3 | 3 | |
3 | 3 | |
||
+-----+</ |
+-----+</syntaxhighlight> |
||
=== Implementation === |
=== Implementation === |
||
< |
<syntaxhighlight lang="stata">void ludec(real matrix a, real matrix l, real matrix u, real vector p) { |
||
real scalar i,j,n,s |
real scalar i,j,n,s |
||
real vector js |
real vector js |
||
Line 4,662: | Line 4,662: | ||
u = uppertriangle(l) |
u = uppertriangle(l) |
||
l = lowertriangle(l, 1) |
l = lowertriangle(l, 1) |
||
}</ |
}</syntaxhighlight> |
||
'''Example''': |
'''Example''': |
||
< |
<syntaxhighlight lang="stata">: ludec(a=(1,3,5\2,4,7\1,1,0),l=.,u=.,p=.) |
||
: a |
: a |
||
Line 4,697: | Line 4,697: | ||
2 | 1 | |
2 | 1 | |
||
3 | 3 | |
3 | 3 | |
||
+-----+</ |
+-----+</syntaxhighlight> |
||
=={{header|Tcl}}== |
=={{header|Tcl}}== |
||
< |
<syntaxhighlight lang="tcl">package require Tcl 8.5 |
||
namespace eval matrix { |
namespace eval matrix { |
||
namespace path {::tcl::mathfunc ::tcl::mathop} |
namespace path {::tcl::mathfunc ::tcl::mathop} |
||
Line 4,767: | Line 4,767: | ||
return $s |
return $s |
||
} |
} |
||
}</ |
}</syntaxhighlight> |
||
Support code: |
Support code: |
||
< |
<syntaxhighlight lang="tcl"># Code adapted from Matrix_multiplication and Matrix_transposition tasks |
||
namespace eval matrix { |
namespace eval matrix { |
||
# Get the size of a matrix; assumes that all rows are the same length, which |
# Get the size of a matrix; assumes that all rows are the same length, which |
||
Line 4,817: | Line 4,817: | ||
return $max |
return $max |
||
} |
} |
||
}</ |
}</syntaxhighlight> |
||
Demonstrating: |
Demonstrating: |
||
< |
<syntaxhighlight lang="tcl"># This does the decomposition and prints it out nicely |
||
proc demo {A} { |
proc demo {A} { |
||
lassign [matrix::luDecompose $A] L U P |
lassign [matrix::luDecompose $A] L U P |
||
Line 4,831: | Line 4,831: | ||
demo {{1 3 5} {2 4 7} {1 1 0}} |
demo {{1 3 5} {2 4 7} {1 1 0}} |
||
puts "=================================" |
puts "=================================" |
||
demo {{11 9 24 2} {1 5 2 6} {3 17 18 1} {2 5 7 1}}</ |
demo {{11 9 24 2} {1 5 2 6} {3 17 18 1} {2 5 7 1}}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 4,881: | Line 4,881: | ||
=={{header|VBA}}== |
=={{header|VBA}}== |
||
{{trans|Phix}} |
{{trans|Phix}} |
||
< |
<syntaxhighlight lang="vb">Option Base 1 |
||
Private Function pivotize(m As Variant) As Variant |
Private Function pivotize(m As Variant) As Variant |
||
Dim n As Integer: n = UBound(m) |
Dim n As Integer: n = UBound(m) |
||
Line 4,975: | Line 4,975: | ||
Debug.Print |
Debug.Print |
||
Next i |
Next i |
||
End Sub</ |
End Sub</syntaxhighlight>{{out}} |
||
<pre>== a,l,u,p: == |
<pre>== a,l,u,p: == |
||
1 3 5 |
1 3 5 |
||
Line 5,017: | Line 5,017: | ||
{{libheader|Wren-matrix}} |
{{libheader|Wren-matrix}} |
||
{{libheader|Wren-fmt}} |
{{libheader|Wren-fmt}} |
||
< |
<syntaxhighlight lang="ecmascript">import "/matrix" for Matrix |
||
import "/fmt" for Fmt |
import "/fmt" for Fmt |
||
Line 5,043: | Line 5,043: | ||
Fmt.mprint(lup[2], 2, 0) |
Fmt.mprint(lup[2], 2, 0) |
||
System.print() |
System.print() |
||
}</ |
}</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
Line 5,102: | Line 5,102: | ||
=={{header|zkl}}== |
=={{header|zkl}}== |
||
Using the GNU Scientific Library, which does the decomposition without returning the permutations: |
Using the GNU Scientific Library, which does the decomposition without returning the permutations: |
||
< |
<syntaxhighlight lang="zkl">var [const] GSL=Import("zklGSL"); // libGSL (GNU Scientific Library) |
||
fcn luTask(A){ |
fcn luTask(A){ |
||
A.LUDecompose(); // in place, contains L & U |
A.LUDecompose(); // in place, contains L & U |
||
Line 5,119: | Line 5,119: | ||
2.0, 5.0, 7.0, 1.0); |
2.0, 5.0, 7.0, 1.0); |
||
L,U:=luTask(A); |
L,U:=luTask(A); |
||
println("L:\n",L.format(8,4),"\nU:\n",U.format(8,4));</ |
println("L:\n",L.format(8,4),"\nU:\n",U.format(8,4));</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 5,146: | Line 5,146: | ||
A matrix is a list of lists, ie list of rows in row major order. |
A matrix is a list of lists, ie list of rows in row major order. |
||
< |
<syntaxhighlight lang="zkl">fcn make_array(n,m,v){ (m).pump(List.createLong(m).write,v)*n } |
||
fcn eye(n){ // Creates a nxn identity matrix. |
fcn eye(n){ // Creates a nxn identity matrix. |
||
I:=make_array(n,n,0.0); |
I:=make_array(n,n,0.0); |
||
Line 5,195: | Line 5,195: | ||
foreach i,j,k in (m,p,n){ ans[i][j]+=a[i][k]*b[k][j]; } |
foreach i,j,k in (m,p,n){ ans[i][j]+=a[i][k]*b[k][j]; } |
||
ans |
ans |
||
}</ |
}</syntaxhighlight> |
||
Example 1 |
Example 1 |
||
< |
<syntaxhighlight lang="zkl">g:=L(L(1.0,3.0,5.0),L(2.0,4.0,7.0),L(1.0,1.0,0.0)); |
||
lu(g).apply2("println");</ |
lu(g).apply2("println");</syntaxhighlight> |
||
{{out}} |
{{out}} |
||
<pre> |
<pre> |
||
Line 5,206: | Line 5,206: | ||
</pre> |
</pre> |
||
Example 2 |
Example 2 |
||
< |
<syntaxhighlight lang="zkl">lu(L( L(11.0, 9.0, 24.0, 2.0), |
||
L( 1.0, 5.0, 2.0, 6.0), |
L( 1.0, 5.0, 2.0, 6.0), |
||
L( 3.0, 17.0, 18.0, 1.0), |
L( 3.0, 17.0, 18.0, 1.0), |
||
Line 5,212: | Line 5,212: | ||
fcn printM(m) { m.pump(Console.println,rowFmt) } |
fcn printM(m) { m.pump(Console.println,rowFmt) } |
||
fcn rowFmt(row){ ("%9.5f "*row.len()).fmt(row.xplode()) }</ |
fcn rowFmt(row){ ("%9.5f "*row.len()).fmt(row.xplode()) }</syntaxhighlight> |
||
The list apply2 method is side effects only, it doesn't aggregate results. When given a list of actions, it applies the action and passes the result to the next action. The fpM method is partial application with a mask, "-" truncates the parameters at that point (in this case, no parameters, ie just print a blank line, not the result of printM). |
The list apply2 method is side effects only, it doesn't aggregate results. When given a list of actions, it applies the action and passes the result to the next action. The fpM method is partial application with a mask, "-" truncates the parameters at that point (in this case, no parameters, ie just print a blank line, not the result of printM). |
||
{{out}} |
{{out}} |