From ab72ad7bb19f08d78218d3558545f9f58e5b36e7 Mon Sep 17 00:00:00 2001 From: Anton Khirnov Date: Thu, 25 Jun 2020 20:51:02 +0200 Subject: Switch to external ndarray library. --- ell_grid_solve.c | 86 ++++++++++++++++++++++++++++---------------------------- 1 file changed, 43 insertions(+), 43 deletions(-) (limited to 'ell_grid_solve.c') diff --git a/ell_grid_solve.c b/ell_grid_solve.c index fd6e455..dad8861 100644 --- a/ell_grid_solve.c +++ b/ell_grid_solve.c @@ -27,6 +27,7 @@ #include #include +#include #include "bicgstab.h" #include "common.h" @@ -37,7 +38,6 @@ #include "mg2d_boundary.h" #include "mg2d_boundary_internal.h" #include "mg2d_constants.h" -#include "ndarray.h" #include "residual_calc.h" #include "timer.h" @@ -169,7 +169,7 @@ static void residual_calc(EGSContext *ctx, int export_res) mg2di_timer_start(&ctx->timer_res_calc); for (int i = 0; i < ARRAY_ELEMS(diff_coeffs); i++) - diff_coeffs[i] = NDPTR2D(priv->diff_coeffs[i], offset[0], offset[1]); + diff_coeffs[i] = NDA_PTR2D(priv->diff_coeffs[i], offset[0], offset[1]); for (int bnd_idx = 0; bnd_idx < 4; bnd_idx++) if (ctx->boundaries[bnd_idx]->type == MG2D_BC_TYPE_REFLECT && @@ -177,9 +177,9 @@ static void residual_calc(EGSContext *ctx, int export_res) reflect_flags |= 1 << bnd_idx; mg2di_residual_calc(priv->rescalc, size, &ctx->residual_max, - NDPTR2D(dst, offset[0], offset[1]), dst->stride[0], - NDPTR2D(ctx->u, offset[0], offset[1]), ctx->u->stride[0], - NDPTR2D(ctx->rhs, offset[0], offset[1]), ctx->rhs->stride[0], + NDA_PTR2D(dst, offset[0], offset[1]), dst->stride[0], + NDA_PTR2D(ctx->u, offset[0], offset[1]), ctx->u->stride[0], + NDA_PTR2D(ctx->rhs, offset[0], offset[1]), ctx->rhs->stride[0], diff_coeffs, priv->diff_coeffs[0]->stride[0], export_res ? 0.0 : 1.0, export_res ? 1.0 : priv->r.relax_factor, reflect_flags, FD_STENCIL_MAX); @@ -580,7 +580,7 @@ static void mat_fill_row(EGSContext *ctx, double *scratch_line, EGSInternal *priv = ctx->priv; EGSExactInternal *e = &priv->e; - const ptrdiff_t idx_src_dc = NDIDX2D(priv->diff_coeffs[0], idx0, idx1); + const ptrdiff_t idx_src_dc = NDA_IDX2D(priv->diff_coeffs[0], idx0, idx1); const ptrdiff_t mat_row_idx = idx1 * ctx->domain_size[0] + idx0; int is_bnd[4], boundary; @@ -785,7 +785,7 @@ static int solve_exact(EGSContext *ctx) for (ptrdiff_t idx1 = 0; idx1 < ctx->domain_size[1]; idx1++) for (ptrdiff_t idx0 = 0; idx0 < ctx->domain_size[0]; idx0++) { const ptrdiff_t mat_row_idx = idx1 * ctx->domain_size[0] + idx0; - e->rhs[mat_row_idx] = e->rhs_factor[mat_row_idx] * (*NDPTR2D(ctx->rhs, idx0, idx1)) + e->rhs_mat[mat_row_idx]; + e->rhs[mat_row_idx] = e->rhs_factor[mat_row_idx] * (*NDA_PTR2D(ctx->rhs, idx0, idx1)) + e->rhs_mat[mat_row_idx]; } mg2di_timer_stop(&ec->timer_mat_construct); @@ -865,7 +865,7 @@ static int init_diff_coeffs_task(void *arg, unsigned int job_idx, unsigned int t const double factor = a->fact; for (ptrdiff_t idx1 = 0; idx1 < dst->shape[1]; idx1++) - *NDPTR2D(dst, idx1, job_idx) = *NDPTR2D(src, idx1, job_idx) * factor; + *NDA_PTR2D(dst, idx1, job_idx) = *NDA_PTR2D(src, idx1, job_idx) * factor; return 0; } @@ -1126,8 +1126,8 @@ static int arrays_alloc(EGSContext *ctx, const size_t domain_size[2]) EGSInternal *priv = ctx->priv; const DomainComponent *dc = &priv->dg->components[priv->local_component]; - Slice slice_interior[2]; - Slice slice_exterior[2]; + NDASlice slice_interior[2]; + NDASlice slice_exterior[2]; size_t size_alloc[2]; size_t size_dc[2]; @@ -1155,73 +1155,73 @@ static int arrays_alloc(EGSContext *ctx, const size_t domain_size[2]) slice_interior[1].end = slice_interior[1].start + dc->interior.size[0]; slice_interior[1].step = 1; - ret = mg2di_ndarray_alloc(&priv->u_base, 2, size_alloc, NDARRAY_ALLOC_ZERO); + ret = ndarray_alloc(&priv->u_base, 2, size_alloc, NDARRAY_ALLOC_ZERO); if (ret < 0) return ret; - ret = mg2di_ndarray_alloc(&priv->u_next_base, 2, size_alloc, NDARRAY_ALLOC_ZERO); + ret = ndarray_alloc(&priv->u_next_base, 2, size_alloc, NDARRAY_ALLOC_ZERO); if (ret < 0) return ret; - ret = mg2di_ndarray_slice(&ctx->u, priv->u_base, slice_interior); + ret = ndarray_slice(&ctx->u, priv->u_base, slice_interior); if (ret < 0) return ret; - ret = mg2di_ndarray_slice(&priv->u_next, priv->u_next_base, slice_interior); + ret = ndarray_slice(&priv->u_next, priv->u_next_base, slice_interior); if (ret < 0) return ret; - ret = mg2di_ndarray_slice(&ctx->u_exterior, priv->u_base, slice_exterior); + ret = ndarray_slice(&ctx->u_exterior, priv->u_base, slice_exterior); if (ret < 0) return ret; - ret = mg2di_ndarray_slice(&priv->u_next_exterior, priv->u_next_base, slice_exterior); + ret = ndarray_slice(&priv->u_next_exterior, priv->u_next_base, slice_exterior); if (ret < 0) return ret; - ret = mg2di_ndarray_alloc(&priv->rhs_base, 2, size_alloc, 0); + ret = ndarray_alloc(&priv->rhs_base, 2, size_alloc, 0); if (ret < 0) return ret; - ret = mg2di_ndarray_slice(&ctx->rhs, priv->rhs_base, slice_interior); + ret = ndarray_slice(&ctx->rhs, priv->rhs_base, slice_interior); if (ret < 0) return ret; - ret = mg2di_ndarray_alloc(&priv->residual_base, 2, size_alloc, + ret = ndarray_alloc(&priv->residual_base, 2, size_alloc, NDARRAY_ALLOC_ZERO); if (ret < 0) return ret; - ret = mg2di_ndarray_slice(&ctx->residual, priv->residual_base, slice_interior); + ret = ndarray_slice(&ctx->residual, priv->residual_base, slice_interior); if (ret < 0) return ret; size_dc[0] = size_alloc[0] * MG2D_DIFF_COEFF_NB; size_dc[1] = size_alloc[1]; - ret = mg2di_ndarray_alloc(&priv->diff_coeffs_public_base, 2, + ret = ndarray_alloc(&priv->diff_coeffs_public_base, 2, size_dc, NDARRAY_ALLOC_ZERO); if (ret < 0) return ret; - ret = mg2di_ndarray_alloc(&priv->diff_coeffs_base, 2, + ret = ndarray_alloc(&priv->diff_coeffs_base, 2, size_dc, NDARRAY_ALLOC_ZERO); if (ret < 0) return ret; for (int i = 0; i < ARRAY_ELEMS(ctx->diff_coeffs); i++) { - const Slice slice_dc[2] = { SLICE(i * size_alloc[0], (i + 1) * size_alloc[0], 1), SLICE_NULL }; + const NDASlice slice_dc[2] = { NDASLICE(i * size_alloc[0], (i + 1) * size_alloc[0], 1), NDASLICE_NULL }; NDArray *tmp; - ret = mg2di_ndarray_slice(&tmp, priv->diff_coeffs_public_base, slice_dc); + ret = ndarray_slice(&tmp, priv->diff_coeffs_public_base, slice_dc); if (ret < 0) return ret; - ret = mg2di_ndarray_slice(&ctx->diff_coeffs[i], tmp, slice_interior); - mg2di_ndarray_free(&tmp); + ret = ndarray_slice(&ctx->diff_coeffs[i], tmp, slice_interior); + ndarray_free(&tmp); if (ret < 0) return ret; - ret = mg2di_ndarray_slice(&tmp, priv->diff_coeffs_base, slice_dc); + ret = ndarray_slice(&tmp, priv->diff_coeffs_base, slice_dc); if (ret < 0) return ret; - ret = mg2di_ndarray_slice(&priv->diff_coeffs[i], tmp, slice_interior); - mg2di_ndarray_free(&tmp); + ret = ndarray_slice(&priv->diff_coeffs[i], tmp, slice_interior); + ndarray_free(&tmp); if (ret < 0) return ret; } @@ -1384,28 +1384,28 @@ void mg2di_egs_free(EGSContext **pctx) exact_arrays_free(ctx); - mg2di_ndarray_free(&ctx->u); - mg2di_ndarray_free(&ctx->u_exterior); - mg2di_ndarray_free(&ctx->priv->u_base); - mg2di_ndarray_free(&ctx->priv->u_next); - mg2di_ndarray_free(&ctx->priv->u_next_exterior); - mg2di_ndarray_free(&ctx->priv->u_next_base); + ndarray_free(&ctx->u); + ndarray_free(&ctx->u_exterior); + ndarray_free(&ctx->priv->u_base); + ndarray_free(&ctx->priv->u_next); + ndarray_free(&ctx->priv->u_next_exterior); + ndarray_free(&ctx->priv->u_next_base); - mg2di_ndarray_free(&ctx->rhs); - mg2di_ndarray_free(&ctx->priv->rhs_base); + ndarray_free(&ctx->rhs); + ndarray_free(&ctx->priv->rhs_base); - mg2di_ndarray_free(&ctx->residual); - mg2di_ndarray_free(&ctx->priv->residual_base); + ndarray_free(&ctx->residual); + ndarray_free(&ctx->priv->residual_base); for (int i = 0; i < ARRAY_ELEMS(ctx->diff_coeffs); i++) { - mg2di_ndarray_free(&ctx->diff_coeffs[i]); + ndarray_free(&ctx->diff_coeffs[i]); } - mg2di_ndarray_free(&ctx->priv->diff_coeffs_public_base); + ndarray_free(&ctx->priv->diff_coeffs_public_base); for (int i = 0; i < ARRAY_ELEMS(ctx->priv->diff_coeffs); i++) { - mg2di_ndarray_free(&ctx->priv->diff_coeffs[i]); + ndarray_free(&ctx->priv->diff_coeffs[i]); } - mg2di_ndarray_free(&ctx->priv->diff_coeffs_base); + ndarray_free(&ctx->priv->diff_coeffs_base); for (int i = 0; i < ARRAY_ELEMS(ctx->boundaries); i++) mg2di_bc_free(&ctx->boundaries[i]); -- cgit v1.2.3