better residuum computation

otvam · Dec 21, 2023 · 1da0e01 · 1da0e01
1 parent abf523e
commit 1da0e01
Show file tree

Hide file tree

Showing 2 changed files with 5 additions and 3 deletions.
diff --git 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
@@ -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