From 4633e6db474e7825979bef26e68628d492a63eb1 Mon Sep 17 00:00:00 2001 From: Anton Khirnov Date: Sun, 28 Aug 2016 10:53:54 +0200 Subject: eval: refactor evaluating the metric --- eval_metric_template.c | 68 ++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 68 insertions(+) create mode 100644 eval_metric_template.c (limited to 'eval_metric_template.c') diff --git a/eval_metric_template.c b/eval_metric_template.c new file mode 100644 index 0000000..8443785 --- /dev/null +++ b/eval_metric_template.c @@ -0,0 +1,68 @@ +/* + * Copyright 2016 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 . + */ + +/** + * @file + * the template for the loop over the grid points +*/ + +#define FUNC3(a, b, c) a ## _ ## b ## _ ## c +#define FUNC2(a, b, c) FUNC3(a, b, c) +#define FUNC(x) FUNC2(x, COMP, DIFF) + +#define DIFF_RHO (DIFF % 3) +#define DIFF_Z (DIFF / 3) + +static void FUNC(eval_metric)(const BDContext *bd, double *out, ptrdiff_t out_stride, + int nb_coords_rho, const double *rho, + int nb_coords_z, const double *z, + double *psi[9], double *q[9]) +{ + BDPriv *s = bd->priv; + + for (int i = 0; i < nb_coords_z; i++) { + const double zz = z[i]; + double *dst = out + i * out_stride; + + for (int j = 0; j < nb_coords_rho; j++) { + const int idx = i * nb_coords_rho + j; + const double rrho = rho[j]; + + const double ppsi = psi[0][idx]; + const double psi2 = SQR(ppsi); + const double psi3 = psi2 * ppsi; + const double psi4 = SQR(psi2); + + const double base = psi4 * exp(2 * q[0][idx]); + const double dqrho = DIFF_RHO > 0 ? q[1][idx] + 2 * psi[1][idx] / ppsi : 0; + const double dqz = DIFF_Z > 0 ? q[3][idx] + 2 * psi[3][idx] / ppsi : 0; + const double dq2rho = DIFF_RHO > 1 ? q[2][idx] + 2 * psi[2][idx] / ppsi - 2 * SQR(psi[1][idx] / ppsi) : 0; + const double dq2z = DIFF_Z > 1 ? q[6][idx] + 2 * psi[6][idx] / ppsi - 2 * SQR(psi[3][idx] / ppsi) : 0; + const double dq2rho_z = DIFF_RHO && DIFF_Z ? q[4][idx] + 2 * psi[4][idx] / ppsi - 2 * psi[1][idx] * psi[3][idx] / psi2 : 0; + + double val; + + DO_EVAL; + + *dst++ = val; + } + } +} + +#undef FUNC3 +#undef FUNC2 +#undef FUNC -- cgit v1.2.3