From 2be5215745b27d2a2d7bf5cab1ff3ebe37be5bef Mon Sep 17 00:00:00 2001 From: Anton Khirnov Date: Wed, 9 Jan 2019 15:07:30 +0100 Subject: mg2d: ignore padding values in findmax() --- mg2d.c | 27 +++++++++++++++------------ 1 file changed, 15 insertions(+), 12 deletions(-) diff --git a/mg2d.c b/mg2d.c index 7f1ce4e..cc75a6e 100644 --- a/mg2d.c +++ b/mg2d.c @@ -66,14 +66,15 @@ struct MG2DInternal { double *boundaries_base[4]; }; -static double findmax(double *arr, size_t len) +static double findmax(double *arr, size_t size[2], ptrdiff_t stride) { double ret = 0.0; - for (size_t i = 0; i < len; i++) { - double val = fabs(*arr++); - if (val > ret) - ret = val; - } + 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; } @@ -449,6 +450,8 @@ static int threadpool_init(MG2DContext *ctx) int mg2d_solve(MG2DContext *ctx) { MG2DInternal *priv = ctx->priv; + MG2DLevel *root = priv->root; + EllRelaxContext *s_root = root->solver; int64_t time_start; double res_prev, res_cur; int ret; @@ -465,18 +468,18 @@ int mg2d_solve(MG2DContext *ctx) if (ret < 0) return ret; - ret = mg2di_ell_relax_init(priv->root->solver); + ret = mg2di_ell_relax_init(s_root); if (ret < 0) return ret; - res_prev = findmax(priv->root->solver->residual, - priv->root->solver->residual_stride * ctx->domain_size); + res_prev = findmax(s_root->residual, s_root->domain_size, + s_root->residual_stride); for (int i = 0; i < ctx->maxiter; i++) { - mg_solve_subgrid(ctx, priv->root); + mg_solve_subgrid(ctx, root); - res_cur = findmax(priv->root->solver->residual, - priv->root->solver->residual_stride * ctx->domain_size); + res_cur = findmax(s_root->residual, s_root->domain_size, + s_root->residual_stride); if (res_cur < ctx->tol) { mg2di_log(&priv->logger, MG2D_LOG_INFO, "converged on iteration %d, residual %g\n", -- cgit v1.2.3