Pasted

flatten [] = []
flatten (x : xs) = x ++ (flatten xs)
 
-- ### I/O ### ----------------------------------------------------------------
-- reads first n lines from stdin and returns them in a list
getLines :: (Ord t, Num t) => t -> IO [String]
getLines n 
    | n <= 0 = return []
    | otherwise = do
        x <- getLine
        xs <- getLines (n - 1)
        return (x:xs)
 
-- transforms a list of lists to a printable matrix string 
vec_to_str [] = "\n"
vec_to_str (x : xs) = (take 6 ((if x < 0 then "" else " ") ++ (show x) ++ [' ' | _ <- [1..]])) ++ "\t" ++ (vec_to_str xs)
mat_to_str [] = ""
mat_to_str (row : rows) = (vec_to_str row) ++ (mat_to_str rows)
 
-- given a string representing vector b, finds vector x s.t. Ax = b
get_solution_string str_b lower upper swaps = 
    let b = map read (words str_b) :: [Double]
        y = forward_sub 1 lower [] (apply_swaps b swaps)
        x = backward_sub upper y
    in "Ax = " ++ (show b) ++ "    =>    x = " ++ (show x) ++ "\n"
 
-- given a list of strings of vectors, returns a string of solutions 
get_all_solutions strings lower upper swaps = 
    flatten [get_solution_string b lower upper swaps | b <- strings]
-- ######################################################################### --
 
-- ### Linear Algebra ### -----------------------------------------------------
mul_vv [] [] = 0
mul_vv (x : xs) (y : ys) = x * y + (mul_vv xs ys)
 
div_vc [] _ = []
div_vc (x : xs) c = (x / c) : (div_vc xs c)
 
div_mc [] _ = []
div_mc (row : rows) c = (div_vc row c) : (div_mc rows c)
 
outer_product [] _ = []
outer_product (x : xs) ys = [x * y | y <- ys] : (outer_product xs ys)
 
sub_vv [] [] = []
sub_vv (x : xs) (y : ys) = (x - y) : (sub_vv xs ys)
 
sub_mm [] [] = []
sub_mm (row1 : mat1) (row2 : mat2) = (sub_vv row1 row2) : (sub_mm mat1 mat2)
-- ######################################################################### --
 
-- ### Mics ### ----------------------------------------------------------------
-- splits the matrix returning the first column and the rest of the matrix 
split [] = ([], [])
split ((x : xs) : rows) = 
    let (col, mat) = split rows
    in ((x : col), (xs : mat))
 
-- reverse of split 
merge [] [] = []
merge (x : xs) (row : mat) = ((x : row) : (merge xs mat))
 
-- swaps two elements of list at positions i1, j1
swap i1 j1 list = 
    if i1 == j1 then
        list
    else
        let i = minimum [i1, j1]
            j = maximum [i1, j1]
            first = take i list
            second = drop (i + 1) (take j list)
            third = drop (j + 1) list
        in first ++ [list !! j] ++ second ++ [list !! i] ++ third
 
apply_swaps list [] = list
apply_swaps list ((i, j) : swaps) = apply_swaps (swap i j list) swaps
-- ######################################################################### --
 
-- ### LUP decomposition ### ---------------------------------------------------
-- returns the index of the maximum element in the list
max_idx _ [] _ res = res
max_idx i (x : xs) cur_max res = 
    if (abs x) > cur_max then
        max_idx (i + 1) xs (abs x) i
    else
        max_idx (i + 1) xs cur_max res
 
-- pivots the matrix bringing the biggest element to the top left corner.
-- returns the new matrix and the permutation represented as a lsit of swaps.
pivot mat =
    let (row, rest) = split mat
        i = max_idx 0 row (-1) (-1)
    in ((swap 0 i mat), (0, i))
 
add_1 [] = []
add_1 ((i, j) : rest) = (((i + 1), (j + 1)) : (add_1 rest))
 
-- computes the LUP decomposition returning (L, U, P).
lup_dec [[x]] = ([[1]], [[x]], [])
lup_dec input_mat =
    let (((a : row) : rows), q) = pivot input_mat
        (col, mat) = split rows
        schur = sub_mm mat (div_mc (outer_product col row) a)
        (l1, u1, p1) = lup_dec schur
        (swap_i, swap_j) = head p1
        l_top = (1 : [0 | _ <- row])
        u_left = (a : [0 | _ <- col])
        l = (l_top : (merge (div_vc (apply_swaps col p1) a) l1))
        u = merge u_left (row : u1)
        p = (q : (add_1 p1))
    in (l, u, p)
    --in ((a : row) : rows, schur, [q])
-- ######################################################################### --
 
-- ### Determinant ### --------------------------------------------------------
determinant [] = 1
determinant ((x : row) : upper) = 
    let (col, mat) = split upper
    in x * (determinant mat)
-- ######################################################################### --
 
-- ### Linear equations ### ---------------------------------------------------
forward_sub _ [] xs _ = xs
forward_sub i (row : lower) xs (b : bs) = 
    let prefix = mul_vv [x | (x, j) <- (zip row [1..(i - 1)])] xs
        x = b - prefix
    in forward_sub (i + 1) lower (xs ++ [x]) bs
 
backward_sub [[x]] [b] = [b / x]
backward_sub ((x : row) : upper) (b : bs) =
    let (col, mat) = split upper
        xs = backward_sub mat bs
        prefix = mul_vv row xs
    in ((b - prefix) / x) : xs
-- ######################################################################### --
 
 
main :: IO()
main = do
    -- reading and parsing the size of the matrix
    n_s <- getLine
    let n = read n_s :: Integer
 
    -- reading rows of the matrix as a list of strings
    rows <- getLines n
 
    -- transforming the list of strings to a double array of Integers
    let mat = [row | row_s <- rows, let row = (map read (words row_s) :: [Double])]
 
    -- calculating the LU-decomposition
    let (lower, upper, p) = lup_dec mat
 
    ---- calculating the determinant using the decomposed matrices
    let det = determinant upper * (if ((length p) `div` 2) == 0 then 1 else -1)
 
    putStrLn $ "Input matrix: \n" ++ (mat_to_str mat) ++ "\n\n"
    putStrLn $ "Unit-Lower Triangular: \n" ++ (mat_to_str lower) ++ "\n"
    putStrLn $ "Upper Triangular: \n" ++ (mat_to_str upper) ++ "\n"
    putStrLn $ "Permutations: " ++ (show p) ++ "\n"
    putStrLn $ "Determinant: " ++ (show det) ++ "\n"
 
    -- reading the amount of queries for linear equation solutions
    k_s <- getLine
    let k = read k_s :: Integer
 
    -- reading queries as a list of strings
    queries <- getLines k
 
    let solutions = get_all_solutions queries lower upper p
 
    putStrLn $ "\n\nSolutions: \n\n" ++ solutions