diff options
| author | Anton Khirnov <anton@khirnov.net> | 2019-06-10 17:41:18 +0200 |
|---|---|---|
| committer | Anton Khirnov <anton@khirnov.net> | 2019-06-10 18:04:28 +0200 |
| commit | 2a22da0b5724e68e4b2644cb2289f2947f9d5d46 (patch) | |
| tree | 74015736bc164827f8828ed81a48679edc188549 | |
| parent | 1afea75597177b736429f8f6368bb5ecd0f0a389 (diff) | |
Update for the new mg2d diff coeffs API.
Treat the PDE coeff non-regularities properly.
| -rw-r--r-- | src/maximal_slicing_axi_mg.c | 22 |
1 files changed, 15 insertions, 7 deletions
diff --git a/src/maximal_slicing_axi_mg.c b/src/maximal_slicing_axi_mg.c index cf56fb9..365c3ce 100644 --- a/src/maximal_slicing_axi_mg.c +++ b/src/maximal_slicing_axi_mg.c @@ -697,10 +697,13 @@ static void fill_eq_coeffs(MSMGContext *ctx, CoordPatch *cp, MG2DContext *solver double *a_Xtx = CCTK_VarDataPtr(gh, 0, "ML_BSSN::Xt1"); double *a_Xtz = CCTK_VarDataPtr(gh, 0, "ML_BSSN::Xt3"); + solver->diff_coeffs[MG2D_DIFF_COEFF_20]->boundaries[MG2D_BOUNDARY_0L].flags = MG2D_DC_FLAG_DISCONT; + solver->diff_coeffs[MG2D_DIFF_COEFF_10]->boundaries[MG2D_BOUNDARY_0L].flags = MG2D_DC_FLAG_POLE; + for (int idx_z = 0; idx_z < solver->local_size[1]; idx_z++) for (int idx_x = 0; idx_x < solver->local_size[0]; idx_x++) { const int idx_src = CCTK_GFINDEX3D(gh, idx_x + cp->offset_left[0], cp->y_idx, idx_z + cp->offset_left[1]); - const int idx_dc = idx_z * solver->diff_coeffs_stride + idx_x; + const int idx_dc = idx_z * solver->diff_coeffs[0]->stride + idx_x; const int idx_rhs = idx_z * solver->rhs_stride + idx_x; const double x = a_x[idx_src]; @@ -770,13 +773,18 @@ static void fill_eq_coeffs(MSMGContext *ctx, CoordPatch *cp, MG2DContext *solver Xx = SQR(phi) * (Xtx + (phi_dx * gtu[0][0] + phi_dz * gtu[0][2]) / phi); Xz = SQR(phi) * (Xtz + (phi_dx * gtu[0][2] + phi_dz * gtu[2][2]) / phi); - solver->diff_coeffs[MG2D_DIFF_COEFF_20][idx_dc] = SQR(phi) * (gtu[0][0] + (on_axis ? gtu[1][1] : 0.0)); - solver->diff_coeffs[MG2D_DIFF_COEFF_02][idx_dc] = SQR(phi) * gtu[2][2]; - solver->diff_coeffs[MG2D_DIFF_COEFF_11][idx_dc] = SQR(phi) * gtu[0][2] * 2; - solver->diff_coeffs[MG2D_DIFF_COEFF_10][idx_dc] = -Xx + (on_axis ? 0.0 : SQR(phi) * gtu[1][1] / x); - solver->diff_coeffs[MG2D_DIFF_COEFF_01][idx_dc] = -Xz; - solver->diff_coeffs[MG2D_DIFF_COEFF_00][idx_dc] = -k2; + solver->diff_coeffs[MG2D_DIFF_COEFF_20]->data[idx_dc] = SQR(phi) * (gtu[0][0] + (on_axis ? gtu[1][1] : 0.0)); + solver->diff_coeffs[MG2D_DIFF_COEFF_02]->data[idx_dc] = SQR(phi) * gtu[2][2]; + solver->diff_coeffs[MG2D_DIFF_COEFF_11]->data[idx_dc] = SQR(phi) * gtu[0][2] * 2; + solver->diff_coeffs[MG2D_DIFF_COEFF_10]->data[idx_dc] = -Xx + (on_axis ? 0.0 : SQR(phi) * gtu[1][1] / x); + solver->diff_coeffs[MG2D_DIFF_COEFF_01]->data[idx_dc] = -Xz; + solver->diff_coeffs[MG2D_DIFF_COEFF_00]->data[idx_dc] = -k2; solver->rhs[idx_rhs] = k2; + + if (on_axis) { + solver->diff_coeffs[MG2D_DIFF_COEFF_20]->boundaries[MG2D_BOUNDARY_0L].val[idx_z] = SQR(phi) * gtu[1][1]; + solver->diff_coeffs[MG2D_DIFF_COEFF_10]->boundaries[MG2D_BOUNDARY_0L].val[idx_z] = SQR(phi) * gtu[1][1]; + } } } |