ell_relax: compute the residual norm in residual_calc()

It is cheap and avoids an extra step in mg2d.

1 files changed, 2 insertions, 16 deletions

diff --git a/mg2d.c b/mg2d.c
index f64e0ae..55b8ed9 100644
--- a/mg2d.c
+++ b/mg2d.c
@@ -66,18 +66,6 @@ struct MG2DInternal {
     double *boundaries_base[4];
 };
 
-static double findmax(double *arr, size_t size[2], ptrdiff_t stride)
-{
-    double ret = 0.0;
-    for (size_t y = 0; y < size[1]; y++)
-        for (size_t x = 0; x < size[0]; x++) {
-            double val = fabs(arr[y * stride + x]);
-            if (val > ret)
-                ret = val;
-        }
-    return ret;
-}
-
 static void log_callback(MG2DLogger *log, int level, const char *fmt, va_list vl)
 {
     MG2DContext *ctx = log->opaque;
@@ -473,14 +461,12 @@ int mg2d_solve(MG2DContext *ctx)
     if (ret < 0)
         return ret;
 
-    res_prev = findmax(s_root->residual, s_root->domain_size,
-                       s_root->residual_stride);
+    res_prev = s_root->residual_max;
 
     for (int i = 0; i < ctx->maxiter; i++) {
         mg_solve_subgrid(ctx, root);
 
-        res_cur = findmax(s_root->residual, s_root->domain_size,
-                          s_root->residual_stride);
+        res_cur = s_root->residual_max;
 
         if (res_cur < ctx->tol) {
             mg2di_log(&priv->logger, MG2D_LOG_INFO, "converged on iteration %d, residual %g\n",