diff --git a/include/bout/interpolation_xz.hxx b/include/bout/interpolation_xz.hxx index def8a60a3e..6c7419f7e4 100644 --- a/include/bout/interpolation_xz.hxx +++ b/include/bout/interpolation_xz.hxx @@ -238,6 +238,9 @@ public: const std::string& region = "RGN_NOBNDRY") const override; }; +/// XZLagrange4pt interpolation class +/// +/// Does not support MPI splitting in X class XZLagrange4pt : public XZInterpolation { Tensor i_corner; // x-index of bottom-left grid point Tensor k_corner; // z-index of bottom-left grid point @@ -271,6 +274,9 @@ public: BoutReal lagrange_4pt(const BoutReal v[], BoutReal offset) const; }; +/// XZBilinear interpolation calss +/// +/// Does not support MPI splitting in X. class XZBilinear : public XZInterpolation { Tensor i_corner; // x-index of bottom-left grid point Tensor k_corner; // z-index of bottom-left grid point diff --git a/src/mesh/interpolation/bilinear_xz.cxx b/src/mesh/interpolation/bilinear_xz.cxx index 8445764a8f..1ed18a7303 100644 --- a/src/mesh/interpolation/bilinear_xz.cxx +++ b/src/mesh/interpolation/bilinear_xz.cxx @@ -31,6 +31,10 @@ XZBilinear::XZBilinear(int y_offset, Mesh* mesh) : XZInterpolation(y_offset, mesh), w0(localmesh), w1(localmesh), w2(localmesh), w3(localmesh) { + if (localmesh->getNXPE() > 1) { + throw BoutException("XZBilinear interpolation does not support MPI splitting in X"); + } + // Index arrays contain guard cells in order to get subscripts right i_corner.reallocate(localmesh->LocalNx, localmesh->LocalNy, localmesh->LocalNz); k_corner.reallocate(localmesh->LocalNx, localmesh->LocalNy, localmesh->LocalNz); diff --git a/src/mesh/interpolation/hermite_spline_xz.cxx b/src/mesh/interpolation/hermite_spline_xz.cxx index c0040d096e..06d5b393a2 100644 --- a/src/mesh/interpolation/hermite_spline_xz.cxx +++ b/src/mesh/interpolation/hermite_spline_xz.cxx @@ -101,8 +101,8 @@ class IndConverter { } }; -XZHermiteSpline::XZHermiteSpline(int y_offset, Mesh* mesh) - : XZInterpolation(y_offset, mesh), h00_x(localmesh), h01_x(localmesh), +XZHermiteSpline::XZHermiteSpline(int y_offset, Mesh* meshin) + : XZInterpolation(y_offset, meshin), h00_x(localmesh), h01_x(localmesh), h10_x(localmesh), h11_x(localmesh), h00_z(localmesh), h01_z(localmesh), h10_z(localmesh), h11_z(localmesh) { @@ -346,6 +346,9 @@ Field3D XZHermiteSpline::interpolate(const Field3D& f, const std::string& region ASSERT1(f.getMesh() == localmesh); Field3D f_interp{emptyFrom(f)}; + const auto region2 = + y_offset == 0 ? "RGN_NOY" : fmt::format("RGN_YPAR_{:+d}", y_offset); + #if USE_NEW_WEIGHTS #ifdef HS_USE_PETSC BoutReal* ptr; @@ -355,7 +358,6 @@ Field3D XZHermiteSpline::interpolate(const Field3D& f, const std::string& region VecRestoreArray(rhs, &ptr); MatMult(petscWeights, rhs, result); VecGetArrayRead(result, &cptr); - const auto region2 = y_offset == 0 ? region : fmt::format("RGN_YPAR_{:+d}", y_offset); BOUT_FOR(i, f.getRegion(region2)) { f_interp[i] = cptr[int(i)]; ASSERT2(std::isfinite(cptr[int(i)])); @@ -375,11 +377,10 @@ Field3D XZHermiteSpline::interpolate(const Field3D& f, const std::string& region } } #endif - return f_interp; #else // Derivatives are used for tension and need to be on dimensionless // coordinates - const auto region2 = fmt::format("RGN_YPAR_{:+d}", y_offset); + // f has been communcated, and thus we can assume that the x-boundaries are // also valid in the y-boundary. Thus the differentiated field needs no // extra comms. @@ -418,8 +419,10 @@ Field3D XZHermiteSpline::interpolate(const Field3D& f, const std::string& region ASSERT2(std::isfinite(f_interp[iyp]) || i.x() < localmesh->xstart || i.x() > localmesh->xend); } - return f_interp; #endif + f_interp.setRegion(region2); + ASSERT2(f_interp.getRegionID()); + return f_interp; } Field3D XZHermiteSpline::interpolate(const Field3D& f, const Field3D& delta_x, diff --git a/src/mesh/interpolation/lagrange_4pt_xz.cxx b/src/mesh/interpolation/lagrange_4pt_xz.cxx index 92c14ecfd5..e16b2699d1 100644 --- a/src/mesh/interpolation/lagrange_4pt_xz.cxx +++ b/src/mesh/interpolation/lagrange_4pt_xz.cxx @@ -29,6 +29,11 @@ XZLagrange4pt::XZLagrange4pt(int y_offset, Mesh* mesh) : XZInterpolation(y_offset, mesh), t_x(localmesh), t_z(localmesh) { + if (localmesh->getNXPE() > 1) { + throw BoutException( + "XZLagrange4pt interpolation does not support MPI splitting in X"); + } + // Index arrays contain guard cells in order to get subscripts right i_corner.reallocate(localmesh->LocalNx, localmesh->LocalNy, localmesh->LocalNz); k_corner.reallocate(localmesh->LocalNx, localmesh->LocalNy, localmesh->LocalNz);