SHARE
TWEET

Untitled

a guest Nov 10th, 2019 99 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. import numpy as np
  2. import copy
  3.  
  4. def gauss(A, b):
  5.     ## преобразуем две матрицы в одну для упрощения подсчетов
  6.     matrix = np.array(copy.deepcopy(A), dtype = np.float64)
  7.     matrix = np.concatenate((matrix, b), axis=1)
  8.     X = np.zeros(matrix.shape[0], dtype = np.float64)
  9.     ## приведение к треугольному виду
  10.     n = matrix.shape[0]
  11.     for i in range(0, n):
  12.         ## находим ненулевой элемент, меняем строки и запоминаем его индекс
  13.         j = i
  14.         while matrix[j][i] == 0:
  15.             j += 1
  16.             if j == n:
  17.                 print("Матрица вырождена!!!\n")
  18.                 return 1
  19.         tmp = copy.deepcopy(matrix[i])
  20.         matrix[i] = matrix[j]
  21.         matrix[j] = tmp
  22.         ## обнуляем оставшиеся элементы в столбце
  23.         matrix[i] = matrix[i] / matrix[i][i]
  24.         for k in range(0, n):
  25.             if k == i:
  26.                 continue
  27.             matrix[k] -= matrix[i] * matrix[k][i]
  28.     ## вычисление решений
  29.     for k in range(matrix.shape[0]):
  30.         X[k] = matrix[k][matrix.shape[1] - 1]
  31.     return X
  32.  
  33. def gauss_main_elem(A, b):
  34.     ## преобразуем две матрицы в одну для упрощения подсчетов
  35.     matrix = np.array(copy.deepcopy(A), dtype = np.float64)
  36.     matrix = np.concatenate((matrix, b), axis=1)
  37.     X = np.zeros(matrix.shape[0], dtype = np.float64)
  38.     indexes = np.arange(matrix.shape[1] - 1, dtype = np.int)
  39.     ## приведение к треугольному виду
  40.     n = matrix.shape[0]
  41.     for i in range(0, n):
  42.         ## находим ненулевой элемент, меняем строки и запоминаем его индекс
  43.         j = i
  44.         while matrix[j][i] == 0:
  45.             j += 1
  46.         tmp = copy.deepcopy(matrix[i])
  47.         matrix[i] = matrix[j]
  48.         matrix[j] = tmp
  49.         ## находим максимальный по модулю элемент в строке и
  50.         ## переставляем его столбец в текущую позицию
  51.         j = np.argmax(np.abs(matrix[i])[:matrix.shape[1] - 1])
  52.         tmp_index = copy.deepcopy(indexes[i])
  53.         indexes[i] = indexes[j]
  54.         indexes[j] = tmp_index
  55.         tmp = copy.deepcopy(matrix.T[i])
  56.         matrix.T[i] = matrix.T[j]
  57.         matrix.T[j] = tmp.T
  58.         ## обнуляем оставшиеся элементы в столбце
  59.         matrix[i] = matrix[i] / matrix[i][i]
  60.         for k in range(0, n):
  61.             if k == i:
  62.                 continue
  63.             matrix[k] -= matrix[i] * matrix[k][i]
  64.     ## вычисление решений
  65.     for k, index in enumerate(indexes):
  66.         X[index] = matrix[k][matrix.shape[1] - 1]
  67.     return X
  68.  
  69. def determinate(A):
  70.     determinate = 1
  71.     ## преобразуем две матрицы в одну для упрощения подсчетов
  72.     matrix = np.array(copy.deepcopy(A), dtype = np.float64)
  73.     ## приведение к треугольному виду
  74.     n = matrix.shape[0]
  75.     for i in range(0, n):
  76.         ## находим ненулевой элемент, меняем строки и запоминаем его индекс
  77.         j = i
  78.         while matrix[j][i] == 0:
  79.             j += 1
  80.         if i != j:
  81.             determinate *= -1
  82.         tmp = copy.deepcopy(matrix[i])
  83.         matrix[i] = matrix[j]
  84.         matrix[j] = tmp
  85.         ## обнуляем оставшиеся элементы в столбце
  86.         determinate *= matrix[i][i]
  87.         matrix[i] = matrix[i] / matrix[i][i]
  88.         for k in range(0, n):
  89.             if k == i:
  90.                 continue
  91.             matrix[k] -= matrix[i] * matrix[k][i]
  92.     return determinate
  93.  
  94. def reverse(A):
  95.     I = np.zeros(A.shape, dtype = np.float64)
  96.     for i in range(I.shape[0]):
  97.         I[i][i] = 1
  98.     ## преобразуем две матрицы в одну для упрощения подсчетов
  99.     matrix = np.array(copy.deepcopy(A), dtype = np.float64)
  100.     matrix = np.concatenate((matrix, I), axis=1)
  101.     ## приведение к треугольному виду
  102.     n = matrix.shape[0]
  103.     for i in range(0, n):
  104.         ## находим ненулевой элемент, меняем строки и запоминаем его индекс
  105.         j = i
  106.         while matrix[j][i] == 0:
  107.             j += 1
  108.         tmp = copy.deepcopy(matrix[i])
  109.         matrix[i] = matrix[j]
  110.         matrix[j] = tmp
  111.         ## обнуляем оставшиеся элементы в столбце
  112.         matrix[i] = matrix[i] / matrix[i][i]
  113.         for k in range(0, n):
  114.             if k == i:
  115.                 continue
  116.             matrix[k] -= matrix[i] * matrix[k][i]
  117.     return np.array(matrix.T[matrix.shape[1] // 2 :].T)
  118.  
  119. def norm(A):
  120.     return np.sum(A ** 2) ** (1 / 2)
  121.  
  122. def cond(A):
  123.     return norm(A) * norm(reverse(A))
  124.  
  125. def check(A, b, oper):
  126. ## Для проверки функций я использую функции библиотеки для анализа данных numpy
  127.     if oper == 1:
  128.         return np.linalg.det(A)
  129.     if oper == 2:
  130.         return np.linalg.inv(A)
  131.     if oper == 3:
  132.         return np.linalg.solve(A, b).flatten()
  133.     if oper == 4:
  134.         return np.linalg.solve(A, b).flatten()
  135.     if oper == 5:
  136.         return np.linalg.cond(A, p = 'fro')
  137.  
  138. def main():
  139.     print("Insert the dimension of matrix")
  140.     n = int(input())
  141.    
  142.     print("\nType of input\n")
  143.     print("\t1 - File\n\t2 - Use formula")
  144.     typ = int(input())
  145.    
  146.     A = np.zeros(shape=(n, n), dtype=np.float64)
  147.     b = np.zeros(shape=(n, 1), dtype=np.float64)
  148.     if typ == 1:
  149.         print("\nEnter filename\n")
  150.         filename = input()
  151.         with open(filename) as f:
  152.             for i in range(n):
  153.                 line = list(f.readline().split())
  154.                 for j, elem in enumerate(line[:n]):
  155.                     A[i][j] = np.float64(elem)
  156.                 b[i][0] = np.float64(line[n])
  157.     elif typ == 2:
  158.         print("\nEnter parametr m\n")
  159.         m = np.float64(input())
  160.         for i in range(n):
  161.             for j in range(n):
  162.                 if i == j:
  163.                     A[i][j] = n + m ** 2 + j / m + i / n
  164.                 else:
  165.                     A[i][j] = (i + j) / (m + n)
  166.             b[i][0] = 200 + 50 * i
  167.     print("\nChoose operation\n")
  168.     print("\t1 - Find determinant\n\t2 - Find reverse matrix\n" +
  169.            "\t3 - Find Gauss solution\n\t4 - Find Gauss solution with main element\n" +
  170.            "\t5 - Find condition number")
  171.     operation = int(input())
  172.     for operation in range(1, 6):
  173.         if operation == 1:
  174.             print("detrminate =\n {}".format(determinate(A)))
  175.         if operation == 2:
  176.             print("reverse(A) =\n {0}".format(reverse(A)))
  177.         if operation == 3:
  178.             print("X =\n {0}".format(gauss(A, b)))
  179.         if operation == 4:
  180.             print("X =\n {0}".format(gauss_main_elem(A, b)))
  181.         if operation == 5:
  182.             print("cond =\n {0}".format(cond(A)))
  183.         print("check:\n {0}\n".format(check(A, b, operation)))
  184.  
  185. main()
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
Top