Um pouco de teoria
Como você sabe, o determinante de uma matriz quadrada n * n é a soma de n! termos, cada um dos quais é um produto contendo exatamente um elemento de matriz de cada coluna e exatamente um de cada linha. O sinal da próxima peça:
determinado pela paridade da substituição:
\ begin {pmatrix} 1 & 2 &… & n \\ {i} _ {1} & {i} _ {2} &… & {i} _ {n} \ end {pmatrix}
Um método direto para calcular o determinante consiste em decompor-o pelos elementos de uma linha ou coluna na soma dos produtos dos elementos de uma linha ou coluna por seus complementos algébricos. Por sua vez, o complemento algébrico do elemento da matriz
Há sim
em que
- há um elemento menor (i, j), ou seja, determinante obtido do determinante original pela exclusão da i-ésima linha e da j-ésima coluna.
Este método gera um processo recursivo que permite que qualquer determinante seja calculado. Mas o desempenho desse algoritmo deixa muito a desejar - O (n!). Portanto, tal cálculo direto é usado, exceto para cálculos simbólicos (e com determinantes de ordem não muito alta).
O método de Gauss acabou sendo muito mais produtivo. Sua essência é baseada nas seguintes disposições:
1. Determinante da matriz triangular superior \ begin {pmatrix} {a} _ {1,1} & {a} _ {1,2} &… & {a} _ {1, n} \\ 0 & {a} _ {2,2} &… & {a} _ {2, n} \\ 0 & 0 & ... & ... \\ 0 & 0 &… & {a} _ {n, n} \\\ end { pmatrix} é igual ao produto de seus elementos diagonais. Este fato decorre imediatamente da decomposição do determinante em termos dos elementos da primeira linha ou da primeira coluna. 2. Se na matriz aos elementos de uma linha somar os elementos de outra linha, multiplicados pelo mesmo número, o valor do determinante não mudará. 3. Se duas linhas (ou duas colunas) forem trocadas na matriz, o valor do determinante mudará seu sinal para o oposto. igual ao produto de seus elementos da diagonal. Este fato decorre imediatamente da decomposição do determinante em termos dos elementos da primeira linha ou da primeira coluna.
Ao escolher os coeficientes, podemos somar a primeira linha da matriz com todas as outras e obter zeros na primeira coluna em todas as posições, exceto na primeira. Para obter zero na segunda linha, você precisa adicionar à segunda linha a primeira, multiplicada por
Para obter zero na terceira linha, você precisa adicionar a primeira linha multiplicada por
etc. Em última análise, a matriz será reduzida a uma forma em que todos os elementos
para n> 1 será igual a zero.
Se na matriz o elemento
acabar sendo igual a zero, então você pode encontrar um elemento diferente de zero na primeira coluna (suponha que ele esteja na ka posição) e trocar a primeira e a ka linha. Com essa transformação, o determinante simplesmente muda de sinal, que pode ser levado em consideração. Se não houver elementos diferentes de zero na primeira coluna, o determinante será igual a zero.
Além disso, agindo de maneira semelhante, você pode obter zeros na segunda coluna, depois na terceira, etc. É importante que os zeros obtidos anteriormente não mudem quando as strings forem adicionadas. Se para qualquer string não for possível encontrar um elemento diferente de zero para o denominador, então o determinante é zero e o processo pode ser interrompido. A conclusão normal do processo gaussiano gera uma matriz em que todos os elementos localizados abaixo da diagonal principal são iguais a zero. Conforme mencionado acima, o determinante de tal matriz é igual ao produto dos elementos diagonais.
Vamos passar para a programação.
Trabalhamos com dados de ponto flutuante. Representamos matrizes como listas de strings. Primeiro, vamos definir dois tipos:
type Row = [Double]
type Matrix = [Row]
Recursão simples
Sem mais hesitações, vamos expandir o determinante pelos elementos da primeira linha (ou seja, zero). Precisamos de um programa para construir o menor, obtido excluindo a primeira linha e a k-ésima coluna.
-- k-
deln :: Matrix -> Int -> Matrix
deln matrix k = map (\ r -> (take (k) r)++(drop (k+1) r)) matrix
E aqui está o menor:
-- k-
minor :: Matrix -> Int -> Double
minor matrix k = det $ deln (drop 1 matrix) k
Observação: menor é um fator determinante. Estamos chamando a função det, que ainda não implementamos. Para implementar det, teremos que formar a soma alternada dos produtos do próximo elemento da primeira linha pelo determinante do próximo menor. Para evitar expressões complicadas, vamos criar uma função separada para formar o sinal de soma:
sgn :: Int -> Double
sgn n = if n `rem` 2 == 0 then 1.0 else (-1.0)
Agora podemos calcular o determinante:
--
det :: Matrix -> Double
det [[a,b],[c,d]] = a*d-b*c
det matrix = sum $ map (\c -> ((matrix !! 0)!!c)*(sgn c)*(minor matrix c)) [0..n]
where n = length matrix - 1
O código é muito simples e não requer nenhum comentário especial. Para testar o desempenho de nossas funções, vamos escrever a função principal:
main = print $ det [[1,2,3],[4,5,6],[7,8,(-9)]]O valor deste determinante é 54, o que pode ser verificado.
Método de Gauss
Precisaremos de algumas funções utilitárias (que podem ser usadas em outro lugar). O primeiro é o intercâmbio de duas linhas na matriz:
--
swap :: Matrix -> Int -> Int -> Matrix
swap matrix n1 n2 = map row [0..n]
where n=length matrix - 1
row k | k==n1 = matrix !! n2
| k==n2 = matrix !! n1
| otherwise = matrix !! k
Como você pode ver no código acima, a função passa linha por linha. Nesse caso, se uma linha com o número n1 for encontrada, a linha n2 será substituída à força (e vice-versa). O resto das linhas permanecem no lugar.
A função a seguir calcula a string r1 adicionada com a string r2 multiplicada elemento por elemento pelo número f:
-- r1+f*r2
comb :: Row -> Row -> Double -> Row
comb r1 r2 f = zipWith (\ x y -> x+f*y) r1 r2
Tudo aqui é extremamente transparente: as ações são realizadas nas linhas da matriz (ou seja, nas listas [Duplas]). Mas a seguinte função realiza essa transformação na matriz (e, naturalmente, obtém uma nova matriz):
-- r1 r2, f
trans :: Matrix -> Int -> Int -> Double -> Matrix
trans matrix n1 n2 f = map row [0..n]
where n=length matrix - 1
row k | k==n1 = comb (matrix !! n1) (matrix !! n2) f
| otherwise = matrix !! k
A função getNz procura o número do primeiro elemento diferente de zero na lista. É necessário quando o próximo elemento da diagonal é igual a zero.
--
getNz :: Row -> Int
getNz xs = if length tmp == 0 then (-1) else snd $ head tmp
where tmp=dropWhile (\ (x,k) -> (abs x) <= 1.0e-10) $ zip xs [0..]
Se todos os elementos da lista forem zero, a função retornará -1.
A função de pesquisa verifica se a matriz é adequada para a próxima transformação (ela deve ter um próximo elemento diagonal diferente de zero). Se não for, a matriz é transformada por permutação de linhas.
--
search :: Matrix -> Int -> Matrix
search matrix k | (abs ((matrix !! k) !! k)) > 1.0e-10 = matrix
| nz < 0 = matrix --
| otherwise = swap matrix k p
where n = length matrix
lst = map (\ r -> r !! k) $ drop k matrix
nz = getNz lst
p = k + nz
Se for impossível encontrar o elemento principal (diferente de zero) (a matriz está degenerada), a função o retornará inalterado. A função mkzero forma zeros na próxima coluna da matriz:
--
mkzero :: Matrix -> Int -> Int -> Matrix
mkzero matrix k p | p>n-1 = matrix
| otherwise = mkzero (trans matrix p k (-f)) k (p+1)
where n = length matrix
f = ((matrix !! p) !! k)/((matrix !! k) !! k)
A função de triângulo forma a forma triangular superior da matriz:
--
triangle :: Matrix -> Int -> Matrix
triangle matrix k | k>=n = matrix
| (abs v) <= 1.0e-10 = [[0.0]] --
| otherwise = triangle (mkzero tmp k k1) k1
where n = length matrix
tmp = search matrix k
v = (tmp !! k) !! k --
k1 = k+1
Se na próxima etapa não for possível encontrar o elemento condutor, a função retorna uma matriz zero de 1ª ordem. Agora podemos compor a função cerimonial de reduzir a matriz à forma triangular superior:
--
gauss :: Matrix -> Matrix
gauss matrix = triangle matrix 0
Para calcular o determinante, precisamos multiplicar os elementos diagonais. Para fazer isso, vamos criar uma função separada:
--
proddiag :: Matrix -> Double
proddiag matrix = product $ map (\ (r,k) -> r !!k) $ zip matrix [0,1..]
Bem, e o "arco" é o cálculo real do determinante:
--
det :: Matrix -> Double
det matrix = proddiag $ triangle matrix 0
Vamos verificar como essa função funciona:
main = print $ det [[1,2,3],[4,5,6],[7,8,-9]]
[1 of 1] Compiling Main ( main.hs, main.o )
Linking a.out ...
54.0
Obrigado a quem leu até o fim!
O código pode ser baixado aqui