import numpy as np
import copy

def gauss(A, b):
    ## преобразуем две матрицы в одну для упрощения подсчетов
    matrix = np.array(copy.deepcopy(A), dtype = np.float64)
    matrix = np.concatenate((matrix, b), axis=1)
    X = np.zeros(matrix.shape[0], dtype = np.float64)
    ## приведение к треугольному виду
    n = matrix.shape[0]
    for i in range(0, n):
        ## находим ненулевой элемент, меняем строки и запоминаем его индекс
        j = i
        while matrix[j][i] == 0:
            j += 1
            if j == n:
                print("Матрица вырождена!!!\n")
                return 1
        tmp = copy.deepcopy(matrix[i])
        matrix[i] = matrix[j]
        matrix[j] = tmp
        ## обнуляем оставшиеся элементы в столбце
        matrix[i] = matrix[i] / matrix[i][i]
        for k in range(0, n):
            if k == i:
                continue
            matrix[k] -= matrix[i] * matrix[k][i]
    ## вычисление решений
    for k in range(matrix.shape[0]):
        X[k] = matrix[k][matrix.shape[1] - 1]
    return X

def gauss_main_elem(A, b):
    ## преобразуем две матрицы в одну для упрощения подсчетов
    matrix = np.array(copy.deepcopy(A), dtype = np.float64)
    matrix = np.concatenate((matrix, b), axis=1)
    X = np.zeros(matrix.shape[0], dtype = np.float64)
    indexes = np.arange(matrix.shape[1] - 1, dtype = np.int)
    ## приведение к треугольному виду
    n = matrix.shape[0]
    for i in range(0, n):
        ## находим ненулевой элемент, меняем строки и запоминаем его индекс
        j = i
        while matrix[j][i] == 0:
            j += 1
        tmp = copy.deepcopy(matrix[i])
        matrix[i] = matrix[j]
        matrix[j] = tmp
        ## находим максимальный по модулю элемент в строке и
        ## переставляем его столбец в текущую позицию
        j = np.argmax(np.abs(matrix[i])[:matrix.shape[1] - 1])
        tmp_index = copy.deepcopy(indexes[i])
        indexes[i] = indexes[j]
        indexes[j] = tmp_index
        tmp = copy.deepcopy(matrix.T[i])
        matrix.T[i] = matrix.T[j]
        matrix.T[j] = tmp.T
        ## обнуляем оставшиеся элементы в столбце
        matrix[i] = matrix[i] / matrix[i][i]
        for k in range(0, n):
            if k == i:
                continue
            matrix[k] -= matrix[i] * matrix[k][i]
    ## вычисление решений
    for k, index in enumerate(indexes):
        X[index] = matrix[k][matrix.shape[1] - 1]
    return X

def determinate(A):
    determinate = 1
    ## преобразуем две матрицы в одну для упрощения подсчетов
    matrix = np.array(copy.deepcopy(A), dtype = np.float64)
    ## приведение к треугольному виду
    n = matrix.shape[0]
    for i in range(0, n):
        ## находим ненулевой элемент, меняем строки и запоминаем его индекс
        j = i
        while matrix[j][i] == 0:
            j += 1
        if i != j:
            determinate *= -1
        tmp = copy.deepcopy(matrix[i])
        matrix[i] = matrix[j]
        matrix[j] = tmp
        ## обнуляем оставшиеся элементы в столбце
        determinate *= matrix[i][i]
        matrix[i] = matrix[i] / matrix[i][i]
        for k in range(0, n):
            if k == i:
                continue
            matrix[k] -= matrix[i] * matrix[k][i]
    return determinate

def reverse(A):
    I = np.zeros(A.shape, dtype = np.float64)
    for i in range(I.shape[0]):
        I[i][i] = 1
    ## преобразуем две матрицы в одну для упрощения подсчетов
    matrix = np.array(copy.deepcopy(A), dtype = np.float64)
    matrix = np.concatenate((matrix, I), axis=1)
    ## приведение к треугольному виду
    n = matrix.shape[0]
    for i in range(0, n):
        ## находим ненулевой элемент, меняем строки и запоминаем его индекс
        j = i
        while matrix[j][i] == 0:
            j += 1
        tmp = copy.deepcopy(matrix[i])
        matrix[i] = matrix[j]
        matrix[j] = tmp
        ## обнуляем оставшиеся элементы в столбце
        matrix[i] = matrix[i] / matrix[i][i]
        for k in range(0, n):
            if k == i:
                continue
            matrix[k] -= matrix[i] * matrix[k][i]
    return np.array(matrix.T[matrix.shape[1] // 2 :].T)

def norm(A):
    return np.sum(A ** 2) ** (1 / 2)

def cond(A):
    return norm(A) * norm(reverse(A))

def check(A, b, oper):
## Для проверки функций я использую функции библиотеки для анализа данных numpy
    if oper == 1:
        return np.linalg.det(A)
    if oper == 2:
        return np.linalg.inv(A)
    if oper == 3:
        return np.linalg.solve(A, b).flatten()
    if oper == 4:
        return np.linalg.solve(A, b).flatten()
    if oper == 5:
        return np.linalg.cond(A, p = 'fro')

def main():
    print("Insert the dimension of matrix")
    n = int(input())

    print("\nType of input\n")
    print("\t1 - File\n\t2 - Use formula")
    typ = int(input())

    A = np.zeros(shape=(n, n), dtype=np.float64)
    b = np.zeros(shape=(n, 1), dtype=np.float64)
    if typ == 1:
        print("\nEnter filename\n")
        filename = input()
        with open(filename) as f:
            for i in range(n):
                line = list(f.readline().split())
                for j, elem in enumerate(line[:n]):
                    A[i][j] = np.float64(elem)
                b[i][0] = np.float64(line[n])
    elif typ == 2:
        print("\nEnter parametr m\n")
        m = np.float64(input())
        for i in range(n):
            for j in range(n):
                if i == j:
                    A[i][j] = n + m ** 2 + j / m + i / n
                else:
                    A[i][j] = (i + j) / (m + n)
            b[i][0] = 200 + 50 * i
    print("\nChoose operation\n")
    print("\t1 - Find determinant\n\t2 - Find reverse matrix\n" +
           "\t3 - Find Gauss solution\n\t4 - Find Gauss solution with main element\n" +
           "\t5 - Find condition number")
    operation = int(input())
    for operation in range(1, 6):
        if operation == 1:
            print("detrminate =\n {}".format(determinate(A)))
        if operation == 2:
            print("reverse(A) =\n {0}".format(reverse(A)))
        if operation == 3:
            print("X =\n {0}".format(gauss(A, b)))
        if operation == 4:
            print("X =\n {0}".format(gauss_main_elem(A, b)))
        if operation == 5:
            print("cond =\n {0}".format(cond(A)))
        print("check:\n {0}\n".format(check(A, b, operation)))

main()