egs: add more timers

3 files changed, 69 insertions, 12 deletions

diff --git a/ell_grid_solve.c b/ell_grid_solve.c
index 4ad671a..5b4fade 100644
--- a/ell_grid_solve.c
+++ b/ell_grid_solve.c
@@ -183,7 +183,7 @@ static void boundaries_apply(EGSContext *ctx, int init)
     static const enum MG2DBCType bnd_type_order[] = { MG2D_BC_TYPE_FIXVAL, MG2D_BC_TYPE_REFLECT, MG2D_BC_TYPE_FALLOFF };
     EGSInternal *priv = ctx->priv;
     const ptrdiff_t strides[2] = { 1, priv->stride };
-    int64_t start;
+    int64_t start, start_bnd;
 
     start = gettime();
     for (int order_idx = 0; order_idx < ARRAY_ELEMS(bnd_type_order); order_idx++) {
@@ -205,23 +205,33 @@ static void boundaries_apply(EGSContext *ctx, int init)
         switch (bnd->type) {
         case MG2D_BC_TYPE_FIXVAL:
             if (init && !priv->bnd_zero[i]) {
+                start_bnd = gettime();
                 boundaries_apply_fixval(dst, dst_strides, bnd->val,
                                         bnd->val_stride, size_boundary);
+                ctx->time_bnd_fixval += gettime() - start_bnd;
+                ctx->count_bnd_fixval++;
             }
             break;
         case MG2D_BC_TYPE_REFLECT:
+            start_bnd = gettime();
             boundaries_apply_reflect(dst, dst_strides, bnd->val,
                                      bnd->val_stride, size_boundary);
+            ctx->time_bnd_reflect += gettime() - start_bnd;
+            ctx->count_bnd_reflect++;
             break;
         case MG2D_BC_TYPE_FALLOFF:
+            start_bnd = gettime();
             boundaries_apply_falloff(dst, dst_strides, bnd->val,
                                      bnd->val_stride, size_boundary, ctx);
+            ctx->time_bnd_falloff += gettime() - start_bnd;
+            ctx->count_bnd_falloff++;
             break;
         }
     }
     }
 
     /* fill in the corner ghosts */
+    start_bnd = gettime();
     for (int pos_y = 0; pos_y < 2; pos_y++) {
         enum MG2DBoundaryLoc loc_y = mg2d_bnd_id(1, pos_y);
         MG2DBoundary        *bnd_y = ctx->boundaries[loc_y];
@@ -252,6 +262,9 @@ static void boundaries_apply(EGSContext *ctx, int init)
                 }
         }
     }
+    ctx->time_bnd_corners += gettime() - start_bnd;
+    ctx->count_bnd_corners++;
+
     ctx->time_boundaries += gettime() - start;
     ctx->count_boundaries++;
 }
@@ -527,6 +540,10 @@ static int solve_exact(EGSContext *ctx)
             }
         }
     }
+    ec->time_mat_construct += gettime() - start;
+    ec->count_mat_construct++;
+
+    start = gettime();
 
     for (int i = 0; i < e->N; i++)
         for (int j = i + 1; j < e->N; j++) {
@@ -534,16 +551,22 @@ static int solve_exact(EGSContext *ctx)
             e->mat[j * e->N + i] = e->mat[i * e->N + j];
             e->mat[i * e->N + j] = tmp;
         }
-    ec->time_mat_construct += gettime() - start;
-    ec->count_mat_construct++;
+    ec->time_mat_transpose += gettime() - start;
+    ec->count_mat_transpose++;
 
     start = gettime();
 
     ret = mg2di_bicgstab_solve(e->bicgstab, e->mat, e->rhs, e->x);
+
+    ec->time_bicgstab_solve += gettime() - start;
+    ec->count_bicgstab_solve++;
+
     if (ret < 0) {
         char equed = 'N';
         double cond, ferr, berr, rpivot;
 
+        start = gettime();
+
         ret = LAPACKE_dgesvx(LAPACK_COL_MAJOR, 'N', 'N', e->N, 1,
                              e->mat, e->N, e->mat_f, e->N, e->ipiv, &equed, NULL, NULL,
                              e->rhs, e->N, e->x, e->N, &cond, &ferr, &berr, &rpivot);
@@ -559,16 +582,21 @@ static int solve_exact(EGSContext *ctx)
                   "condition number %16.16g; forward error %16.16g backward error %16.16g\n",
                   e->N, e->N, cond, ferr, berr);
 
+        ec->time_lu_solve += gettime() - start;
+        ec->count_lu_solve++;
+
         ret = mg2di_bicgstab_init(e->bicgstab, e->mat_f, e->x);
         if (ret < 0)
             return ret;
     }
 
+    start = gettime();
+
     for (size_t idx1 = 0; idx1 < ctx->domain_size[1]; idx1++)
         memcpy(ctx->u->data + idx1 * ctx->u->stride[0], e->x + idx1 * ctx->domain_size[0], ctx->domain_size[0] * sizeof(*e->x));
 
-    ec->time_lin_solve += gettime() - start;
-    ec->count_lin_solve++;
+    ec->time_export += gettime() - start;
+    ec->count_export++;
 
     boundaries_apply(ctx, 0);
     residual_calc(ctx);
@@ -597,11 +625,13 @@ int mg2di_egs_solve(EGSContext *ctx)
 int mg2di_egs_init(EGSContext *ctx)
 {
     EGSInternal *priv = ctx->priv;
-    int64_t start;
+    int64_t start, start_init;
     int ret;
 
     start = gettime();
 
+    start_init = gettime();
+
     if (ctx->solver_type == EGS_SOLVER_EXACT) {
         switch (ctx->fd_stencil) {
         case 1: priv->e.fill_mat = fill_mat_s1; break;
@@ -680,6 +710,9 @@ int mg2di_egs_init(EGSContext *ctx)
     if (ret < 0)
         return ret;
 
+    ctx->time_init += gettime() - start_init;
+    ctx->count_init++;
+
     boundaries_apply(ctx, 1);
     residual_calc(ctx);
 
diff --git a/ell_grid_solve.h b/ell_grid_solve.h
index d56616b..abe3958 100644
--- a/ell_grid_solve.h
+++ b/ell_grid_solve.h
@@ -91,8 +91,14 @@ typedef struct EGSRelaxContext {
 typedef struct EGSExactContext {
     int64_t count_mat_construct;
     int64_t time_mat_construct;
-    int64_t count_lin_solve;
-    int64_t time_lin_solve;
+    int64_t count_mat_transpose;
+    int64_t time_mat_transpose;
+    int64_t count_bicgstab_solve;
+    int64_t time_bicgstab_solve;
+    int64_t count_lu_solve;
+    int64_t time_lu_solve;
+    int64_t count_export;
+    int64_t time_export;
 } EGSExactContext;
 
 typedef struct EGSContext {
@@ -192,8 +198,18 @@ typedef struct EGSContext {
     /* timings */
     int64_t time_boundaries;
     int64_t count_boundaries;
+    int64_t time_bnd_fixval;
+    int64_t count_bnd_fixval;
+    int64_t time_bnd_falloff;
+    int64_t count_bnd_falloff;
+    int64_t time_bnd_reflect;
+    int64_t count_bnd_reflect;
+    int64_t time_bnd_corners;
+    int64_t count_bnd_corners;
     int64_t time_res_calc;
     int64_t count_res;
+    int64_t time_init;
+    int64_t count_init;
     int64_t time_total;
 } EGSContext;
 
diff --git a/mg2d.c b/mg2d.c
index 3e02ac1..0ae651b 100644
--- a/mg2d.c
+++ b/mg2d.c
@@ -792,9 +792,14 @@ void mg2d_print_stats(MG2DContext *ctx, const char *prefix)
         }
 
         ret = snprintf(p, sizeof(buf) - (p - buf),
-                       "||%2.2f%% residual %2.2f%% boundaries",
+                       "||%2.2f%% init %2.2f%% residual %2.2f%% boundaries (%2.2f%% fixval %2.2f%% reflect %2.2f%% falloff %2.2f%% corners)",
+                       level->solver->time_init       * 100.0 / level->solver->time_total,
                        level->solver->time_res_calc   * 100.0 / level->solver->time_total,
-                       level->solver->time_boundaries * 100.0 / level->solver->time_total);
+                       level->solver->time_boundaries * 100.0 / level->solver->time_total,
+                       level->solver->time_bnd_fixval * 100.0 / level->solver->time_total,
+                       level->solver->time_bnd_reflect * 100.0 / level->solver->time_total,
+                       level->solver->time_bnd_falloff * 100.0 / level->solver->time_total,
+                       level->solver->time_bnd_corners * 100.0 / level->solver->time_total);
         if (ret > 0)
             p += ret;
 
@@ -806,9 +811,12 @@ void mg2d_print_stats(MG2DContext *ctx, const char *prefix)
                 p += ret;
         } else if (e) {
             ret = snprintf(p, sizeof(buf) - (p - buf),
-                           " %2.2f%% matrix construct %2.2f%% linear solve",
+                           " %2.2f%% const %2.2f%% transp %2.2f%% bicgstab %2.2f%% lu %2.2f%% export",
                            e->time_mat_construct * 100.0 / level->solver->time_total,
-                           e->time_lin_solve     * 100.0 / level->solver->time_total);
+                           e->time_mat_transpose     * 100.0 / level->solver->time_total,
+                           e->time_bicgstab_solve     * 100.0 / level->solver->time_total,
+                           e->time_lu_solve     * 100.0 / level->solver->time_total,
+                           e->time_export     * 100.0 / level->solver->time_total);
             if (ret > 0)
                 p += ret;
         }