Rosetta Code

Cholesky decomposition: Difference between revisions

Page Discussion

Cholesky decomposition (view source)

Revision as of 19:41, 26 August 2022

2,515 bytes added ,  1 year ago
m
syntax highlighting fixup automation
No edit summary
m (syntax highlighting fixup automation)
Line 82:
{{trans|Python}}
 
<langsyntaxhighlight lang=11l>F cholesky(A)
V l = [[0.0] * A.len] * A.len
L(i) 0 .< A.len
Line 106:
[54, 86, 174, 134],
[42, 62, 134, 106]]
pprint(cholesky(m2))</langsyntaxhighlight>
 
{{out}}
Line 123:
{{works with|Ada 2005}}
decomposition.ads:
<langsyntaxhighlight lang=Ada>with Ada.Numerics.Generic_Real_Arrays;
generic
with package Matrix is new Ada.Numerics.Generic_Real_Arrays (<>);
Line 131:
procedure Decompose (A : Matrix.Real_Matrix; L : out Matrix.Real_Matrix);
 
end Decomposition;</langsyntaxhighlight>
 
decomposition.adb:
<langsyntaxhighlight lang=Ada>with Ada.Numerics.Generic_Elementary_Functions;
 
package body Decomposition is
Line 166:
end loop;
end Decompose;
end Decomposition;</langsyntaxhighlight>
 
Example usage:
<langsyntaxhighlight lang=Ada>with Ada.Numerics.Real_Arrays;
with Ada.Text_IO;
with Decomposition;
Line 213:
Ada.Text_IO.Put_Line ("A:"); Print (Example_2);
Ada.Text_IO.Put_Line ("L:"); Print (L_2);
end Decompose_Example;</langsyntaxhighlight>
{{out}}
<pre>Example 1:
Line 242:
{{works with|ALGOL 68G|Any - tested with release [http://sourceforge.net/projects/algol68/files/algol68g/algol68g-1.18.0/algol68g-1.18.0-9h.tiny.el5.centos.fc11.i386.rpm/download 1.18.0-9h.tiny].}}
{{wont work with|ELLA ALGOL 68|Any (with appropriate job cards) - tested with release [http://sourceforge.net/projects/algol68/files/algol68toc/algol68toc-1.8.8d/algol68toc-1.8-8d.fc9.i386.rpm/download 1.8-8d] - due to extensive use of '''format'''[ted] ''transput''.}}
<langsyntaxhighlight lang=algol68>#!/usr/local/bin/a68g --script #
 
MODE FIELD=LONG REAL;
Line 300:
MAT c2 = cholesky(m2);
print matrix(c2)
)</langsyntaxhighlight>
{{out}}
<pre>
Line 313:
 
=={{header|AutoHotkey}}==
<langsyntaxhighlight lang=AutoHotkey>Cholesky_Decomposition(A){
L := [], n := A.Count()
L[1,1] := Sqrt(A[1,1])
Line 347:
}
return "[" Trim(output, "`n,") "]"
}</langsyntaxhighlight>
Examples:<langsyntaxhighlight lang=AutoHotkey>A := [[25, 15, -5]
, [15, 18, 0]
, [-5, 0 , 11]]
Line 360:
 
MsgBox % Result := ShowMatrix(L1) "`n----`n" ShowMatrix(L2) "`n----"
return</langsyntaxhighlight>
{{out}}
<pre>[[5.000, 0.000, 0.000]
Line 374:
=={{header|BBC BASIC}}==
{{works with|BBC BASIC for Windows}}
<langsyntaxhighlight lang=bbcbasic> DIM m1(2,2)
m1() = 25, 15, -5, \
\ 15, 18, 0, \
Line 419:
PRINT
NEXT row%
ENDPROC</langsyntaxhighlight>
'''Output:'''
<pre>
Line 433:
 
=={{header|C}}==
<langsyntaxhighlight lang=c>#include <stdio.h>
#include <stdlib.h>
#include <math.h>
Line 483:
 
return 0;
}</langsyntaxhighlight>
{{out}}
<pre>5.00000 0.00000 0.00000
Line 495:
 
=={{header|C sharp|C#}}==
<langsyntaxhighlight lang=csharp>using System;
using System.Collections.Generic;
using System.Linq;
Line 597:
}
}
}</langsyntaxhighlight>
 
{{out}}
Line 626:
 
=={{header|C++}}==
<langsyntaxhighlight lang=cpp>#include <cassert>
#include <cmath>
#include <iomanip>
Line 727:
 
return 0;
}</langsyntaxhighlight>
 
{{out}}
Line 753:
=={{header|Clojure}}==
{{trans|Python}}
<langsyntaxhighlight lang=clojure>(defn cholesky
[matrix]
(let [n (count matrix)
Line 763:
(Math/sqrt (- (aget A i i) s))
(* (/ 1.0 (aget L j j)) (- (aget A i j) s))))))
(vec (map vec L))))</langsyntaxhighlight>
Example:
<langsyntaxhighlight lang=clojure>(cholesky [[25 15 -5] [15 18 0] [-5 0 11]])
;=> [[ 5.0 0.0 0.0]
; [ 3.0 3.0 0.0]
Line 774:
; [ 5.185449728701349 6.565905201197403 0.0 0.0 ]
; [12.727922061357857 3.0460384954008553 1.6497422479090704 0.0 ]
; [ 9.899494936611667 1.624553864213788 1.8497110052313648 1.3926212476456026]]</langsyntaxhighlight>
 
=={{header|Common Lisp}}==
<langsyntaxhighlight lang=lisp>;; Calculates the Cholesky decomposition matrix L
;; for a positive-definite, symmetric nxn matrix A.
(defun chol (A)
Line 805:
 
;; Return the calculated matrix L.
L))</langsyntaxhighlight>
 
<langsyntaxhighlight lang=lisp>;; Example 1:
(setf A (make-array '(3 3) :initial-contents '((25 15 -5) (15 18 0) (-5 0 11))))
(chol A)
#2A((5.0 0 0)
(3.0 3.0 0)
(-1.0 1.0 3.0))</langsyntaxhighlight>
 
<langsyntaxhighlight lang=lisp>;; Example 2:
(setf B (make-array '(4 4) :initial-contents '((18 22 54 42) (22 70 86 62) (54 86 174 134) (42 62 134 106))))
(chol B)
Line 820:
(5.18545 6.565905 0 0)
(12.727922 3.0460374 1.6497375 0)
(9.899495 1.6245536 1.849715 1.3926151))</langsyntaxhighlight>
 
<langsyntaxhighlight lang=lisp>;; case of matrix stored as a list of lists (inner lists are rows of matrix)
;; as above, returns the Cholesky decomposition matrix of a square positive-definite, symmetric matrix
(defun cholesky (m)
Line 832:
(setf (cdr (last l)) (list (nconc x (list (sqrt (- (elt cm (incf j)) (*v x x))))))))))
;; where *v is the scalar product defined as
(defun *v (v1 v2) (reduce #'+ (mapcar #'* v1 v2)))</langsyntaxhighlight>
 
<langsyntaxhighlight lang=lisp>;; example 1
CL-USER> (setf a '((25 15 -5) (15 18 0) (-5 0 11)))
((25 15 -5) (15 18 0) (-5 0 11))
Line 843:
3 3 0
-1 1 3
NIL</langsyntaxhighlight>
 
<langsyntaxhighlight lang=lisp>;; example 2
CL-USER> (setf a '((18 22 54 42) (22 70 86 62) (54 86 174 134) (42 62 134 106)))
((18 22 54 42) (22 70 86 62) (54 86 174 134) (42 62 134 106))
Line 855:
12.72792 3.04604 1.64974 0.00000
9.89950 1.62455 1.84971 1.39262
NIL</langsyntaxhighlight>
 
=={{header|D}}==
<langsyntaxhighlight lang=d>import std.stdio, std.math, std.numeric;
 
T[][] cholesky(T)(in T[][] A) pure nothrow /*@safe*/ {
Line 884:
[42, 62, 134, 106]];
writefln("%(%(%2.3f %)\n%)", m2.cholesky);
}</langsyntaxhighlight>
{{out}}
<pre> 5 0 0
Line 898:
=={{header|DWScript}}==
{{Trans|C}}
<langsyntaxhighlight lang=delphi>function Cholesky(a : array of Float) : array of Float;
var
i, j, k, n : Integer;
Line 943:
42.0, 62.0, 134.0, 106.0 ];
var c2 := Cholesky(m2);
ShowMatrix(c2);</langsyntaxhighlight>
 
=={{header|F_Sharp|F#}}==
 
<langsyntaxhighlight lang=fsharp>open Microsoft.FSharp.Collections
 
let cholesky a =
Line 992:
printfn "ex2:"
printMat ex2
</syntaxhighlight>
</lang>
{{out}}
<pre>ex1:
Line 1,006:
 
=={{header|Fantom}}==
<langsyntaxhighlight lang=fantom>**
** Cholesky decomposition
**
Line 1,060:
runTest ([[18f,22f,54f,42f],[22f,70f,86f,62f],[54f,86f,174f,134f],[42f,62f,134f,106f]])
}
}</langsyntaxhighlight>
{{out}}
<pre>
Line 1,070:
 
=={{header|Fortran}}==
<langsyntaxhighlight lang=Fortran>Program Cholesky_decomp
! *************************************************!
! LBH @ ULPGC 06/03/2014
Line 1,139:
end do
 
End program Cholesky_decomp</langsyntaxhighlight>
{{out}}
<pre>
Line 1,149:
=={{header|FreeBASIC}}==
{{trans|BBC BASIC}}
<langsyntaxhighlight lang=freebasic>' version 18-01-2017
' compile with: fbc -s console
 
Line 1,214:
Print : Print "hit any key to end program"
Sleep
End</langsyntaxhighlight>
{{out}}
<pre> 5.00000 0.00000 0.00000
Line 1,228:
===Real===
This version works with real matrices, like most other solutions on the page. The representation is packed, however, storing only the lower triange of the input symetric matrix and the output lower matrix. The decomposition algorithm computes rows in order from top to bottom but is a little different thatn Cholesky–Banachiewicz.
<langsyntaxhighlight lang=go>package main
 
import (
Line 1,332:
fmt.Println("L:")
a.choleskyLower().print()
}</langsyntaxhighlight>
{{out}}
<pre>
Line 1,356:
===Hermitian===
This version handles complex Hermitian matricies as described on the WP page. The matrix representation is flat, and storage is allocated for all elements, not just the lower triangles. The decomposition algorithm is Cholesky–Banachiewicz.
<langsyntaxhighlight lang=go>package main
 
import (
Line 1,430:
a.print(heading)
a.choleskyDecomp().print("Cholesky factor L:")
}</langsyntaxhighlight>
{{out}}
<pre>
Line 1,465:
 
===Library gonum/mat===
<langsyntaxhighlight lang=go>package main
 
import (
Line 1,491:
42, 62, 134, 106,
}))
}</langsyntaxhighlight>
{{out}}
<pre>
Line 1,505:
 
===Library go.matrix===
<langsyntaxhighlight lang=go>package main
 
import (
Line 1,537:
fmt.Println("L:")
fmt.Println(l)
}</langsyntaxhighlight>
Output:
<pre>
Line 1,562:
=={{header|Groovy}}==
{{Trans|Java}}
<langsyntaxhighlight lang=groovy>def decompose = { a ->
assert a.size > 0 && a[0].size == a.size
def m = a.size
Line 1,575:
}
l
}</langsyntaxhighlight>
Test:
<langsyntaxhighlight lang=groovy>def test1 = [[25, 15, -5],
[15, 18, 0],
[-5, 0, 11]]
Line 1,589:
println()
decompose(test).each { println it[0..<(test.size)] }
}</langsyntaxhighlight>
{{out}}
<pre>[5.0, 0, 0]
Line 1,606:
 
The Cholesky module:
<langsyntaxhighlight lang=haskell>module Cholesky (Arr, cholesky) where
 
import Data.Array.IArray
Line 1,634:
return l
where maxBnd = fst . snd . bounds
update a l i = unsafeFreeze l >>= \l' -> writeArray l i (get a l' i)</langsyntaxhighlight>
The main module:
<langsyntaxhighlight lang=haskell>import Data.Array.IArray
import Data.List
import Cholesky
Line 1,665:
main = do
putStrLn $ showMatrice 3 $ elems $ cholesky ex1
putStrLn $ showMatrice 4 $ elems $ cholesky ex2</langsyntaxhighlight>
<b>output:</b>
<pre>
Line 1,678:
 
===With Numeric.LinearAlgebra===
<langsyntaxhighlight lang=haskell>import Numeric.LinearAlgebra
 
a,b :: Matrix R
Line 1,701:
print $ tr $ chol sb
 
</syntaxhighlight>
</lang>
{{out}}
<pre>Herm (3><3)
Line 1,728:
 
=={{header|Icon}} and {{header|Unicon}}==
<langsyntaxhighlight lang=Icon>procedure cholesky (array)
result := make_square_array (*array)
every (i := 1 to *array) do {
Line 1,768:
do_cholesky ([[25,15,-5],[15,18,0],[-5,0,11]])
do_cholesky ([[18,22,54,42],[22,70,86,62],[54,86,174,134],[42,62,134,106]])
end</langsyntaxhighlight>
{{out}}
<pre>
Line 1,795:
 
'''Solution:'''
<langsyntaxhighlight lang=Idris>module Main
 
import Data.Vect
Line 1,862:
print ex2
putStrLn "\n"
</syntaxhighlight>
</lang>
 
{{out}}
Line 1,873:
=={{header|J}}==
'''Solution:'''
<langsyntaxhighlight lang=j>mp=: +/ . * NB. matrix product
h =: +@|: NB. conjugate transpose
 
Line 1,887:
L0,(T mp L0),.L1
end.
)</langsyntaxhighlight>
See [[j:Essays/Cholesky Decomposition|Cholesky Decomposition essay]] on the J Wiki.
{{out|Examples}}
<langsyntaxhighlight lang=j> eg1=: 25 15 _5 , 15 18 0 ,: _5 0 11
eg2=: 18 22 54 42 , 22 70 86 62 , 54 86 174 134 ,: 42 62 134 106
cholesky eg1
Line 1,900:
5.18545 6.56591 0 0
12.7279 3.04604 1.64974 0
9.89949 1.62455 1.84971 1.39262</langsyntaxhighlight>
'''Using `math/lapack` addon'''
<langsyntaxhighlight lang=j> load 'math/lapack'
load 'math/lapack/potrf'
potrf_jlapack_ eg1
Line 1,912:
5.18545 6.56591 0 0
12.7279 3.04604 1.64974 0
9.89949 1.62455 1.84971 1.39262</langsyntaxhighlight>
 
=={{header|Java}}==
{{works with|Java|1.5+}}
<langsyntaxhighlight lang=java5>import java.util.Arrays;
 
public class Cholesky {
Line 1,946:
System.out.println(Arrays.deepToString(chol(test2)));
}
}</langsyntaxhighlight>
{{out}}
<pre>[[5.0, 0.0, 0.0], [3.0, 3.0, 0.0], [-1.0, 1.0, 3.0]]
Line 1,952:
 
=={{header|JavaScript}}==
<langsyntaxhighlight lang=javascript>
const cholesky = function (array) {
const zeros = [...Array(array.length)].map( _ => Array(array.length).fill(0));
Line 1,966:
let arr4 = [[18, 22, 54, 42], [22, 70, 86, 62], [54, 86, 174, 134], [42, 62, 134, 106]];
console.log(cholesky(arr4));
</syntaxhighlight>
</lang>
 
{{output}}
Line 1,983:
{{Works with|jq|1.4}}
'''Infrastructure''':
<langsyntaxhighlight lang=jq># Create an m x n matrix
def matrix(m; n; init):
if m == 0 then []
Line 2,021:
else [ ., 0, 0, length ] | test
end ;
</langsyntaxhighlight>'''Cholesky Decomposition''':<syntaxhighlight lang =jq>def cholesky_factor:
if is_symmetric then
length as $length
Line 2,042:
end ))
else error( "cholesky_factor: matrix is not symmetric" )
end ;</langsyntaxhighlight>
'''Task 1''':
[[25,15,-5],[15,18,0],[-5,0,11]] | cholesky_factor
Line 2,053:
[42, 62, 134, 106]] | cholesky_factor | neatly(20)
{{Out}}
<langsyntaxhighlight lang=jq> 4.242640687119285 0 0 0
5.185449728701349 6.565905201197403 0 0
12.727922061357857 3.0460384954008553 1.6497422479090704 0
9.899494936611665 1.6245538642137891 1.849711005231382 1.3926212476455924</langsyntaxhighlight>
 
=={{header|Julia}}==
Julia's strong linear algebra support includes Cholesky decomposition.
<langsyntaxhighlight lang=Julia>
a = [25 15 5; 15 18 0; -5 0 11]
b = [18 22 54 22; 22 70 86 62; 54 86 174 134; 42 62 134 106]
Line 2,066:
println(a, "\n => \n", chol(a, :L))
println(b, "\n => \n", chol(b, :L))
</syntaxhighlight>
</lang>
 
{{out}}
Line 2,090:
=={{header|Kotlin}}==
{{trans|C}}
<langsyntaxhighlight lang=scala>// version 1.0.6
 
fun cholesky(a: DoubleArray): DoubleArray {
Line 2,129:
val c2 = cholesky(m2)
showMatrix(c2)
}</langsyntaxhighlight>
 
{{out}}
Line 2,146:
{{trans|Go}}
Translated from the Go Real version: This version works with real matrices, like most other solutions on the page. The representation is packed, however, storing only the lower triange of the input symetric matrix and the output lower matrix. The decomposition algorithm computes rows in order from top to bottom but is a little different than Cholesky–Banachiewicz.
<langsyntaxhighlight lang=Lobster>import std
 
// choleskyLower returns the cholesky decomposition of a symmetric real
Line 2,225:
54.0, 86.0, 174.0,
42.0, 62.0, 134.0, 106.0])
</syntaxhighlight>
</lang>
{{out}}
<pre>
Line 2,250:
=={{header|Maple}}==
The Cholesky decomposition is obtained by passing the `method = Cholesky' option to the LUDecomposition procedure in the LinearAlgebra pacakge. This is illustrated below for the two requested examples. The first is computed exactly; the second is also, but the subsequent application of `evalf' to the result produces a matrix with floating point entries which can be compared with the expected output in the problem statement.
<langsyntaxhighlight lang=Maple>> A := << 25, 15, -5; 15, 18, 0; -5, 0, 11 >>;
[25 15 -5]
[ ]
Line 2,301:
[12.72792206 3.046038495 1.649742248 0. ]
[ ]
[9.899494934 1.624553864 1.849711006 1.392621248]</langsyntaxhighlight>
 
=={{header|Mathematica}} / {{header|Wolfram Language}}==
<langsyntaxhighlight lang=Mathematica>CholeskyDecomposition[{{25, 15, -5}, {15, 18, 0}, {-5, 0, 11}}]</langsyntaxhighlight>
Without the use of built-in functions, making use of memoization:
<langsyntaxhighlight lang=Mathematica>chol[A_] :=
Module[{L},
L[k_, k_] := L[k, k] = Sqrt[A[[k, k]] - Sum[L[k, j]^2, {j, 1, k-1}]];
L[i_, k_] := L[i, k] = L[k, k]^-1 (A[[i, k]] - Sum[L[i, j] L[k, j], {j, 1, k-1}]);
PadRight[Table[L[i, j], {i, Length[A]}, {j, i}]]
]</langsyntaxhighlight>
 
=={{header|MATLAB}} / {{header|Octave}}==
The cholesky decomposition chol() is an internal function
<langsyntaxhighlight lang=Matlab> A = [
25 15 -5
15 18 0
Line 2,328:
[L] = chol(A,'lower')
[L] = chol(B,'lower')
</syntaxhighlight>
</lang>
{{out}}
<pre> > [L] = chol(A,'lower')
Line 2,346:
 
=={{header|Maxima}}==
<langsyntaxhighlight lang=maxima>/* Cholesky decomposition is built-in */
 
a: hilbert_matrix(4)$
Line 2,357:
b . transpose(b) - a;
matrix([0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 0, 0])</langsyntaxhighlight>
 
=={{header|Nim}}==
{{trans|C}}
<langsyntaxhighlight lang=nim>import math, strutils, strformat
 
type Matrix[N: static int, T: SomeFloat] = array[N, array[N, T]]
Line 2,392:
[54.0, 86.0, 174.0, 134.0],
[42.0, 62.0, 134.0, 106.0]]
echo cholesky(m2)</langsyntaxhighlight>
 
{{out}}
Line 2,407:
{{trans|C}}
 
<langsyntaxhighlight lang=objeck>
class Cholesky {
function : Main(args : String[]) ~ Nil {
Line 2,451:
}
}
</syntaxhighlight>
</lang>
 
<pre>
Line 2,465:
 
=={{header|OCaml}}==
<langsyntaxhighlight lang=OCaml>let cholesky inp =
let n = Array.length inp in
let res = Array.make_matrix n n 0.0 in
Line 2,494:
[|22.0; 70.0; 86.0; 62.0|];
[|54.0; 86.0; 174.0; 134.0|];
[|42.0; 62.0; 134.0; 106.0|] |];</langsyntaxhighlight>
{{out}}
<pre>in:
Line 2,518:
=={{header|ooRexx}}==
{{trans|REXX}}
<langsyntaxhighlight lang=oorexx>/*REXX program performs the Cholesky decomposition on a square matrix. */
niner = '25 15 -5' , /*define a 3x3 matrix. */
'15 18 0' ,
Line 2,568:
do j=0 while h>9; m.j=h; h=h%2+1; end /*j*/
do k=j+5 to 0 by -1; numeric digits m.k; g=(g+x/g)*.5; end /*k*/
numeric digits d; return (g/1)i /*make complex if X < 0.*/</langsyntaxhighlight>
 
=={{header|PARI/GP}}==
 
<langsyntaxhighlight lang=parigp>cholesky(M) =
{
my (L = matrix(#M,#M));
Line 2,583:
);
L
}</langsyntaxhighlight>
 
Output: (set displayed digits with: \p 5)
Line 2,607:
 
=={{header|Pascal}}==
<langsyntaxhighlight lang=pascal>program CholeskyApp;
 
type
Line 2,674:
cOut := cholesky(cIn);
printM(cOut);
end.</langsyntaxhighlight>
{{out}}
<pre>
Line 2,688:
 
=={{header|Perl}}==
<langsyntaxhighlight lang=perl>sub cholesky {
my $matrix = shift;
my $chol = [ map { [(0) x @$matrix ] } @$matrix ];
Line 2,713:
print "\nExample 2:\n";
print +(map { sprintf "%7.4f\t", $_ } @$_), "\n" for @{ cholesky $example2 };
</syntaxhighlight>
</lang>
{{out}}
<pre>
Line 2,730:
=={{header|Phix}}==
{{trans|Sidef}}
<!--<langsyntaxhighlight lang=Phix>(phixonline)-->
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
<span style="color: #008080;">function</span> <span style="color: #000000;">cholesky</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">matrix</span><span style="color: #0000FF;">)</span>
Line 2,755:
<span style="color: #0000FF;">{</span> <span style="color: #000000;">54</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">86</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">174</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">134</span><span style="color: #0000FF;">},</span>
<span style="color: #0000FF;">{</span> <span style="color: #000000;">42</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">62</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">134</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">106</span><span style="color: #0000FF;">}}))</span>
<!--</langsyntaxhighlight>-->
{{out}}
<pre>
Line 2,768:
 
=={{header|PicoLisp}}==
<langsyntaxhighlight lang=PicoLisp>(scl 9)
(load "@lib/math.l")
 
Line 2,784:
(for R L
(for N R (prin (align 9 (round N 5))))
(prinl) ) ) )</langsyntaxhighlight>
Test:
<langsyntaxhighlight lang=PicoLisp>(cholesky
'((25.0 15.0 -5.0) (15.0 18.0 0) (-5.0 0 11.0)) )
 
Line 2,796:
(22.0 70.0 86.0 62.0)
(54.0 86.0 174.0 134.0)
(42.0 62.0 134.0 106.0) ) )</langsyntaxhighlight>
{{out}}
<pre> 5.00000 0.00000 0.00000
Line 2,808:
 
=={{header|PL/I}}==
<langsyntaxhighlight lang=PL/I>(subscriptrange):
decompose: procedure options (main); /* 31 October 2013 */
declare a(*,*) float controlled;
Line 2,854:
end cholesky;
 
end decompose;</langsyntaxhighlight>
ACTUAL RESULTS:-
<pre>Original matrix:
Line 2,877:
 
=={{header|PowerShell}}==
<langsyntaxhighlight lang=PowerShell>
function cholesky ($a) {
$l = @()
Line 2,929:
"l2 ="
show (cholesky $a2)
</syntaxhighlight>
</lang>
<b>Output:</b>
<pre>
Line 2,957:
=={{header|Python}}==
===Python2.X version===
<langsyntaxhighlight lang=python>from __future__ import print_function
 
from pprint import pprint
Line 2,983:
[54, 86, 174, 134],
[42, 62, 134, 106]]
pprint(cholesky(m2), width=120)</langsyntaxhighlight>
 
{{out}}
Line 2,997:
Factors out accesses to <code>A[i], L[i], and L[j]</code> by creating <code>Ai, Li and Lj</code> respectively as well as using <code>enumerate</code> instead of <code>range(len(some_array))</code>.
 
<langsyntaxhighlight lang=python>def cholesky(A):
L = [[0.0] * len(A) for _ in range(len(A))]
for i, (Ai, Li) in enumerate(zip(A, L)):
Line 3,004:
Li[j] = sqrt(Ai[i] - s) if (i == j) else \
(1.0 / Lj[j] * (Ai[j] - s))
return L</langsyntaxhighlight>
 
{{out}}
Line 3,011:
=={{header|q}}==
 
<langsyntaxhighlight lang=q>solve:{[A;B] $[0h>type A;B%A;inv[A] mmu B]}
ak:{[m;k] (),/:m[;k]til k:k-1}
akk:{[m;k] m[k;k:k-1]}
Line 3,026:
show cholesky (25 15 -5f;15 18 0f;-5 0 11f)
-1"";
show cholesky (18 22 54 42f;22 70 86 62f;54 86 174 134f;42 62 134 106f)</langsyntaxhighlight>
 
{{out}}
Line 3,040:
 
=={{header|R}}==
<langsyntaxhighlight lang=r>t(chol(matrix(c(25, 15, -5, 15, 18, 0, -5, 0, 11), nrow=3, ncol=3)))
# [,1] [,2] [,3]
# [1,] 5 0 0
Line 3,051:
# [2,] 5.185450 6.565905 0.000000 0.000000
# [3,] 12.727922 3.046038 1.649742 0.000000
# [4,] 9.899495 1.624554 1.849711 1.392621</langsyntaxhighlight>
 
=={{header|Racket}}==
<langsyntaxhighlight lang=racket>
#lang racket
(require math)
Line 3,083:
[54 86 174 134]
[42 62 134 106]]))
</syntaxhighlight>
</lang>
Output:
<langsyntaxhighlight lang=racket>
'#(#(5 0 0)
#(3 3 0)
Line 3,094:
#( 9.899494936611665 1.6245538642137891 1.849711005231382 1.3926212476455924))
 
</syntaxhighlight>
</lang>
 
=={{header|Raku}}==
(formerly Perl 6)
{{works with|Rakudo|2015.12}}
<syntaxhighlight lang=raku perl6line>sub cholesky(@A) {
my @L = @A »*» 0;
for ^@A -> $i {
Line 3,121:
[54, 86, 174, 134],
[42, 62, 134, 106],
];</langsyntaxhighlight>
 
=={{header|REXX}}==
Line 3,127:
<br>REXX number normalization) &nbsp; can be removed from the line &nbsp; (number 40) &nbsp; which contains the source statement:
::::: &nbsp; <b> z=z &nbsp; right( format(@.row.col, , &nbsp; dPlaces) / 1, &nbsp; &nbsp; width) </b>
<langsyntaxhighlight lang=rexx>/*REXX program performs the Cholesky decomposition on a square matrix & displays results*/
niner = '25 15 -5' , /*define a 3x3 matrix with elements. */
'15 18 0' ,
Line 3,175:
numeric form; m.=9; parse value format(x,2,1,,0) 'E0' with g 'E' _ .; g=g*.5'e'_ %2
do j=0 while h>9; m.j=h; h=h%2+1; end /*j*/
do k=j+5 to 0 by -1; numeric digits m.k; g=(g+x/g)*.5; end /*k*/; return g/1</langsyntaxhighlight>
{{out|output}}
<pre>
Line 3,208:
 
=={{header|Ring}}==
<langsyntaxhighlight lang=ring>
# Project : Cholesky decomposition
 
Line 3,251:
see nl
next
</syntaxhighlight>
</lang>
Output:
<pre>
Line 3,265:
 
=={{header|Ruby}}==
<langsyntaxhighlight lang=ruby>require 'matrix'
 
class Matrix
Line 3,302:
[22, 70, 86, 62],
[54, 86, 174, 134],
[42, 62, 134, 106]].cholesky_factor</langsyntaxhighlight>
{{out}}
<pre>
Line 3,315:
 
{{trans|C}}
<langsyntaxhighlight lang=rust>fn cholesky(mat: Vec<f64>, n: usize) -> Vec<f64> {
let mut res = vec![0.0; mat.len()];
for i in 0..n {
Line 3,355:
show_matrix(res2, dimension);
}
</syntaxhighlight>
</lang>
 
{{out}}
Line 3,370:
 
=={{header|Scala}}==
<langsyntaxhighlight lang=scala>case class Matrix( val matrix:Array[Array[Double]] ) {
 
// Assuming matrix is positive-definite, symmetric and not empty...
Line 3,430:
(l2.matrix.flatten.zip(r2)).foreach{ case (result,test) =>
assert(math.round( result * 100000 ) * 0.00001 == math.round( test * 100000 ) * 0.00001)
}</langsyntaxhighlight>
 
=={{header|Scilab}}==
Line 3,436:
The Cholesky decomposition is builtin, and an upper triangular matrix is returned, such that $A=L^TL$.
 
<langsyntaxhighlight lang=scilab>a = [25 15 -5; 15 18 0; -5 0 11];
chol(a)
ans =
Line 3,454:
0. 6.5659052 3.0460385 1.6245539
0. 0. 1.6497422 1.849711
0. 0. 0. 1.3926212</langsyntaxhighlight>
 
=={{header|Seed7}}==
<langsyntaxhighlight lang=seed7>$ include "seed7_05.s7i";
include "float.s7i";
include "math.s7i";
Line 3,516:
writeln;
writeMat(cholesky(m2));
end func;</langsyntaxhighlight>
 
Output:
Line 3,532:
=={{header|Sidef}}==
{{trans|Perl}}
<langsyntaxhighlight lang=ruby>func cholesky(matrix) {
var chol = matrix.len.of { matrix.len.of(0) }
for row in ^matrix {
Line 3,544:
}
return chol
}</langsyntaxhighlight>
 
Examples:
<langsyntaxhighlight lang=ruby>var example1 = [ [ 25, 15, -5 ],
[ 15, 18, 0 ],
[ -5, 0, 11 ] ];
Line 3,564:
cholesky(example2).each { |row|
say row.map {'%7.4f' % _}.join(' ');
}</langsyntaxhighlight>
{{out}}
<pre>
Line 3,580:
 
=={{header|Smalltalk}}==
<langsyntaxhighlight lang=Smalltalk>
FloatMatrix>>#cholesky
| l |
Line 3,603:
l at: k @ i put: aki - partialSum * factor reciprocal]]].
^l
</syntaxhighlight>
</lang>
 
=={{header|Stata}}==
See [http://www.stata.com/help.cgi?mf_cholesky Cholesky square-root decomposition] in Stata help.
<langsyntaxhighlight lang=stata>mata
: a=25,15,-5\15,18,0\-5,0,11
 
Line 3,646:
3 | 12.72792206 3.046038495 1.649742248 0 |
4 | 9.899494937 1.624553864 1.849711005 1.392621248 |
+---------------------------------------------------------+</langsyntaxhighlight>
 
=={{header|Swift}}==
{{trans|Rust}}
 
<langsyntaxhighlight lang=swift>func cholesky(matrix: [Double], n: Int) -> [Double] {
var res = [Double](repeating: 0, count: matrix.count)
 
Line 3,700:
printMatrix(res1, n: 3)
print()
printMatrix(res2, n: 4)</langsyntaxhighlight>
 
{{out}}
Line 3,714:
=={{header|Tcl}}==
{{trans|Java}}
<langsyntaxhighlight lang=tcl>proc cholesky a {
set m [llength $a]
set n [llength [lindex $a 0]]
Line 3,732:
}
return $l
}</langsyntaxhighlight>
Demonstration code:
<langsyntaxhighlight lang=tcl>set test1 {
{25 15 -5}
{15 18 0}
Line 3,746:
{42 62 134 106}
}
puts [cholesky $test2]</langsyntaxhighlight>
{{out}}
<pre>
Line 3,756:
This function returns the lower Cholesky decomposition of a square matrix fed to it. It does not check for positive semi-definiteness, although it does check for squareness. It assumes that <code>Option Base 0</code> is set, and thus the matrix entry indices need to be adjusted if Base is set to 1. It also assumes a matrix of size less than 256x256. To handle larger matrices, change all <code>Byte</code>-type variables to <code>Long</code>. It takes the square matrix range as an input, and can be implemented as an array function on the same sized square range of cells as output. For example, if the matrix is in cells A1:E5, highlighting cells A10:E14, typing "<code>=Cholesky(A1:E5)</code>" and htting <code>Ctrl-Shift-Enter</code> will populate the target cells with the lower Cholesky decomposition.
 
<langsyntaxhighlight lang=vb>Function Cholesky(Mat As Range) As Variant
 
Dim A() As Double, L() As Double, sum As Double, sum2 As Double
Line 3,808:
Cholesky = L
End Function
</syntaxhighlight>
</lang>
 
=={{header|Vlang}}==
{{trans|go}}
<langsyntaxhighlight lang=vlang>import math
 
// Symmetric and Lower use a packed representation that stores only
Line 3,910:
println("L:")
a.cholesky_lower().print()
}</langsyntaxhighlight>
 
{{out}}
Line 3,937:
{{libheader|Wren-matrix}}
{{libheader|Wren-fmt}}
<langsyntaxhighlight lang=ecmascript>import "/matrix" for Matrix
import "/fmt" for Fmt
 
Line 3,957:
Fmt.mprint(Matrix.new(array).cholesky(), 8, 5)
System.print()
}</langsyntaxhighlight>
 
{{out}}
Line 3,986:
=={{header|zkl}}==
Using the GNU Scientific Library:
<langsyntaxhighlight lang=zkl>var [const] GSL=Import("zklGSL"); // libGSL (GNU Scientific Library)
fcn lowerCholesky(m){ // trans: C
rows:=m.rows;
Line 3,996:
}
lcm
}</langsyntaxhighlight>
{{out}}
<pre>
Line 4,022:
Or, using lists:
{{trans|C}}
<langsyntaxhighlight lang=zkl>fcn cholesky(mat){
rows:=mat.len();
r:=(0).pump(rows,List().write, (0).pump(rows,List,0.0).copy); // matrix of zeros
Line 4,031:
}
r
}</langsyntaxhighlight>
<langsyntaxhighlight lang=zkl>ex1:=L( L(25.0,15.0,-5.0), L(15.0,18.0,0.0), L(-5.0,0.0,11.0) );
printM(cholesky(ex1));
println("-----------------");
Line 4,039:
L(54.0, 86.0, 174.0, 134.0,),
L(42.0, 62.0, 134.0, 106.0,) );
printM(cholesky(ex2));</langsyntaxhighlight>
<langsyntaxhighlight lang=zkl>fcn printM(m){ m.pump(Console.println,rowFmt) }
fcn rowFmt(row){ ("%9.5f "*row.len()).fmt(row.xplode()) }</langsyntaxhighlight>
{{out}}
<pre>
Line 4,056:
=={{header|ZX Spectrum Basic}}==
{{trans|BBC_BASIC}}
<langsyntaxhighlight lang=zxbasic>10 LET d=2000: GO SUB 1000: GO SUB 4000: GO SUB 5000
20 LET d=3000: GO SUB 1000: GO SUB 4000: GO SUB 5000
30 STOP
Line 4,090:
5050 PRINT
5060 NEXT r
5070 RETURN</langsyntaxhighlight>
Thundergnat
10,333

edits

Retrieved from "https://rosettacode.org/wiki/Special:MobileDiff/327217"