From 0007a0b0c11fa7c12b228883453368f105a4324b Mon Sep 17 00:00:00 2001 From: Anton Khirnov Date: Thu, 16 Nov 2017 13:11:07 +0100 Subject: Initial commit. The following code is present: * the basis API * the BiCGSTAB solver * the pseudospectral linear system solver * helper APIs: - threadpool - logging - cpuid --- tests/pssolve.c | 139 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 139 insertions(+) create mode 100644 tests/pssolve.c (limited to 'tests') diff --git a/tests/pssolve.c b/tests/pssolve.c new file mode 100644 index 0000000..94afadd --- /dev/null +++ b/tests/pssolve.c @@ -0,0 +1,139 @@ +/* + * Copyright 2017 Anton Khirnov + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ + +#include +#include +#include + +#include "basis.h" +#include "common.h" +#include "log.h" +#include "pssolve.h" + +#define N_X 8 +#define N_Z 4 +#define FUNC_X 4 +#define FUNC_Z 2 +#define TOL 5e-15 +static int scalar0(void) +{ + BasisSetContext *basis = NULL; + PSSolveContext *ctx = NULL; + + double eq_coeffs00[PSSOLVE_DIFF_ORDER_NB][N_X * N_Z]; + const double *eq_coeffs0[PSSOLVE_DIFF_ORDER_NB]; + const double *(*eq_coeffs)[PSSOLVE_DIFF_ORDER_NB]; + double rhs[N_X * N_Z]; + double coeffs[N_X * N_Z]; + int ret; + + ret = tdi_basis_init(&basis, BASIS_FAMILY_SB_EVEN, 1.0); + if (ret < 0) + goto finish; + + ret = tdi_pssolve_context_alloc(&ctx, 1); + if (ret < 0) + return ret; + + ctx->basis[0][0] = basis; + ctx->basis[0][1] = basis; + + ctx->solve_order[0][0] = N_X; + ctx->solve_order[0][1] = N_Z; + + ctx->logger.log = tdi_log_default_callback; + + ret = tdi_pssolve_context_init(ctx); + if (ret < 0) + goto finish; + + for (int j = 0; j < N_Z; j++) { + const double z = ctx->colloc_grid[0][1][j]; + + for (int i = 0; i < N_X; i++) { + const double x = ctx->colloc_grid[0][0][i]; + + const int idx = j * N_X + i; + + eq_coeffs00[PSSOLVE_DIFF_ORDER_20][idx] = 1.0; + eq_coeffs00[PSSOLVE_DIFF_ORDER_02][idx] = 1.0; + eq_coeffs00[PSSOLVE_DIFF_ORDER_10][idx] = 0.0; + eq_coeffs00[PSSOLVE_DIFF_ORDER_01][idx] = 0.0; + eq_coeffs00[PSSOLVE_DIFF_ORDER_11][idx] = 0.0; + eq_coeffs00[PSSOLVE_DIFF_ORDER_00][idx] = 0.0; + + rhs[idx] = tdi_basis_eval(basis, BS_EVAL_TYPE_DIFF2, x, FUNC_X) * + tdi_basis_eval(basis, BS_EVAL_TYPE_VALUE, z, FUNC_Z) + + tdi_basis_eval(basis, BS_EVAL_TYPE_VALUE, x, FUNC_X) * + tdi_basis_eval(basis, BS_EVAL_TYPE_DIFF2, z, FUNC_Z); + } + } + for (int i = 0; i < PSSOLVE_DIFF_ORDER_NB; i++) + eq_coeffs0[i] = eq_coeffs00[i]; + eq_coeffs = &eq_coeffs0; + + ret = tdi_pssolve_solve(ctx, &eq_coeffs, rhs, coeffs); + if (ret < 0) + goto finish; + + for (int j = 0; j < N_Z; j++) + for (int i = 0; i < N_X; i++) { + const int idx = j * N_X + i; + const double val = (i == FUNC_X && j == FUNC_Z) ? 1.0 : 0.0; + + if (fabs(coeffs[idx] - val) > TOL) { + fprintf(stderr, "unexpected value %g at %d/%d\n", + coeffs[idx], i, j); + ret = -1; + goto finish; + } + } + +finish: + tdi_pssolve_context_free(&ctx); + tdi_basis_free(&basis); + + return ret; +} +#undef N +#undef FUNC_X +#undef FUNC_Z +#undef TOL + +static const struct { + const char *name; + int (*test)(void); +} tests[] = { + { "scalar0", scalar0 }, +}; + +int main(void) +{ + int ret = 0; + + for (int i = 0; i < ARRAY_ELEMS(tests); i++) { + fprintf(stderr, "executing test '%s'\n", tests[i].name); + ret = tests[i].test(); + if (ret < 0) { + fprintf(stderr, "test '%s' failed\n", tests[i].name); + return -1; + } + fprintf(stderr, "test '%s' succeeded\n", tests[i].name); + } + + return 0; +} -- cgit v1.2.3