diff --git a/code/richardson_method.py b/code/richardson_method.py index 88ae9890..2d963d11 100644 --- a/code/richardson_method.py +++ b/code/richardson_method.py @@ -33,11 +33,11 @@ class RichardsonMethod: def solve(self): x = self.x0[:] if RichardsonMethod.convergence_norm(self.A, self.omega, self.I) >= 1: - return "Richardson method for those values will NOT converge" + return RichardsonMethod.convergence_norm(self.A, self.omega, self.I), "Richardson method for those values will NOT converge", for iteration in range(self.max_iterations): Ax = LinearAlgebraUtils.matrix_vector_multiply(self.A, x) residual = LinearAlgebraUtils.vector_vector_subtraction(self.b, Ax) x = LinearAlgebraUtils.vector_vector_addition(x, LinearAlgebraUtils.scalar_matrix_multiply(self.omega, residual)) - return x + return x, 0 diff --git a/code/tests.py b/code/tests.py index 6e16b4d4..1bdf7961 100644 --- a/code/tests.py +++ b/code/tests.py @@ -23,13 +23,13 @@ def calculate_eigenvalues(A): return lambda_min, lambda_max -def numpy_additional_richardson_calculations(): +def calcualte_norm_from_matrix_numpy(A, n, max_iterations): lambda_min, lambda_max = calculate_eigenvalues(A) - print("eigenvalues: ", lambda_min, lambda_max, RichardsonMethod.calculate_eigenvalues(A, max_iterations)) omega = 2 / (lambda_min + lambda_max) - print("omega: ", omega, RichardsonMethod.calculate_omega(lambda_min, lambda_max)) I = np.eye(n) - print("norms: ", calculate_norm_numpy(I, omega, A), RichardsonMethod.convergence_norm(A, omega, I)) + return calculate_norm_numpy(I, omega, A) + + @pytest.mark.parametrize("n", [2, 3, 4, 5, 10, 20, 50, 100]) def test_richardson_vs_cg(n: int): @@ -38,20 +38,20 @@ def test_richardson_vs_cg(n: int): max_iterations=1000 A, b = MatrixGenerator.generate_random_matrix_and_vector(n) richardson_solver = RichardsonMethod(A, b, max_iterations, size=n, tol=1e-7) - solution_richardson = richardson_solver.solve() + solution_richardson, info_richardson = richardson_solver.solve() solution_cg, info = cg(A, b) if info == 0: # SciPy CG converged - assert_scipy_converged(solution_richardson, solution_cg, tolerance, A, b) + assert_scipy_converged(solution_richardson, info_richardson, solution_cg, tolerance, A, b, max_iterations, n) else: # SciPy CG did not converge - assert_scipy_not_converged(solution_richardson, A, b) + assert_scipy_not_converged(solution_richardson, info_richardson, A, b) -def assert_scipy_converged(solution_richardson, solution_cg, tolerance, A, b): - if solution_richardson == "Richardson method for those values will NOT converge": +def assert_scipy_converged(solution_richardson, info_richardson, solution_cg, tolerance, A, b, max_iterations, n): + if info_richardson == "Richardson method for those values will NOT converge": print("Richardson did not converge, while SciPy did") - print("Matrix A:\n", A) - print("Vector b:\n", b) + numpy_norm = calcualte_norm_from_matrix_numpy(A, n, max_iterations) + print("Numpy norm: ", numpy_norm, " Richardson norm: ", solution_richardson) assert False, "Richardson did not converge, while SciPy did" else: difference = np.linalg.norm(solution_richardson - solution_cg) @@ -65,8 +65,8 @@ def assert_scipy_converged(solution_richardson, solution_cg, tolerance, A, b): print("Vector b:\n", b) assert difference < tolerance, f"The solutions are different! Difference: {difference:.8f}" -def assert_scipy_not_converged(solution_richardson, A, b): - if solution_richardson == "Richardson method for those values will NOT converge": +def assert_scipy_not_converged(solution_richardson, info_richardson, A, b): + if info_richardson == "Richardson method for those values will NOT converge": print("Richardson and SciPy did not converge") else: print("Richardson converged while SciPy did not:", solution_richardson)