Fix eighth-order C derivative and lopsided stencils

This commit is contained in:
2026-05-14 20:32:53 +08:00
parent 207efdbbb2
commit 925889fbe9
5 changed files with 140 additions and 222 deletions

View File

@@ -398,7 +398,7 @@ void fderivs(const int ex[3],
#elif (ghost_width == 5)
/* ---- 8th-order ----------------------------------------------------- */
{
const int ord = 4;
const int ord = 5;
int iminF = 1, jminF = 1, kminF = 1;
if (Symmetry > NO_SYMM && fabs(Z[0]) < dZ) kminF = -3;
@@ -484,14 +484,14 @@ void fderivs(const int ex[3],
const size_t p = idx_ex(i0, j0, k0, ex);
fx[p] = d2dx * (
-fh[idx_fh_F_ord4(iF - 1, jF, kF, ex)] +
fh[idx_fh_F_ord4(iF + 1, jF, kF, ex)]);
-fh[idx_fh_F_ord5(iF - 1, jF, kF, ex)] +
fh[idx_fh_F_ord5(iF + 1, jF, kF, ex)]);
fy[p] = d2dy * (
-fh[idx_fh_F_ord4(iF, jF - 1, kF, ex)] +
fh[idx_fh_F_ord4(iF, jF + 1, kF, ex)]);
-fh[idx_fh_F_ord5(iF, jF - 1, kF, ex)] +
fh[idx_fh_F_ord5(iF, jF + 1, kF, ex)]);
fz[p] = d2dz * (
-fh[idx_fh_F_ord4(iF, jF, kF - 1, ex)] +
fh[idx_fh_F_ord4(iF, jF, kF + 1, ex)]);
-fh[idx_fh_F_ord5(iF, jF, kF - 1, ex)] +
fh[idx_fh_F_ord5(iF, jF, kF + 1, ex)]);
}
}
}
@@ -507,22 +507,22 @@ void fderivs(const int ex[3],
const size_t p = idx_ex(i0, j0, k0, ex);
fx[p] = d12dx * (
fh[idx_fh_F_ord4(iF - 2, jF, kF, ex)] -
EIT * fh[idx_fh_F_ord4(iF - 1, jF, kF, ex)] +
EIT * fh[idx_fh_F_ord4(iF + 1, jF, kF, ex)] -
fh[idx_fh_F_ord4(iF + 2, jF, kF, ex)]);
fh[idx_fh_F_ord5(iF - 2, jF, kF, ex)] -
EIT * fh[idx_fh_F_ord5(iF - 1, jF, kF, ex)] +
EIT * fh[idx_fh_F_ord5(iF + 1, jF, kF, ex)] -
fh[idx_fh_F_ord5(iF + 2, jF, kF, ex)]);
fy[p] = d12dy * (
fh[idx_fh_F_ord4(iF, jF - 2, kF, ex)] -
EIT * fh[idx_fh_F_ord4(iF, jF - 1, kF, ex)] +
EIT * fh[idx_fh_F_ord4(iF, jF + 1, kF, ex)] -
fh[idx_fh_F_ord4(iF, jF + 2, kF, ex)]);
fh[idx_fh_F_ord5(iF, jF - 2, kF, ex)] -
EIT * fh[idx_fh_F_ord5(iF, jF - 1, kF, ex)] +
EIT * fh[idx_fh_F_ord5(iF, jF + 1, kF, ex)] -
fh[idx_fh_F_ord5(iF, jF + 2, kF, ex)]);
fz[p] = d12dz * (
fh[idx_fh_F_ord4(iF, jF, kF - 2, ex)] -
EIT * fh[idx_fh_F_ord4(iF, jF, kF - 1, ex)] +
EIT * fh[idx_fh_F_ord4(iF, jF, kF + 1, ex)] -
fh[idx_fh_F_ord4(iF, jF, kF + 2, ex)]);
fh[idx_fh_F_ord5(iF, jF, kF - 2, ex)] -
EIT * fh[idx_fh_F_ord5(iF, jF, kF - 1, ex)] +
EIT * fh[idx_fh_F_ord5(iF, jF, kF + 1, ex)] -
fh[idx_fh_F_ord5(iF, jF, kF + 2, ex)]);
}
}
}
@@ -538,28 +538,28 @@ void fderivs(const int ex[3],
const size_t p = idx_ex(i0, j0, k0, ex);
fx[p] = d60dx * (
-fh[idx_fh_F_ord4(iF - 3, jF, kF, ex)] +
F9 * fh[idx_fh_F_ord4(iF - 2, jF, kF, ex)] -
F45 * fh[idx_fh_F_ord4(iF - 1, jF, kF, ex)] +
F45 * fh[idx_fh_F_ord4(iF + 1, jF, kF, ex)] -
F9 * fh[idx_fh_F_ord4(iF + 2, jF, kF, ex)] +
fh[idx_fh_F_ord4(iF + 3, jF, kF, ex)]);
-fh[idx_fh_F_ord5(iF - 3, jF, kF, ex)] +
F9 * fh[idx_fh_F_ord5(iF - 2, jF, kF, ex)] -
F45 * fh[idx_fh_F_ord5(iF - 1, jF, kF, ex)] +
F45 * fh[idx_fh_F_ord5(iF + 1, jF, kF, ex)] -
F9 * fh[idx_fh_F_ord5(iF + 2, jF, kF, ex)] +
fh[idx_fh_F_ord5(iF + 3, jF, kF, ex)]);
fy[p] = d60dy * (
-fh[idx_fh_F_ord4(iF, jF - 3, kF, ex)] +
F9 * fh[idx_fh_F_ord4(iF, jF - 2, kF, ex)] -
F45 * fh[idx_fh_F_ord4(iF, jF - 1, kF, ex)] +
F45 * fh[idx_fh_F_ord4(iF, jF + 1, kF, ex)] -
F9 * fh[idx_fh_F_ord4(iF, jF + 2, kF, ex)] +
fh[idx_fh_F_ord4(iF, jF + 3, kF, ex)]);
-fh[idx_fh_F_ord5(iF, jF - 3, kF, ex)] +
F9 * fh[idx_fh_F_ord5(iF, jF - 2, kF, ex)] -
F45 * fh[idx_fh_F_ord5(iF, jF - 1, kF, ex)] +
F45 * fh[idx_fh_F_ord5(iF, jF + 1, kF, ex)] -
F9 * fh[idx_fh_F_ord5(iF, jF + 2, kF, ex)] +
fh[idx_fh_F_ord5(iF, jF + 3, kF, ex)]);
fz[p] = d60dz * (
-fh[idx_fh_F_ord4(iF, jF, kF - 3, ex)] +
F9 * fh[idx_fh_F_ord4(iF, jF, kF - 2, ex)] -
F45 * fh[idx_fh_F_ord4(iF, jF, kF - 1, ex)] +
F45 * fh[idx_fh_F_ord4(iF, jF, kF + 1, ex)] -
F9 * fh[idx_fh_F_ord4(iF, jF, kF + 2, ex)] +
fh[idx_fh_F_ord4(iF, jF, kF + 3, ex)]);
-fh[idx_fh_F_ord5(iF, jF, kF - 3, ex)] +
F9 * fh[idx_fh_F_ord5(iF, jF, kF - 2, ex)] -
F45 * fh[idx_fh_F_ord5(iF, jF, kF - 1, ex)] +
F45 * fh[idx_fh_F_ord5(iF, jF, kF + 1, ex)] -
F9 * fh[idx_fh_F_ord5(iF, jF, kF + 2, ex)] +
fh[idx_fh_F_ord5(iF, jF, kF + 3, ex)]);
}
}
}
@@ -576,34 +576,34 @@ void fderivs(const int ex[3],
const size_t p = idx_ex(i0, j0, k0, ex);
fx[p] = d840dx * (
+(double)3 * fh[idx_fh_F_ord4(iF - 4, jF, kF, ex)] -
F32 * fh[idx_fh_F_ord4(iF - 3, jF, kF, ex)] +
F168 * fh[idx_fh_F_ord4(iF - 2, jF, kF, ex)] -
F672 * fh[idx_fh_F_ord4(iF - 1, jF, kF, ex)] +
F672 * fh[idx_fh_F_ord4(iF + 1, jF, kF, ex)] -
F168 * fh[idx_fh_F_ord4(iF + 2, jF, kF, ex)] +
F32 * fh[idx_fh_F_ord4(iF + 3, jF, kF, ex)] -
(double)3 * fh[idx_fh_F_ord4(iF + 4, jF, kF, ex)]);
+(double)3 * fh[idx_fh_F_ord5(iF - 4, jF, kF, ex)] -
F32 * fh[idx_fh_F_ord5(iF - 3, jF, kF, ex)] +
F168 * fh[idx_fh_F_ord5(iF - 2, jF, kF, ex)] -
F672 * fh[idx_fh_F_ord5(iF - 1, jF, kF, ex)] +
F672 * fh[idx_fh_F_ord5(iF + 1, jF, kF, ex)] -
F168 * fh[idx_fh_F_ord5(iF + 2, jF, kF, ex)] +
F32 * fh[idx_fh_F_ord5(iF + 3, jF, kF, ex)] -
(double)3 * fh[idx_fh_F_ord5(iF + 4, jF, kF, ex)]);
fy[p] = d840dy * (
+(double)3 * fh[idx_fh_F_ord4(iF, jF - 4, kF, ex)] -
F32 * fh[idx_fh_F_ord4(iF, jF - 3, kF, ex)] +
F168 * fh[idx_fh_F_ord4(iF, jF - 2, kF, ex)] -
F672 * fh[idx_fh_F_ord4(iF, jF - 1, kF, ex)] +
F672 * fh[idx_fh_F_ord4(iF, jF + 1, kF, ex)] -
F168 * fh[idx_fh_F_ord4(iF, jF + 2, kF, ex)] +
F32 * fh[idx_fh_F_ord4(iF, jF + 3, kF, ex)] -
(double)3 * fh[idx_fh_F_ord4(iF, jF + 4, kF, ex)]);
+(double)3 * fh[idx_fh_F_ord5(iF, jF - 4, kF, ex)] -
F32 * fh[idx_fh_F_ord5(iF, jF - 3, kF, ex)] +
F168 * fh[idx_fh_F_ord5(iF, jF - 2, kF, ex)] -
F672 * fh[idx_fh_F_ord5(iF, jF - 1, kF, ex)] +
F672 * fh[idx_fh_F_ord5(iF, jF + 1, kF, ex)] -
F168 * fh[idx_fh_F_ord5(iF, jF + 2, kF, ex)] +
F32 * fh[idx_fh_F_ord5(iF, jF + 3, kF, ex)] -
(double)3 * fh[idx_fh_F_ord5(iF, jF + 4, kF, ex)]);
fz[p] = d840dz * (
+(double)3 * fh[idx_fh_F_ord4(iF, jF, kF - 4, ex)] -
F32 * fh[idx_fh_F_ord4(iF, jF, kF - 3, ex)] +
F168 * fh[idx_fh_F_ord4(iF, jF, kF - 2, ex)] -
F672 * fh[idx_fh_F_ord4(iF, jF, kF - 1, ex)] +
F672 * fh[idx_fh_F_ord4(iF, jF, kF + 1, ex)] -
F168 * fh[idx_fh_F_ord4(iF, jF, kF + 2, ex)] +
F32 * fh[idx_fh_F_ord4(iF, jF, kF + 3, ex)] -
(double)3 * fh[idx_fh_F_ord4(iF, jF, kF + 4, ex)]);
+(double)3 * fh[idx_fh_F_ord5(iF, jF, kF - 4, ex)] -
F32 * fh[idx_fh_F_ord5(iF, jF, kF - 3, ex)] +
F168 * fh[idx_fh_F_ord5(iF, jF, kF - 2, ex)] -
F672 * fh[idx_fh_F_ord5(iF, jF, kF - 1, ex)] +
F672 * fh[idx_fh_F_ord5(iF, jF, kF + 1, ex)] -
F168 * fh[idx_fh_F_ord5(iF, jF, kF + 2, ex)] +
F32 * fh[idx_fh_F_ord5(iF, jF, kF + 3, ex)] -
(double)3 * fh[idx_fh_F_ord5(iF, jF, kF + 4, ex)]);
}
}
}