From 1da0e0109e88e38e61cb5719607c8398ea8d1365 Mon Sep 17 00:00:00 2001 From: Thomas Guillod Date: Thu, 21 Dec 2023 12:30:56 -0500 Subject: [PATCH] better residuum computation --- pypeec/lib_matrix/matrix_iter.py | 5 ++++- pypeec/lib_solver/equation_solver.py | 3 +-- 2 files changed, 5 insertions(+), 3 deletions(-) diff --git a/pypeec/lib_matrix/matrix_iter.py b/pypeec/lib_matrix/matrix_iter.py index a65a7758..63c2fbbc 100644 --- a/pypeec/lib_matrix/matrix_iter.py +++ b/pypeec/lib_matrix/matrix_iter.py @@ -168,6 +168,9 @@ def fct_callback(sol): # check for convergence status = flag == 0 + # get residuum + res = sys_op(sol)-rhs + # get the number of iterations n_sys_eval = sys_obj.get_n_eval() n_pcd_eval = pcd_obj.get_n_eval() @@ -177,7 +180,7 @@ def fct_callback(sol): # assign results alg = {"n_sys_eval": n_sys_eval, "n_pcd_eval": n_pcd_eval, "n_iter": n_iter, "conv": conv} - return status, alg, sol + return status, alg, sol, res def get_solve(sol_init, sys_op, pcd_op, rhs, fct_conv, iter_options): diff --git a/pypeec/lib_solver/equation_solver.py b/pypeec/lib_solver/equation_solver.py index a3176f7b..2cbe7b50 100644 --- a/pypeec/lib_solver/equation_solver.py +++ b/pypeec/lib_solver/equation_solver.py @@ -83,10 +83,9 @@ def fct_sys(sol): sys_op = sla.LinearOperator((n_dof, n_dof), matvec=fct_sys, dtype=NP_TYPES.COMPLEX) # call the solver - (status_solver, alg, sol) = matrix_iter.get_solve(sol_init, sys_op, pcd_op, rhs, fct_conv, iter_options) + (status_solver, alg, sol, res) = matrix_iter.get_solve(sol_init, sys_op, pcd_op, rhs, fct_conv, iter_options) # compute and check the residuum - res = sys_op(sol)-rhs res_rms = np.sqrt(np.mean(np.abs(res)**2)) # get status