From 19df41ffde95df208c198ce3bf3041460280108e Mon Sep 17 00:00:00 2001 From: Anton Khirnov Date: Thu, 24 Mar 2016 10:53:52 +0100 Subject: Cartoonize, convert to C. --- src/apply_dissipation.c | 40 ++++++++++++++++++++++++++++++++++++++++ src/dissipation.c | 17 +++++------------ src/make.code.defn | 2 +- 3 files changed, 46 insertions(+), 13 deletions(-) create mode 100644 src/apply_dissipation.c diff --git a/src/apply_dissipation.c b/src/apply_dissipation.c new file mode 100644 index 0000000..6953757 --- /dev/null +++ b/src/apply_dissipation.c @@ -0,0 +1,40 @@ +#include + +#include "cctk.h" + +void apply_dissipation(const cGH *gh, const double *var, double *rhs, double dx[3], + int order, double *epsdis) +{ + double inv_dx[3] = { 1.0 / dx[0], 1.0 / dx[1], 1.0 / dx[2] }; + + int ni = gh->cctk_lsh[0]; + int nj = gh->cctk_lsh[1]; + int nk = gh->cctk_lsh[2]; + + double *y; + int y_idx = -1; + + y = CCTK_VarDataPtr(gh, 0, "grid::y"); + for (int j = 0; j < nj; j++) { + if (fabs(y[CCTK_GFINDEX3D(gh, 0, j, 0)]) < 1e-8) { + y_idx = j; + break; + } + } + + switch (order) { + case 5: +#pragma omp parallel for + for (int k = 3; k < nk - 3; k++) + for (int i = 3; i < nk - 3; i++) { + int idx = CCTK_GFINDEX3D(gh, i, y_idx, k); + rhs[idx] += epsdis[idx] / 64 * + ((var[idx - 3] - 6 * var[idx - 2] + 15 * var[idx - 1] - 20 * var[idx] + 15 * var[idx + 1] - 6 * var[idx + 2] + var[idx + 3]) * inv_dx[0] + + (var[idx - 3 * ni * nj] - 6 * var[idx - 2 * ni * nj] + 15 * var[idx - 1 * ni * nj] - 20 * var[idx] + 15 * var[idx + 1 * ni * nj] - 6 * var[idx + 2 * ni * nj] + var[idx + 3 * ni * nj]) * inv_dx[2]); +// + (var(i,j-3,k) - 6*var(i,j-2,k) + 15*var(i,j-1,k) - 20*var(i,j,k) + 15*var(i,j+1,k) - 6*var(i,j+2,k) + var(i,j+3,k) ) * idx(2) + } + break; + default: + CCTK_WARN(0, "internal error"); + } +} diff --git a/src/dissipation.c b/src/dissipation.c index a0ff84f..f4f3a15 100644 --- a/src/dissipation.c +++ b/src/dissipation.c @@ -8,15 +8,8 @@ #include "cctk_Arguments.h" #include "cctk_Parameters.h" -void CCTK_FCALL -CCTK_FNAME(apply_dissipation) (CCTK_REAL const * const var, - CCTK_REAL * const rhs, - int const * const ni, - int const * const nj, - int const * const nk, - CCTK_REAL const * const dx, - CCTK_INT const * const order, - CCTK_REAL const * const epsdis); +void apply_dissipation(const cGH *gh, const double *var, double *rhs, double dx[3], + int order, double *epsdis); static void apply (int const varindex, char const * const optstring, void * const arg); @@ -48,6 +41,8 @@ apply (int const varindex, char const * const optstring, void * const arg) assert (varindex >= 0); for (d=0; d