diff --git a/matrix_calculus/matrix_norms.py b/matrix_calculus/matrix_norms.py index 40d774b..300a64c 100644 --- a/matrix_calculus/matrix_norms.py +++ b/matrix_calculus/matrix_norms.py @@ -47,3 +47,37 @@ def p_norm(matrix: np.array, p: int, n_samples: int = 1000, cond_number=True): res = (max_val, condidion_number) if cond_number else max_val return res + + +def first_norm(matrix: np.array, cond_number=False): + columns_sums = np.sum(np.abs(matrix), axis=0) + m_1 = np.max(columns_sums) + + if cond_number: + try: + matrix_inv = np.linalg.inv(matrix) + norm_M_inv = first_norm(matrix_inv, cond_number=False) + condition_number = m_1 * norm_M_inv + except np.linalg.LinAlgError: + condition_number = np.inf + + return m_1, condition_number + + return m_1 + + +def inf_norm(matrix: np.array, cond_number=False): + rows_sums = np.sum(np.abs(matrix), axis=1) + m_inf = np.max(rows_sums) + + if cond_number: + try: + matrix_inv = np.linalg.inv(matrix) + norm_M_inv = inf_norm(matrix_inv, cond_number=False) + condition_number = m_inf * norm_M_inv + except np.linalg.LinAlgError: + condition_number = np.inf + + return m_inf, condition_number + + return m_inf diff --git a/reports/matrix_norms_report.md b/reports/matrix_norms_report.md index e651ddc..872fe88 100644 --- a/reports/matrix_norms_report.md +++ b/reports/matrix_norms_report.md @@ -2,27 +2,121 @@ *Autorzy: Maja Byrecka, Michał Kowalczyk* +## Norma 1 + +Norma indukowana ||M||_1 wskazuje na maksymalne wzmocnienie sumy bezwzględnych wartości wektora po przekształceniu przez macierz. + +Złożoność obliczeniowa: Zaledwie O(n^2) dla wyznaczenia samej normy dla macierzy n x n. Obliczenie współczynnika uwarunkowania podnosi złożoność do O(n^3). + +### Algorytm + +``` +1. Oblicz sumy wartości bezwzględnych elementów dla każdej kolumny macierzy M. +2. Znajdź maksymalną wartość z obliczonych sum i zapisz w m_1. +3. Jeśli wymagane jest obliczenie współczynnika uwarunkowania: + 3.1. Oblicz macierz odwrotną M_inv. + 3.2. Wywołaj first_norm dla M_inv (bez liczenia współczynnika). + 3.3. Współczynnik uwarunkowania = m_1 * norm_M_inv. +4. Zwróć m_1 (oraz współczynnik uwarunkowania jeśli obliczono). +``` + +### Kod + +```python +def first_norm(matrix: np.array, cond_number=False): + columns_sums = np.sum(np.abs(matrix), axis=0) + m_1 = np.max(columns_sums) + + if cond_number: + try: + matrix_inv = np.linalg.inv(matrix) + norm_M_inv = first_norm(matrix_inv, cond_number=False) + condition_number = m_1 * norm_M_inv + except np.linalg.LinAlgError: + condition_number = np.inf + + return m_1, condition_number + + return m_1 +``` + +### Wynik dla Macierzy M + +Dla macierzy M = [[4, 9, 2], [3, 5, 7], [8, 1, 6]], norma macierzowa M_1 wynosi `15`, a współczynnik uwarunkowania wynosi `5.333333333333333`. + ## Norma 2 +Norma indukowana ||M||_2 odpowiada największej wartości własnej macierzy. + +Złożoność obliczeniowa: Jest to operacja wymagająca wyznaczenia rozkładu według wartości własnych macierzy. Złożoność wynosi O(n^3). + ### Algorytm ``` -1. Znajdź wartości własne macierzy M -2. Zwróć maksymalną wartość własną +1. Znajdź wartości własne macierzy M. +2. Norma m_2 to największa wartość własna macierzy M. +3. Współczynnik uwarunkowania to stosunek największej wartości własnej do najmniejszej. +4. Zwróć m_2 oraz obliczony współczynnik uwarunkowania. ``` ### Kod -```python3 -def second_norm(matrix: np.array) -> float: +```python +def second_norm(matrix: np.array): singular_values = np.linalg.svd(matrix, compute_uv=False) - return np.max(singular_values) + s = np.linalg.svd(matrix, compute_uv=False) + m_2 = np.max(singular_values) + condition_number = s.max() / s.min() + + return m_2, condition_number ``` ### Wynik dla Macierzy M Dla macierzy M = [[4, 9, 2], [3, 5, 7], [8, 1, 6]], norma macierzowa wyznaczona w powyższy sposób, z parametrem `p=2` wynosiła `15.000000000000002`. +## Norma inf + +Norma indukowana nieskończoność ||M||_inf definiuje najgorszy możliwy przypadek wzmocnienia sygnału dla pojedynczej współrzędnej wyjściowej. + +Złożoność obliczeniowa: Podobnie jak w przypadku normy 1, wynosi O(n^2) dla sumowania i wyznaczenia normy oraz O(n^3) dla obliczenia współczynnika uwarunkowania. + +### Algorytm + +``` +1. Oblicz sumy wartości bezwzględnych elementów dla każdego wiersza macierzy M. +2. Znajdź maksymalną wartość z obliczonych sum i przypisz do m_inf. +3. Jeśli wymagane jest obliczenie współczynnika uwarunkowania: + 3.1. Oblicz macierz odwrotną M_inv. + 3.2. Wywołaj inf_norm dla M_inv (bez liczenia współczynnika) i wynik zapisz w norm_M_inv. + 3.3. Współczynnik uwarunkowania = m_inf * norm_M_inv. +4. Zwróć m_inf (oraz współczynnik uwarunkowania jeśli obliczono). +``` + +### Kod + +```python +def inf_norm(matrix: np.array, cond_number=False): + rows_sums = np.sum(np.abs(matrix), axis=1) + m_inf = np.max(rows_sums) + + if cond_number: + try: + matrix_inv = np.linalg.inv(matrix) + norm_M_inv = inf_norm(matrix_inv, cond_number=False) + condition_number = m_inf * norm_M_inv + except np.linalg.LinAlgError: + condition_number = np.inf + + return m_inf, condition_number + + return m_inf +``` + +### Wynik dla Macierzy M + +Dla przykładowej macierzy M = [[4, 4, 2], [3, 2, 7], [8, 8, 3]], norma macierzowa M_inf wynosi `19`, a współczynnik uwarunkowania wynosi `370.5`. + ## Dowolna Norma P Dowolną normę *p* można obliczyć jako: @@ -30,6 +124,8 @@ Dowolną normę *p* można obliczyć jako: ![Źródło: Wykład - Wartości i wektory własne & Dekompozycja na wartości osobliwe 2026, Maciej Paszyński](imgs/matrix_p_norm.png) *Źródło: Wykład - Wartości i wektory własne & Dekompozycja na wartości osobliwe 2026, Maciej Paszyński* +W ten sposób możliwe jest obliczenie przybliżonej wartości dowolnej normy p macierzy. Złożoność algorytmu to O(k * n^2), gdzie k to liczba próbek. + ### Pseudokod ``` @@ -42,16 +138,19 @@ Dowolną normę *p* można obliczyć jako: 2.5. matrix_p <- p-norma wektora matrix_x 2.6. Jeśli matrix_p>max_val: 2.6.1. max_val <- matrix_p -3. Zwróć max_val + 2.6.2. Jeśli liczymy uwarunkowanie: oblicz odwrotność M_inv, wywołaj p_norm dla M_inv i zaktualizuj condition_number. +3. Zwróć max_val lub krotkę kax_val, condition_number jeśli obliczono uwarunkowanie ``` ### Algorytm ```python -def p_norm(matrix: np.array, p: int, n_samples: int =1000): +def p_norm(matrix: np.array, p: int, n_samples: int = 1000, cond_number=True): n = matrix.shape[1] max_val = 0.0 + if cond_number: + condidion_number = np.inf # Search for max value of the norm calculated from matrix multiplied by a randomly generated vector for _ in range(n_samples): @@ -71,7 +170,16 @@ def p_norm(matrix: np.array, p: int, n_samples: int =1000): if val > max_val: max_val = val - return max_val + if cond_number: + try: + matrix_inv = np.linalg.inv(matrix) + norm_M_inv = p_norm(matrix_inv, p, n_samples, cond_number=False) + condidion_number = max_val * norm_M_inv + except np.linalg.LinAlgError: + condidion_number = np.inf + + res = (max_val, condidion_number) if cond_number else max_val + return res ``` ### Wynik dla Macierzy M diff --git a/reports/resources/matrix-norms/code.txt b/reports/resources/matrix-norms/code.txt index 3f78b2e..0886e53 100644 --- a/reports/resources/matrix-norms/code.txt +++ b/reports/resources/matrix-norms/code.txt @@ -4,9 +4,11 @@ from __future__ import annotations import numpy as np import pytest # noqa: F401 -from scipy.linalg import norm -from matrix_calculus.matrix_norms import second_norm, p_norm +from matrix_calculus.matrix_norms import first_norm +from matrix_calculus.matrix_norms import inf_norm +from matrix_calculus.matrix_norms import p_norm +from matrix_calculus.matrix_norms import second_norm def test_m_2_norm(): @@ -18,6 +20,7 @@ def test_m_2_norm(): m_2, condition_number = second_norm(matrix) assert np.allclose(m_2, 27.343, rtol=0.01) + def test_p_norm_2(): matrix = np.array([[4, 9, 2], [3, 5, 7], [8, 1, 6]]) m_2, condition_number = p_norm(matrix, p=2) @@ -28,11 +31,57 @@ def test_p_norm_2(): assert np.allclose(m_2, 27.343, rtol=0.01) +def test_m_1_norm(): + matrix1 = np.array([[4, 9, 2], [3, 5, 7], [8, 1, 6]]) + m_1, cond = first_norm(matrix1, cond_number=True) + + # Porównanie z NumPy + expected_m1 = np.linalg.norm(matrix1, ord=1) + expected_cond = np.linalg.cond(matrix1, p=1) + + assert np.allclose(m_1, expected_m1), f"Błąd normy m1: {m_1} != {expected_m1}" + assert np.allclose(cond, expected_cond), f"Błąd cond m1: {cond} != {expected_cond}" + + matrix2 = np.array([[1, 2, 3], [3, 2, 1], [10, 15, 20]]) + m_1, cond = first_norm(matrix2, cond_number=True) + + expected_m1 = np.linalg.norm(matrix2, ord=1) + expected_cond = np.linalg.cond(matrix2, p=1) + + assert np.allclose(m_1, expected_m1) + assert np.allclose(cond, expected_cond) + + +def test_m_inf_norm(): + matrix1 = np.array([[4, 9, 2], [3, 5, 7], [8, 1, 6]]) + m_inf, cond = inf_norm(matrix1, cond_number=True) + + expected_minf = np.linalg.norm(matrix1, ord=np.inf) + expected_cond = np.linalg.cond(matrix1, p=np.inf) + + assert np.allclose( + m_inf, + expected_minf, + ), f"Błąd normy inf: {m_inf} != {expected_minf}" + assert np.allclose(cond, expected_cond), f"Błąd cond inf: {cond} != {expected_cond}" + + matrix2 = np.array([[1, 2, 3], [3, 2, 1], [10, 15, 20]]) + m_inf, cond = inf_norm(matrix2, cond_number=True) + + expected_minf = np.linalg.norm(matrix2, ord=np.inf) + expected_cond = np.linalg.cond(matrix2, p=np.inf) + + assert np.allclose(m_inf, expected_minf) + assert np.allclose(cond, expected_cond) + + ========== matrix_norms.py ========== +from __future__ import annotations + import numpy as np -def second_norm(matrix: np.array) -> float: +def second_norm(matrix: np.array): singular_values = np.linalg.svd(matrix, compute_uv=False) s = np.linalg.svd(matrix, compute_uv=False) m_2 = np.max(singular_values) @@ -40,7 +89,8 @@ def second_norm(matrix: np.array) -> float: return m_2, condition_number -def p_norm(matrix: np.array, p: int, n_samples: int =1000, cond_number=True): + +def p_norm(matrix: np.array, p: int, n_samples: int = 1000, cond_number=True): n = matrix.shape[1] max_val = 0.0 @@ -68,10 +118,44 @@ def p_norm(matrix: np.array, p: int, n_samples: int =1000, cond_number=True): if cond_number: try: matrix_inv = np.linalg.inv(matrix) - norm_M_inv = p_norm(matrix_inv, p, n_samples, cond_number = False) + norm_M_inv = p_norm(matrix_inv, p, n_samples, cond_number=False) condidion_number = max_val * norm_M_inv except np.linalg.LinAlgError: condidion_number = np.inf res = (max_val, condidion_number) if cond_number else max_val return res + + +def first_norm(matrix: np.array, cond_number=False): + columns_sums = np.sum(np.abs(matrix), axis=0) + m_1 = np.max(columns_sums) + + if cond_number: + try: + matrix_inv = np.linalg.inv(matrix) + norm_M_inv = first_norm(matrix_inv, cond_number=False) + condition_number = m_1 * norm_M_inv + except np.linalg.LinAlgError: + condition_number = np.inf + + return m_1, condition_number + + return m_1 + + +def inf_norm(matrix: np.array, cond_number=False): + rows_sums = np.sum(np.abs(matrix), axis=1) + m_inf = np.max(rows_sums) + + if cond_number: + try: + matrix_inv = np.linalg.inv(matrix) + norm_M_inv = inf_norm(matrix_inv, cond_number=False) + condition_number = m_inf * norm_M_inv + except np.linalg.LinAlgError: + condition_number = np.inf + + return m_inf, condition_number + + return m_inf diff --git a/reports/resources/matrix-norms/report.pdf b/reports/resources/matrix-norms/report.pdf index 4b9f1dd..3f394c3 100644 Binary files a/reports/resources/matrix-norms/report.pdf and b/reports/resources/matrix-norms/report.pdf differ diff --git a/tests/matrix_norms_test.py b/tests/matrix_norms_test.py index ea0f039..c21d6b2 100644 --- a/tests/matrix_norms_test.py +++ b/tests/matrix_norms_test.py @@ -3,6 +3,8 @@ import numpy as np import pytest # noqa: F401 +from matrix_calculus.matrix_norms import first_norm +from matrix_calculus.matrix_norms import inf_norm from matrix_calculus.matrix_norms import p_norm from matrix_calculus.matrix_norms import second_norm @@ -25,3 +27,47 @@ def test_p_norm_2(): matrix = np.array([[1, 2, 3], [3, 2, 1], [10, 15, 20]]) m_2, condition_number = p_norm(matrix, p=2) assert np.allclose(m_2, 27.343, rtol=0.01) + + +def test_m_1_norm(): + matrix1 = np.array([[4, 9, 2], [3, 5, 7], [8, 1, 6]]) + m_1, cond = first_norm(matrix1, cond_number=True) + + # Porównanie z NumPy + expected_m1 = np.linalg.norm(matrix1, ord=1) + expected_cond = np.linalg.cond(matrix1, p=1) + + assert np.allclose(m_1, expected_m1), f"Błąd normy m1: {m_1} != {expected_m1}" + assert np.allclose(cond, expected_cond), f"Błąd cond m1: {cond} != {expected_cond}" + + matrix2 = np.array([[1, 2, 3], [3, 2, 1], [10, 15, 20]]) + m_1, cond = first_norm(matrix2, cond_number=True) + + expected_m1 = np.linalg.norm(matrix2, ord=1) + expected_cond = np.linalg.cond(matrix2, p=1) + + assert np.allclose(m_1, expected_m1) + assert np.allclose(cond, expected_cond) + + +def test_m_inf_norm(): + matrix1 = np.array([[4, 9, 2], [3, 5, 7], [8, 1, 6]]) + m_inf, cond = inf_norm(matrix1, cond_number=True) + + expected_minf = np.linalg.norm(matrix1, ord=np.inf) + expected_cond = np.linalg.cond(matrix1, p=np.inf) + + assert np.allclose( + m_inf, + expected_minf, + ), f"Błąd normy inf: {m_inf} != {expected_minf}" + assert np.allclose(cond, expected_cond), f"Błąd cond inf: {cond} != {expected_cond}" + + matrix2 = np.array([[1, 2, 3], [3, 2, 1], [10, 15, 20]]) + m_inf, cond = inf_norm(matrix2, cond_number=True) + + expected_minf = np.linalg.norm(matrix2, ord=np.inf) + expected_cond = np.linalg.cond(matrix2, p=np.inf) + + assert np.allclose(m_inf, expected_minf) + assert np.allclose(cond, expected_cond)