From 472c95da5fbea917407ecc70709dbd65a2e22507 Mon Sep 17 00:00:00 2001 From: Anton Khirnov Date: Tue, 5 Feb 2019 13:27:46 +0100 Subject: Change the meaning of boundary location to make more sense. API bump. --- ell_grid_solve.c | 74 ++++++++++++++++++++++++++++---------------------------- libmg2d.v | 2 +- 2 files changed, 38 insertions(+), 38 deletions(-) diff --git a/ell_grid_solve.c b/ell_grid_solve.c index 13d1854..75960e7 100644 --- a/ell_grid_solve.c +++ b/ell_grid_solve.c @@ -82,19 +82,19 @@ static const struct { unsigned int is_upper; } boundary_def[] = { [MG2D_BOUNDARY_0L] = { - .stride_idx = 0, - .is_upper = 0 + .stride_idx = 1, + .is_upper = 0, }, [MG2D_BOUNDARY_0U] = { - .stride_idx = 0, + .stride_idx = 1, .is_upper = 1, }, [MG2D_BOUNDARY_1L] = { - .stride_idx = 1, - .is_upper = 0, + .stride_idx = 0, + .is_upper = 0 }, [MG2D_BOUNDARY_1U] = { - .stride_idx = 1, + .stride_idx = 0, .is_upper = 1, }, }; @@ -198,7 +198,7 @@ static void boundaries_apply(EGSContext *ctx) double *dst = ctx->u; int fact_x = -1, fact_z = -1; - if (ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXDIFF) + if (ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXDIFF) fact_z *= -1; else fact_x *= -1; @@ -210,12 +210,12 @@ static void boundaries_apply(EGSContext *ctx) dst[idx_dst] = dst[idx_src]; } } - if (ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXDIFF || - ctx->boundaries[MG2D_BOUNDARY_1U]->type == MG2D_BC_TYPE_FIXDIFF) { + if (ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXDIFF || + ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXDIFF) { double *dst = ctx->u + ctx->domain_size[0] - 1; int fact_x = 1, fact_z = -1; - if (ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXDIFF) + if (ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXDIFF) fact_z *= -1; else fact_x *= -1; @@ -227,12 +227,12 @@ static void boundaries_apply(EGSContext *ctx) dst[idx_dst] = dst[idx_src]; } } - if (ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXDIFF || - ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXDIFF) { + if (ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXDIFF || + ctx->boundaries[MG2D_BOUNDARY_1U]->type == MG2D_BC_TYPE_FIXDIFF) { double *dst = ctx->u + strides[1] * (ctx->domain_size[1] - 1); int fact_x = -1, fact_z = 1; - if (ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXDIFF) + if (ctx->boundaries[MG2D_BOUNDARY_1U]->type == MG2D_BC_TYPE_FIXDIFF) fact_z *= -1; else fact_x *= -1; @@ -249,7 +249,7 @@ static void boundaries_apply(EGSContext *ctx) double *dst = ctx->u + strides[1] * (ctx->domain_size[1] - 1) + ctx->domain_size[0] - 1; int fact_x = 1, fact_z = 1; - if (ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXDIFF) + if (ctx->boundaries[MG2D_BOUNDARY_1U]->type == MG2D_BC_TYPE_FIXDIFF) fact_z *= -1; else fact_x *= -1; @@ -397,12 +397,12 @@ static int solve_exact(EGSContext *ctx) memset(e->mat, 0, SQR(e->N) * sizeof(*e->mat)); for (ptrdiff_t idx1 = 0; idx1 < ctx->domain_size[1]; idx1++) { - const int boundary_x_lower = idx1 < ctx->fd_stencil; - const int boundary_x_upper = idx1 >= ctx->domain_size[1] - ctx->fd_stencil; + const int boundary_y_lower = idx1 < ctx->fd_stencil; + const int boundary_y_upper = idx1 >= ctx->domain_size[1] - ctx->fd_stencil; for (ptrdiff_t idx0 = 0; idx0 < ctx->domain_size[0]; idx0++) { - const int boundary_y_lower = idx0 < ctx->fd_stencil; - const int boundary_y_upper = idx0 >= ctx->domain_size[0] - ctx->fd_stencil; + const int boundary_x_lower = idx0 < ctx->fd_stencil; + const int boundary_x_upper = idx0 >= ctx->domain_size[0] - ctx->fd_stencil; const int boundary = boundary_y_lower || boundary_y_upper || boundary_x_lower || boundary_x_upper; @@ -430,10 +430,10 @@ static int solve_exact(EGSContext *ctx) e->rhs[mat_row_idx] = ctx->rhs[idx_src]; if (boundary) { - const int bnd_fixval = (idx1 == 0 && ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXVAL) || - (idx0 == 0 && ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXVAL) || - (idx1 == ctx->domain_size[1] - 1 && ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXVAL) || - (idx0 == ctx->domain_size[0] - 1 && ctx->boundaries[MG2D_BOUNDARY_1U]->type == MG2D_BC_TYPE_FIXVAL); + const int bnd_fixval = (idx0 == 0 && ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXVAL) || + (idx1 == 0 && ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXVAL) || + (idx0 == ctx->domain_size[1] - 1 && ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXVAL) || + (idx1 == ctx->domain_size[0] - 1 && ctx->boundaries[MG2D_BOUNDARY_1U]->type == MG2D_BC_TYPE_FIXVAL); double *mat_dst = e->mat + e->N * mat_row_idx; ptrdiff_t bnd_idx = (boundary_x_lower || boundary_x_upper) ? idx0 : idx1; @@ -447,8 +447,8 @@ static int solve_exact(EGSContext *ctx) memset(e->scratch_line, 0, e->N_ghosts * sizeof(*e->scratch_line)); } else { /* apply the boundary conditions to eliminate the ghostpoint values */ - // fixval - if (boundary_x_lower && ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXDIFF) { + // fixdiff + if (boundary_y_lower && ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXDIFF) { for (ptrdiff_t idx1_bnd = 0; idx1_bnd < ctx->fd_stencil; idx1_bnd++) for (ptrdiff_t idx0_col = -(ptrdiff_t)ctx->fd_stencil; idx0_col < (ptrdiff_t)(ctx->domain_size[0] + ctx->fd_stencil); idx0_col++) { const ptrdiff_t idx_ghost = -(idx1_bnd + 1) * row_stride + idx0_col; @@ -456,7 +456,7 @@ static int solve_exact(EGSContext *ctx) mat_row[idx_ghost] = 0.0; } } - if (boundary_x_upper && ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXDIFF) { + if (boundary_y_upper && ctx->boundaries[MG2D_BOUNDARY_1U]->type == MG2D_BC_TYPE_FIXDIFF) { for (ptrdiff_t idx1_bnd = 0; idx1_bnd < ctx->fd_stencil; idx1_bnd++) for (ptrdiff_t idx0_col = -(ptrdiff_t)ctx->fd_stencil; idx0_col < (ptrdiff_t)(ctx->domain_size[0] + ctx->fd_stencil); idx0_col++) { const ptrdiff_t idx_ghost = (ctx->domain_size[1] - 1 + idx1_bnd + 1) * row_stride + idx0_col; @@ -464,7 +464,7 @@ static int solve_exact(EGSContext *ctx) mat_row[idx_ghost] = 0.0; } } - if (boundary_y_lower && ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXDIFF) { + if (boundary_x_lower && ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXDIFF) { for (ptrdiff_t idx1_col = -(ptrdiff_t)ctx->fd_stencil; idx1_col < (ptrdiff_t)(ctx->domain_size[1] + ctx->fd_stencil); idx1_col++) for (ptrdiff_t idx0_bnd = 0; idx0_bnd < ctx->fd_stencil; idx0_bnd++) { const ptrdiff_t idx_ghost = idx1_col * row_stride - (idx0_bnd + 1); @@ -472,7 +472,7 @@ static int solve_exact(EGSContext *ctx) mat_row[idx_ghost] = 0.0; } } - if (boundary_y_upper && ctx->boundaries[MG2D_BOUNDARY_1U]->type == MG2D_BC_TYPE_FIXDIFF) { + if (boundary_x_upper && ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXDIFF) { for (ptrdiff_t idx1_col = -(ptrdiff_t)ctx->fd_stencil; idx1_col < (ptrdiff_t)(ctx->domain_size[1] + ctx->fd_stencil); idx1_col++) for (ptrdiff_t idx0_bnd = 0; idx0_bnd < ctx->fd_stencil; idx0_bnd++) { const ptrdiff_t idx_ghost = idx1_col * row_stride + ctx->domain_size[0] - 1 + idx0_bnd + 1; @@ -482,8 +482,8 @@ static int solve_exact(EGSContext *ctx) } // fixval - if (boundary_x_lower && ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXVAL) { - const MG2DBoundary *bnd = ctx->boundaries[MG2D_BOUNDARY_0L]; + if (boundary_y_lower && ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXVAL) { + const MG2DBoundary *bnd = ctx->boundaries[MG2D_BOUNDARY_1L]; for (ptrdiff_t idx1_bnd = 1; idx1_bnd < ctx->fd_stencil; idx1_bnd++) for (ptrdiff_t idx0_col = -((ptrdiff_t)ctx->fd_stencil - 1); idx0_col < (ptrdiff_t)(ctx->domain_size[0] + ctx->fd_stencil - 1); idx0_col++) { const ptrdiff_t idx_ghost = -idx1_bnd * row_stride + idx0_col; @@ -491,8 +491,8 @@ static int solve_exact(EGSContext *ctx) mat_row[idx_ghost] = 0.0; } } - if (boundary_x_upper && ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXVAL) { - const MG2DBoundary *bnd = ctx->boundaries[MG2D_BOUNDARY_0U]; + if (boundary_y_upper && ctx->boundaries[MG2D_BOUNDARY_1U]->type == MG2D_BC_TYPE_FIXVAL) { + const MG2DBoundary *bnd = ctx->boundaries[MG2D_BOUNDARY_1U]; for (ptrdiff_t idx1_bnd = 1; idx1_bnd < ctx->fd_stencil; idx1_bnd++) for (ptrdiff_t idx0_col = -((ptrdiff_t)ctx->fd_stencil - 1); idx0_col < (ptrdiff_t)(ctx->domain_size[0] + ctx->fd_stencil - 1); idx0_col++) { const ptrdiff_t idx_ghost = (ctx->domain_size[1] - 1 + idx1_bnd) * row_stride + idx0_col; @@ -500,8 +500,8 @@ static int solve_exact(EGSContext *ctx) mat_row[idx_ghost] = 0.0; } } - if (boundary_y_lower && ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXVAL) { - const MG2DBoundary *bnd = ctx->boundaries[MG2D_BOUNDARY_1L]; + if (boundary_x_lower && ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXVAL) { + const MG2DBoundary *bnd = ctx->boundaries[MG2D_BOUNDARY_0L]; for (ptrdiff_t idx1_col = -((ptrdiff_t)ctx->fd_stencil - 1); idx1_col < (ptrdiff_t)(ctx->domain_size[1] + ctx->fd_stencil - 1); idx1_col++) for (ptrdiff_t idx0_bnd = 1; idx0_bnd < ctx->fd_stencil; idx0_bnd++) { const ptrdiff_t idx_ghost = idx1_col * row_stride - idx0_bnd; @@ -509,8 +509,8 @@ static int solve_exact(EGSContext *ctx) mat_row[idx1_col * row_stride - idx0_bnd] = 0.0; } } - if (boundary_y_upper && ctx->boundaries[MG2D_BOUNDARY_1U]->type == MG2D_BC_TYPE_FIXVAL) { - const MG2DBoundary *bnd = ctx->boundaries[MG2D_BOUNDARY_1U]; + if (boundary_x_upper && ctx->boundaries[MG2D_BOUNDARY_0U]->type == MG2D_BC_TYPE_FIXVAL) { + const MG2DBoundary *bnd = ctx->boundaries[MG2D_BOUNDARY_0U]; for (ptrdiff_t idx1_col = -((ptrdiff_t)ctx->fd_stencil - 1); idx1_col < (ptrdiff_t)(ctx->domain_size[1] + ctx->fd_stencil - 1); idx1_col++) for (ptrdiff_t idx0_bnd = 1; idx0_bnd < ctx->fd_stencil; idx0_bnd++) { const ptrdiff_t idx_ghost = idx1_col * row_stride + ctx->domain_size[0] - 1 + idx0_bnd; @@ -630,9 +630,9 @@ int mg2di_egs_init(EGSContext *ctx) priv->residual_calc_offset = 0; priv->residual_calc_size[0] = ctx->domain_size[0]; priv->residual_calc_size[1] = ctx->domain_size[1]; - if (ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXVAL) - priv->residual_calc_offset += ctx->residual_stride; if (ctx->boundaries[MG2D_BOUNDARY_1L]->type == MG2D_BC_TYPE_FIXVAL) + priv->residual_calc_offset += ctx->residual_stride; + if (ctx->boundaries[MG2D_BOUNDARY_0L]->type == MG2D_BC_TYPE_FIXVAL) priv->residual_calc_offset++; for (int i = 0; i < ARRAY_ELEMS(ctx->boundaries); i++) { diff --git a/libmg2d.v b/libmg2d.v index f30680f..400c46f 100644 --- a/libmg2d.v +++ b/libmg2d.v @@ -1,4 +1,4 @@ -LIBMG2D_6 { +LIBMG2D_7 { global: mg2d_*; local: *; }; -- cgit v1.2.3