Voici une des 257 implémentations de l'algorithme de Gauss... Petite originalité : on utilise des listes dynamiques et on ne les parcourt que par «tête et queue».
C'est assez efficace...Quelques microsecondes pour le calcul du déterminant d'une matrice de taille 1000.
Haskell
det matrice taille 1000 : 592 microsec warming up estimating clock resolution... mean is 1.373373 us (640001 iterations) found 3269 outliers among 639999 samples (0.5%) 2863 (0.4%) high severe estimating cost of a clock call... mean is 33.83300 ns (9 iterations) found 1 outliers among 9 samples (11.1%) 1 (11.1%) high severe benchmarking detMat/1000 collecting 100 samples, 1 iterations each, in estimated 5.147409 s mean: 25.38584 ms, lb 25.30688 ms, ub 25.58190 ms, ci 0.950 std dev: 592.3872 us, lb 196.0864 us, ub 1.124442 ms, ci 0.950 found 8 outliers among 100 samples (8.0%) 5 (5.0%) high mild 3 (3.0%) high severe variance introduced by outliers: 17.070% variance is moderately inflated by outliers
Voici le code avec des noms les plus explicites possibles:
Haskell
-- déterminant rotatePivotDet :: (Eq a, Num a) => Matrice a -> Matrice a rotatePivotDet (Mat (ligne:lignes)) | ((head ligne) /= 0) || (all (== 0) (head (List.transpose (ligne:lignes)))) = Mat (ligne:lignes) | otherwise = rotatePivotDet (Mat (lignes ++ [map (* ((-1) ^ (length lignes))) ligne])) determinant :: (Eq a, Fractional a) => Matrice a -> a determinant (Mat []) = error "Matrice vide" determinant (Mat [[a]]) = a determinant m = if pivot == 0 then 0 else pivot * (determinant (Mat lignes')) where Mat ((pivot:l):lignes) = rotatePivotDet m lignes' = map transVec lignes transVec (t:q) | t == 0 = q | otherwise = fmap (/ pivot) (tail $ zipWith (-) (map (* pivot) (t:q)) (map (* t) (pivot:l))) -- inverse rotatePivot :: (Eq a, Num a) => [[a]] -> [[ a]] rotatePivot (ligne:lignes) | ((head ligne) /= 0) || (all (== 0) (head (List.transpose (ligne:lignes)))) = ligne:lignes | otherwise = rotatePivot (lignes ++ [ligne]) triangle :: (Eq a, Num a) => [[a]] -> [[a]] triangle [] = [] triangle m = (pivot:l):(triangle lignes') where (pivot:l):lignes = rotatePivot m lignes' = map transVec lignes transVec (t:q) | t == 0 = q | otherwise = tail $ zipWith (-) (map (* pivot) (t:q)) (map (* t) (pivot:l)) completeLigne :: (Num a) => Int -> [a] -> [a] completeLigne t l = [0 | x <- [0..(t - (length l) - 1)]] ++ l complete :: (Num a) => [[a]] -> [[a]] complete m = map (completeLigne (2 * (length m))) m retourne_complet :: [[a]] -> [[a]] retourne_complet m = reverse (map reverse m) unitet :: (Num a) => Int -> [[a]] unitet n = [[if i == j then 1 else 0 | j <- [0..n-1]] | i <- [0..n-1]] collet :: (Num a) => [[a]] -> [[a]] -> [[a]] collet m n = zipWith (++) m n tab :: (Num a) => [[a]] -> [[a]] tab m = collet m (unitet (length m)) carre_g :: [[a]] -> [[a]] carre_g m = map (take long) m where long = length m carre_d :: [[a]] -> [[a]] carre_d m = map (drop long) m where long = length m nm :: (Num a, Eq a) => [[a]] -> [[a]] nm a = collet (retourne_complet m1) (reverse m2) where (m1,m2) = (carre_g (m a),carre_d (m a)) where m a = complete $ triangle $ tab a nnm :: (Fractional a, Eq a) => [[a]] -> [[a]] nnm a = map (\ligne -> (map (\x -> x / (head ligne))) ligne) (triangle (nm a)) inverse :: (Fractional a, Eq a) => Matrice a -> Matrice a inverse (Mat a) = if (determinant (Mat a) == 0) then error "Matrice non inversible" else Mat (reverse $ carre_d $ complete (nnm a))