diff options
Diffstat (limited to 'eval.c')
| -rw-r--r-- | eval.c | 32 |
1 files changed, 16 insertions, 16 deletions
@@ -35,7 +35,7 @@ int bd_eval_psi(const BDContext *bd, double *psi, unsigned int psi_stride) { BDPriv *s = bd->priv; - const double sfr = bd->basis_scale_factor_rho, sfz = bd->basis_scale_factor_z; + const double *sf = bd->basis_scale_factor; double (*eval)(double coord, int idx, double sf); double *basis_val_rho = NULL, *basis_val_z = NULL; @@ -53,34 +53,34 @@ int bd_eval_psi(const BDContext *bd, add = 1.0; /* precompute the basis values on the grid points */ - basis_val_rho = malloc(sizeof(*basis_val_rho) * bd->nb_coeffs_rho * nb_coords_rho); - basis_val_z = malloc(sizeof(*basis_val_z) * bd->nb_coeffs_z * nb_coords_z); + basis_val_rho = malloc(sizeof(*basis_val_rho) * s->nb_coeffs[0] * nb_coords_rho); + basis_val_z = malloc(sizeof(*basis_val_z) * s->nb_coeffs[1] * nb_coords_z); if (!basis_val_rho || !basis_val_z) { ret = -ENOMEM; goto fail; } switch (diff_order[0]) { - case 0: eval = s->basis->eval; break; - case 1: eval = s->basis->eval_diff1; break; - case 2: eval = s->basis->eval_diff2; break; + case 0: eval = s->basis[0]->eval; break; + case 1: eval = s->basis[0]->eval_diff1; break; + case 2: eval = s->basis[0]->eval_diff2; break; } for (int i = 0; i < nb_coords_rho; i++) { double rrho = rho[i]; - for (int j = 0; j < bd->nb_coeffs_rho; j++) - basis_val_rho[i * bd->nb_coeffs_rho + j] = eval(rrho, j, sfr); + for (int j = 0; j < s->nb_coeffs[0]; j++) + basis_val_rho[i * s->nb_coeffs[0] + j] = eval(rrho, j, sf[0]); } switch (diff_order[1]) { - case 0: eval = s->basis->eval; break; - case 1: eval = s->basis->eval_diff1; break; - case 2: eval = s->basis->eval_diff2; break; + case 0: eval = s->basis[1]->eval; break; + case 1: eval = s->basis[1]->eval_diff1; break; + case 2: eval = s->basis[1]->eval_diff2; break; } for (int i = 0; i < nb_coords_z; i++) { double zz = z[i]; - for (int j = 0; j < bd->nb_coeffs_z; j++) - basis_val_z[i * bd->nb_coeffs_z + j] = eval(zz, j, sfz); + for (int j = 0; j < s->nb_coeffs[1]; j++) + basis_val_z[i * s->nb_coeffs[1] + j] = eval(zz, j, sf[1]); } for (int i = 0; i < nb_coords_z; i++) { @@ -89,9 +89,9 @@ int bd_eval_psi(const BDContext *bd, for (int j = 0; j < nb_coords_rho; j++) { double ppsi = add; - for (int m = 0; m < bd->nb_coeffs_z; m++) - for (int n = 0; n < bd->nb_coeffs_rho; n++) - ppsi += s->coeffs[m * bd->nb_coeffs_rho + n] * basis_val_rho[j * bd->nb_coeffs_rho + n] * basis_val_z[i * bd->nb_coeffs_z + m]; + for (int m = 0; m < s->nb_coeffs[1]; m++) + for (int n = 0; n < s->nb_coeffs[0]; n++) + ppsi += s->coeffs[m * s->nb_coeffs[0] + n] * basis_val_rho[j * s->nb_coeffs[0] + n] * basis_val_z[i * s->nb_coeffs[1] + m]; *dst++ = ppsi; } |