# LU decomposition (Eiffel)

(Redirected from LU Decomposition (Eiffel))

This class makes use of Real matrix (Eiffel).

```<<lu_decomposition.e>>=
class LU_DECOMPOSITION[X -> REAL_MATRIX] inherit MATRIX_DECOMPOSITION[X]
rename
matrix1 as lower,
matrix2 as upper
redefine
recompose
end
creation
make
feature {ANY} -- extra properties
permutation : X
parity : REAL
feature {ANY} -- creation
make(a : X) is
require
a.rows = a.columns
do
create lower.make(a.rows, a.columns)
create upper.make_identity(a.rows)
create permutation.make_identity(a.rows)
parity := 1.0
decompose(a.deep_twin)
end
feature {ANY}
decompose(a : X) is
local
i, j, n : INTEGER
sum : REAL
element_type : REAL
do
pivot(a)
-- first column of lower triangular, first row of upper triangular
from
i := 1
until
i > a.rows
loop
lower.set_element(i, 1, a.element(i,1))
upper.set_element(1, i, a.element(1,i)/lower.element(1,1))
i := i + 1
end
-- caluculate the remaining portion
from
i := 2
until
i > a.rows
loop
from
j := i
until
j > a.rows
loop
from
n := 1
sum := element_type.zero
until
n > i - 1
loop
sum := sum + lower.element(j,n) * upper.element(n,i)
n := n + 1
end
lower.set_element(j, i, a.element(j,i) - sum)
j := j + 1
end
from
j := i + 1
until
j > a.rows
loop
from
n := 1
sum := element_type.zero
until
n > i - 1
loop
sum := sum + lower.element(i,n) * upper.element(n, j)
n := n + 1
end
upper.set_element(i, j, ((a.element(i,j) - sum)/lower.element(i,i)))
j := j + 1
end
i := i + 1
end
ensure then
permutation /= Void
end -- decompose
recompose : X is
do
Result := permutation * lower * upper
end
feature {NONE}
pivot(a : X) is
local
i, j, max_index : INTEGER
col_max : REAL
do
from
j := 1
until
j > a.columns
loop
from
i := j
col_max := 0
until
i > a.rows
loop
if a.element(i,j).abs >= col_max then
col_max := a.element(i,j).abs
max_index := i
end
i := i + 1
end
if j /= max_index then
a.swap_rows(j, max_index)
permutation.swap_rows(j, max_index)
parity := parity * -1.0
end
j := j + 1
end
end -- pivot
end -- class LU_DECOMPOSITION
```
```<<matrix_decomposition.e>>=
deferred class MATRIX_DECOMPOSITION[X -> MATRIX[NUMERIC]]
feature {ANY} -- storage
matrix1 : X
matrix2 : X
feature {ANY} -- (de|re)-composition functions
decompose(a : X) is
require
a.rows = a.columns
deferred
ensure
matrix1 /= Void
matrix2 /= Void
end -- decompose
recompose : X is
require
matrix1 /= Void
matrix2 /= Void
do
Result := matrix1 * matrix2
end -- recompose
end -- class MATRIX_DECOMPOSITION
```