Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
12 changes: 0 additions & 12 deletions examples/performance/bracket/bracket.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -68,9 +68,6 @@ int main(int argc, char** argv) {
ITERATOR_TEST_BLOCK("Bracket [2D,3D] ARAKAWA",
result = bracket(a, c, BRACKET_ARAKAWA););

ITERATOR_TEST_BLOCK("Bracket [2D,3D] ARAKAWA_OLD",
result = bracket(a, c, BRACKET_ARAKAWA_OLD););

ITERATOR_TEST_BLOCK("Bracket [2D,3D] SIMPLE",
result = bracket(a, c, BRACKET_SIMPLE););

Expand All @@ -81,21 +78,12 @@ int main(int argc, char** argv) {
ITERATOR_TEST_BLOCK("Bracket [3D,3D] ARAKAWA",
result = bracket(a, b, BRACKET_ARAKAWA););

ITERATOR_TEST_BLOCK("Bracket [3D,3D] ARAKAWA_OLD",
result = bracket(a, b, BRACKET_ARAKAWA_OLD););

ITERATOR_TEST_BLOCK("Bracket [3D,3D] SIMPLE",
result = bracket(a, b, BRACKET_SIMPLE););

ITERATOR_TEST_BLOCK("Bracket [3D,3D] DEFAULT", result = bracket(a, b, BRACKET_STD););
}

// Uncomment below for a "correctness" check
// Field3D resNew = bracket(a, b, BRACKET_ARAKAWA); mesh->communicate(resNew);
// Field3D resOld = bracket(a, b, BRACKET_ARAKAWA_OLD); mesh->communicate(resOld);
// time_output << "Max abs diff is
// "<<max(abs(resNew-resOld),true)/max(abs(resOld),true)<<std::endl;

Comment thread
ZedThree marked this conversation as resolved.
if (profileMode) {
int nthreads = 0;
#if BOUT_USE_OPENMP
Expand Down
12 changes: 5 additions & 7 deletions include/bout/difops.hxx
Original file line number Diff line number Diff line change
Expand Up @@ -271,18 +271,16 @@ Field3D b0xGrad_dot_Grad(const Field3D& phi, const Field3D& A,
* Poisson bracket methods
*/
enum class BRACKET_METHOD {
standard, ///< Use b0xGrad_dot_Grad
simple, ///< Keep only terms in X-Z
arakawa, ///< Arakawa method in X-Z (optimised)
ctu, ///< Corner Transport Upwind (CTU) method. Explicit method only, needs the
/// timestep from the solver
arakawa_old ///< Older version, for regression testing of optimised version.
standard, ///< Use b0xGrad_dot_Grad
simple, ///< Keep only terms in X-Z
arakawa, ///< Arakawa method in X-Z
ctu, ///< Corner Transport Upwind (CTU) method. Explicit method only, needs the
/// timestep from the solver
};
constexpr BRACKET_METHOD BRACKET_STD = BRACKET_METHOD::standard;
constexpr BRACKET_METHOD BRACKET_SIMPLE = BRACKET_METHOD::simple;
constexpr BRACKET_METHOD BRACKET_ARAKAWA = BRACKET_METHOD::arakawa;
constexpr BRACKET_METHOD BRACKET_CTU = BRACKET_METHOD::ctu;
constexpr BRACKET_METHOD BRACKET_ARAKAWA_OLD = BRACKET_METHOD::arakawa_old;

/*!
* Compute advection operator terms, which can be cast as
Expand Down
2 changes: 1 addition & 1 deletion manual/sphinx/user_docs/python_boutpp.rst
Original file line number Diff line number Diff line change
Expand Up @@ -160,7 +160,7 @@ A real example - check derivative contributions:
phi[-1,:,z]=phi_arr
with open(path+"/equilibrium/phi_eq.dat","rb") as inf:
phi_arr=np.fromfile(inf,dtype=np.double)
bm="BRACKET_ARAKAWA_OLD"
bm="BRACKET_ARAKAWA"

for tind in range(len(times)):
vort = Field3D.fromCollect("vort" ,path=path,tind=tind,info=False)
Expand Down
97 changes: 0 additions & 97 deletions src/mesh/difops.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -771,46 +771,6 @@ Field3D bracket(const Field3D& f, const Field2D& g, BRACKET_METHOD method,

break;
}
case BRACKET_ARAKAWA_OLD: {
#if not(BOUT_USE_METRIC_3D)
const int ncz = mesh->LocalNz;
BOUT_OMP_PERF(parallel for)
for (int jx = mesh->xstart; jx <= mesh->xend; jx++) {
for (int jy = mesh->ystart; jy <= mesh->yend; jy++) {
const BoutReal partialFactor = 1.0 / (12 * metric->dz(jx, jy));
const BoutReal spacingFactor = partialFactor / metric->dx(jx, jy);
for (int jz = 0; jz < mesh->LocalNz; jz++) {
const int jzp = jz + 1 < ncz ? jz + 1 : 0;
// Above is alternative to const int jzp = (jz + 1) % ncz;
const int jzm = jz - 1 >= 0 ? jz - 1 : ncz - 1;
// Above is alternative to const int jzmTmp = (jz - 1 + ncz) % ncz;

// J++ = DDZ(f)*DDX(g) - DDX(f)*DDZ(g)
BoutReal Jpp =
((f(jx, jy, jzp) - f(jx, jy, jzm)) * (g(jx + 1, jy) - g(jx - 1, jy))
- (f(jx + 1, jy, jz) - f(jx - 1, jy, jz)) * (g(jx, jy) - g(jx, jy)));

// J+x
BoutReal Jpx = (g(jx + 1, jy) * (f(jx + 1, jy, jzp) - f(jx + 1, jy, jzm))
- g(jx - 1, jy) * (f(jx - 1, jy, jzp) - f(jx - 1, jy, jzm))
- g(jx, jy) * (f(jx + 1, jy, jzp) - f(jx - 1, jy, jzp))
+ g(jx, jy) * (f(jx + 1, jy, jzm) - f(jx - 1, jy, jzm)));

// Jx+
BoutReal Jxp = (g(jx + 1, jy) * (f(jx, jy, jzp) - f(jx + 1, jy, jz))
- g(jx - 1, jy) * (f(jx - 1, jy, jz) - f(jx, jy, jzm))
- g(jx - 1, jy) * (f(jx, jy, jzp) - f(jx - 1, jy, jz))
+ g(jx + 1, jy) * (f(jx + 1, jy, jz) - f(jx, jy, jzm)));

result(jx, jy, jz) = (Jpp + Jpx + Jxp) * spacingFactor;
}
}
}
#else
throw BoutException("BRACKET_ARAKAWA_OLD not valid with 3D metrics yet.");
#endif
break;
}
case BRACKET_SIMPLE: {
// Use a subset of terms for comparison to BOUT-06
result = VDDX(DDZ(f, outloc), g, outloc);
Expand Down Expand Up @@ -1090,63 +1050,6 @@ Field3D bracket(const Field3D& f, const Field3D& g, BRACKET_METHOD method,
}
break;
}
case BRACKET_ARAKAWA_OLD: {
// Arakawa scheme for perpendicular flow

const int ncz = mesh->LocalNz;

// We need to discard const qualifier in order to manipulate
// storage array directly
Field3D f_temp = f;
Field3D g_temp = g;

BOUT_OMP_PERF(parallel for)
for (int jx = mesh->xstart; jx <= mesh->xend; jx++) {
for (int jy = mesh->ystart; jy <= mesh->yend; jy++) {
#if not(BOUT_USE_METRIC_3D)
const BoutReal spacingFactor =
1.0 / (12 * metric->dz(jx, jy) * metric->dx(jx, jy));
#endif
const BoutReal* Fxm = f_temp(jx - 1, jy);
const BoutReal* Fx = f_temp(jx, jy);
const BoutReal* Fxp = f_temp(jx + 1, jy);
const BoutReal* Gxm = g_temp(jx - 1, jy);
const BoutReal* Gx = g_temp(jx, jy);
const BoutReal* Gxp = g_temp(jx + 1, jy);
for (int jz = 0; jz < mesh->LocalNz; jz++) {
#if BOUT_USE_METRIC_3D
const BoutReal spacingFactor =
1.0 / (12 * metric->dz(jx, jy, jz) * metric->dx(jx, jy, jz));
#endif
const int jzp = jz + 1 < ncz ? jz + 1 : 0;
// Above is alternative to const int jzp = (jz + 1) % ncz;
const int jzm = jz - 1 >= 0 ? jz - 1 : ncz - 1;
// Above is alternative to const int jzm = (jz - 1 + ncz) % ncz;

// J++ = DDZ(f)*DDX(g) - DDX(f)*DDZ(g)
// NOLINTNEXTLINE
BoutReal Jpp = ((Fx[jzp] - Fx[jzm]) * (Gxp[jz] - Gxm[jz])
- (Fxp[jz] - Fxm[jz]) * (Gx[jzp] - Gx[jzm]));

// J+x
// NOLINTNEXTLINE
BoutReal Jpx =
(Gxp[jz] * (Fxp[jzp] - Fxp[jzm]) - Gxm[jz] * (Fxm[jzp] - Fxm[jzm])
- Gx[jzp] * (Fxp[jzp] - Fxm[jzp]) + Gx[jzm] * (Fxp[jzm] - Fxm[jzm]));

// Jx+
// NOLINTNEXTLINE
BoutReal Jxp =
(Gxp[jzp] * (Fx[jzp] - Fxp[jz]) - Gxm[jzm] * (Fxm[jz] - Fx[jzm])
- Gxm[jzp] * (Fx[jzp] - Fxm[jz]) + Gxp[jzm] * (Fxp[jz] - Fx[jzm]));

result(jx, jy, jz) = (Jpp + Jpx + Jxp) * spacingFactor;
}
}
}

break;
}
case BRACKET_SIMPLE: {
// Use a subset of terms for comparison to BOUT-06
result = VDDX(DDZ(f, outloc), g, outloc) + VDDZ(-DDX(f, outloc), g, outloc);
Expand Down
1 change: 0 additions & 1 deletion tests/MMS/spatial/advection/runtest
Original file line number Diff line number Diff line change
Expand Up @@ -43,7 +43,6 @@ build_and_log("MMS steady-state advection test")
# List of options to be passed for each test
options = [
("method=2", "Arakawa", "-^", 2),
("method=4", "Arakawa-old", "-.", 2),
(
"method=1 mesh:ddx:upwind=u1 mesh:ddz:upwind=u1",
"SIMPLE: 1st order upwind",
Expand Down
2 changes: 1 addition & 1 deletion tools/pylib/_boutpp_build/other_enums.hxx
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
#ifdef YOU_SHOULDNT_READ_THIS
#error YOU_SHOULDNT_BE_HERE
enum class BRACKET_METHOD { standard, simple, arakawa, ctu, arakawa_old };
enum class BRACKET_METHOD { standard, simple, arakawa, ctu };
#endif
2 changes: 0 additions & 2 deletions tools/pylib/_boutpp_build/scan_enums.py
Original file line number Diff line number Diff line change
Expand Up @@ -78,8 +78,6 @@ def __str__(self):
"ARAKAWA": "arakawa",
"BRACKET_CTU": "ctu",
"CTU": "ctu",
"BRACKET_ARAKAWA_OLD": "arakawa_old",
"ARAKAWA_OLD": "arakawa_old",
}
enums["CELL_LOC"].extra = {
"CELL_DEFAULT": "deflt",
Expand Down
Loading