вычислить бесконечную норму разности двух решений

Question

В следующем коде я смог:

1. Реализация исключения Гаусса без поворота для общей квадратной линейной системы.
2. Я проверил это, решив Ax = bгде A - это случайная матрица 100x100, а b - случайный вектор 100x1.
3. Я сравнил свое решение с решением, полученным с использованием numpy.linalg.solve

Однако в конечном итогеЗадача мне нужно вычислить бесконечность нормы разности между двумя решениями.Я знаю, что норма бесконечности - это наибольшая абсолютная сумма строк матрицы.Но как я могу сделать это, чтобы вычислить бесконечную норму разницы между двумя решениями, моим решением и numpy.linalg.solve.Если вам нужна помощь с этим!

import numpy as np
def GENP(A, b):
    '''
    Gaussian elimination with no pivoting.
    % input: A is an n x n nonsingular matrix
    %        b is an n x 1 vector
    % output: x is the solution of Ax=b.
    % post-condition: A and b have been modified. 
    '''
    n =  len(A)
    if b.size != n:
        raise ValueError("Invalid argument: incompatible sizes between A & b.", b.size, n)
    for pivot_row in range(n-1):
        for row in range(pivot_row+1, n):
            multiplier = A[row][pivot_row]/A[pivot_row][pivot_row]
            #the only one in this column since the rest are zero
            A[row][pivot_row] = multiplier
            for col in range(pivot_row + 1, n):
                A[row][col] = A[row][col] - multiplier*A[pivot_row][col]
            #Equation solution column
            b[row] = b[row] - multiplier*b[pivot_row]
    x = np.zeros(n)
    k = n-1
    x[k] = b[k]/A[k,k]
    while k >= 0:
        x[k] = (b[k] - np.dot(A[k,k+1:],x[k+1:]))/A[k,k]
        k = k-1
    return x
if __name__ == "__main__":
    A = np.round(np.random.rand(100, 100)*10)
    b =  np.round(np.random.rand(100)*10)
    print (GENP(np.copy(A), np.copy(b)))

Например, этот код дает следующий вывод для задачи 1, указанной выше:

[-6.61537666  0.95704368  1.30101768 -3.69577873 -2.51427519 -4.56927017
 -1.61201589  2.88242622  1.67836096  2.18145556  2.60831672  0.08055869
 -2.39347903  2.19672137 -0.91609732 -1.17994959 -3.87309152 -2.53330865
  5.97476318  3.74687301  5.38585146 -2.71597978  2.0034079  -0.35045844
  0.43988439 -2.2623829  -1.82137544  3.20545721 -4.98871738 -6.94378666
 -6.5076601   3.28448129  3.42318453 -1.63900434  4.70352047 -4.12289961
 -0.79514656  3.09744616  2.96397264  2.60408589  2.38707091  8.72909353
 -1.33584905  1.30879264 -0.28008339  0.93560728 -1.40591226  1.31004142
 -1.43422946  0.41875924  3.28412668  3.82169545  1.96675247  2.76094378
 -0.90069455  1.3641636  -0.60520103  3.4814196  -1.43076816  5.01222382
  0.19160657  2.23163261  2.42183726 -0.52941262 -7.35597457 -3.41685057
 -0.24359225 -5.33856181 -1.41741354 -0.35654736 -1.71158503 -2.24469314
 -3.26453092  1.0932765   1.58333208  0.15567584  0.02793548  1.59561909
  0.31732915 -1.00695954  3.41663177 -4.06869021  3.74388762 -0.82868155
  1.49789582 -1.63559124  0.2741194  -1.11709237  1.97177449  0.66410154
  0.48397714 -1.96241854  0.34975886  1.3317751   2.25763568 -6.80055066
 -0.65903682 -1.07105965 -0.40211347 -0.30507635]

, тогда для задачи два мой код дает следующее:

my_solution = GENP(np.copy(A), np.copy(b))
numpy_solution = np.linalg.solve(A, b)
print(numpy_solution)

в результате:

[-6.61537666  0.95704368  1.30101768 -3.69577873 -2.51427519 -4.56927017
-1.61201589  2.88242622  1.67836096  2.18145556  2.60831672  0.08055869
-2.39347903  2.19672137 -0.91609732 -1.17994959 -3.87309152 -2.53330865
 5.97476318  3.74687301  5.38585146 -2.71597978  2.0034079  -0.35045844
 0.43988439 -2.2623829  -1.82137544  3.20545721 -4.98871738 -6.94378666
-6.5076601   3.28448129  3.42318453 -1.63900434  4.70352047 -4.12289961
-0.79514656  3.09744616  2.96397264  2.60408589  2.38707091  8.72909353
-1.33584905  1.30879264 -0.28008339  0.93560728 -1.40591226  1.31004142
-1.43422946  0.41875924  3.28412668  3.82169545  1.96675247  2.76094378
 -0.90069455  1.3641636  -0.60520103  3.4814196  -1.43076816  5.01222382
 0.19160657  2.23163261  2.42183726 -0.52941262 -7.35597457 -3.41685057
-0.24359225 -5.33856181 -1.41741354 -0.35654736 -1.71158503 -2.24469314
-3.26453092  1.0932765   1.58333208  0.15567584  0.02793548  1.59561909
 0.31732915 -1.00695954  3.41663177 -4.06869021  3.74388762 -0.82868155
 1.49789582 -1.63559124  0.2741194  -1.11709237  1.97177449  0.66410154
 0.48397714 -1.96241854  0.34975886  1.3317751   2.25763568 -6.80055066
-0.65903682 -1.07105965 -0.40211347 -0.30507635]

наконец для задачи 3:

if np.allclose(my_solution, numpy_solution):
    print("These solutions agree")
else:
   print("These solutions do not agree")

, в результате:

These solutions agree

Answer 1

Если то, что вы хотите, является только нормой бесконечности для матрицы, она обычно должна выглядеть примерно так:

def inf_norm(matrix):
    return max(abs(row.sum()) for row in matrix)

Но так как ваши my_solution и numpy_solution являются только одномерными векторами, выможет изменить их (я предполагаю, что 100x1, что у вас в вашем примере) для использования с вышеуказанной функцией:

альтернатива 1:

def inf_norm(matrix):
    return max(abs(row.sum()) for row in matrix)

diff = my_solution - numpy_solution
inf_norm_result = inf_norm(diff.reshape((100, 1))

альтернатива 2:

Или, если вы знаете, что они всегда будут одномерными векторами, вы можете опустить сумму (поскольку все строки будут иметь длину 1) и вычислить ее напрямую:

abs(my_solution - numpy_solution).max()

альтернатива 3:

или как написано в numpy.linalg.norm (см. Ниже) документации:

max(sum(abs(my_solution - numpy_solution), axis=1))

альтернатива 4:

или используйте numpy.linalg.norm() (см .: https://docs.scipy.org/doc/numpy-1.15.1/reference/generated/numpy.linalg.norm.html):

np.linalg.norm(my_solution - numpy_solution, np.inf)