From b489564311f5f2c83be3b3e3f0311010877f1757 Mon Sep 17 00:00:00 2001 From: CGH0S7 <776459475@qq.com> Date: Thu, 14 May 2026 16:00:11 +0800 Subject: [PATCH] Fix C derivative ghost-buffer indexing across FD orders --- AMSS_NCKU_source/fdderivs_c.C | 90 ++++++++++++++--------------- AMSS_NCKU_source/fdderivs_sh_c.C | 20 +++---- AMSS_NCKU_source/fderivs_c.C | 90 ++++++++++++++--------------- AMSS_NCKU_source/fderivs_sh_c.C | 20 +++---- AMSS_NCKU_source/kodiss_c.C | 12 ++-- AMSS_NCKU_source/kodiss_sh_c.C | 8 +-- AMSS_NCKU_source/lopsided_c.C | 12 ++-- AMSS_NCKU_source/lopsided_kodis_c.C | 24 ++++---- AMSS_NCKU_source/share_func.h | 4 +- 9 files changed, 140 insertions(+), 140 deletions(-) diff --git a/AMSS_NCKU_source/fdderivs_c.C b/AMSS_NCKU_source/fdderivs_c.C index ea9d069..05db174 100644 --- a/AMSS_NCKU_source/fdderivs_c.C +++ b/AMSS_NCKU_source/fdderivs_c.C @@ -93,11 +93,11 @@ void fdderivs(const int ex[3], if (i2_lo <= i2_hi && j2_lo <= j2_hi && k2_lo <= k2_hi) { for (int k0 = k2_lo; k0 <= k2_hi; ++k0) { - const int kF = k0; + const int kF = k0 + 1; for (int j0 = j2_lo; j0 <= j2_hi; ++j0) { - const int jF = j0; + const int jF = j0 + 1; for (int i0 = i2_lo; i0 <= i2_hi; ++i0) { - const int iF = i0; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fxx[p] = Sdxdx * ( @@ -471,19 +471,19 @@ void fdderivs(const int ex[3], } /* loop bounds for each pass (from widest to narrowest) */ - const int i2_lo = (iminF + 1 > 0) ? (iminF + 1) : 0; - const int j2_lo = (jminF + 1 > 0) ? (jminF + 1) : 0; - const int k2_lo = (kminF + 1 > 0) ? (kminF + 1) : 0; + const int i2_lo = (iminF > 0) ? iminF : 0; + const int j2_lo = (jminF > 0) ? jminF : 0; + const int k2_lo = (kminF > 0) ? kminF : 0; const int i2_hi = ex1 - 2, j2_hi = ex2 - 2, k2_hi = ex3 - 2; - const int i4_lo = (iminF + 2 > 0) ? (iminF + 2) : 0; - const int j4_lo = (jminF + 2 > 0) ? (jminF + 2) : 0; - const int k4_lo = (kminF + 2 > 0) ? (kminF + 2) : 0; + const int i4_lo = (iminF + 1 > 0) ? (iminF + 1) : 0; + const int j4_lo = (jminF + 1 > 0) ? (jminF + 1) : 0; + const int k4_lo = (kminF + 1 > 0) ? (kminF + 1) : 0; const int i4_hi = ex1 - 3, j4_hi = ex2 - 3, k4_hi = ex3 - 3; - const int i6_lo = (iminF + 3 > 0) ? (iminF + 3) : 0; - const int j6_lo = (jminF + 3 > 0) ? (jminF + 3) : 0; - const int k6_lo = (kminF + 3 > 0) ? (kminF + 3) : 0; + const int i6_lo = (iminF + 2 > 0) ? (iminF + 2) : 0; + const int j6_lo = (jminF + 2 > 0) ? (jminF + 2) : 0; + const int k6_lo = (kminF + 2 > 0) ? (kminF + 2) : 0; const int i6_hi = ex1 - 4, j6_hi = ex2 - 4, k6_hi = ex3 - 4; const int has4 = (i4_lo <= i4_hi && j4_lo <= j4_hi && k4_lo <= k4_hi); @@ -492,13 +492,13 @@ void fdderivs(const int ex[3], /* 2nd-order: skip 4th+6th overlap */ if (i2_lo <= i2_hi && j2_lo <= j2_hi && k2_lo <= k2_hi) { for (int k0 = k2_lo; k0 <= k2_hi; ++k0) { - const int kF = k0 + 2; + const int kF = k0 + 1; for (int j0 = j2_lo; j0 <= j2_hi; ++j0) { - const int jF = j0 + 2; + const int jF = j0 + 1; for (int i0 = i2_lo; i0 <= i2_hi; ++i0) { _Bool in4 = has4 && i0>=i4_lo && i0<=i4_hi && j0>=j4_lo && j0<=j4_hi && k0>=k4_lo && k0<=k4_hi; if (in4) continue; - const int iF = i0 + 2; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fxx[p] = Sdxdx * (fh[idx_fh_F(iF - 1, jF, kF, ex)] - TWO*fh[idx_fh_F(iF,jF,kF,ex)] + fh[idx_fh_F(iF + 1, jF, kF, ex)]); @@ -516,12 +516,12 @@ void fdderivs(const int ex[3], /* 4th-order: skip 6th overlap */ if (has4) { for (int k0 = k4_lo; k0 <= k4_hi; ++k0) { - const int kF = k0 + 2; + const int kF = k0 + 1; for (int j0 = j4_lo; j0 <= j4_hi; ++j0) { - const int jF = j0 + 2; + const int jF = j0 + 1; for (int i0 = i4_lo; i0 <= i4_hi; ++i0) { if (has6 && i0>=i6_lo && i0<=i6_hi && j0>=j6_lo && j0<=j6_hi && k0>=k6_lo && k0<=k6_hi) continue; - const int iF = i0 + 2; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fxx[p] = Fdxdx * (-fh[idx_fh_F(iF - 2, jF, kF, ex)] + F16*fh[idx_fh_F(iF-1,jF,kF,ex)] - F30*fh[idx_fh_F(iF,jF,kF,ex)] - fh[idx_fh_F(iF+2,jF,kF,ex)] + F16*fh[idx_fh_F(iF+1,jF,kF,ex)]); @@ -557,11 +557,11 @@ void fdderivs(const int ex[3], /* 6th-order: interior only */ if (has6) { for (int k0 = k6_lo; k0 <= k6_hi; ++k0) { - const int kF = k0 + 2; + const int kF = k0 + 1; for (int j0 = j6_lo; j0 <= j6_hi; ++j0) { - const int jF = j0 + 2; + const int jF = j0 + 1; for (int i0 = i6_lo; i0 <= i6_hi; ++i0) { - const int iF = i0 + 2; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); /* Diagonal: [+2,-27,+270,-490,+270,-27,+2] / (180*dx^2) */ @@ -685,24 +685,24 @@ void fdderivs(const int ex[3], } /* Loop bounds for each pass */ - const int i2_lo = (iminF + 2 > 0) ? (iminF + 2) : 0; - const int j2_lo = (jminF + 2 > 0) ? (jminF + 2) : 0; - const int k2_lo = (kminF + 2 > 0) ? (kminF + 2) : 0; + const int i2_lo = (iminF > 0) ? iminF : 0; + const int j2_lo = (jminF > 0) ? jminF : 0; + const int k2_lo = (kminF > 0) ? kminF : 0; const int i2_hi = ex1 - 2, j2_hi = ex2 - 2, k2_hi = ex3 - 2; - const int i4_lo = (iminF + 3 > 0) ? (iminF + 3) : 0; - const int j4_lo = (jminF + 3 > 0) ? (jminF + 3) : 0; - const int k4_lo = (kminF + 3 > 0) ? (kminF + 3) : 0; + const int i4_lo = (iminF + 1 > 0) ? (iminF + 1) : 0; + const int j4_lo = (jminF + 1 > 0) ? (jminF + 1) : 0; + const int k4_lo = (kminF + 1 > 0) ? (kminF + 1) : 0; const int i4_hi = ex1 - 3, j4_hi = ex2 - 3, k4_hi = ex3 - 3; - const int i6_lo = (iminF + 4 > 0) ? (iminF + 4) : 0; - const int j6_lo = (jminF + 4 > 0) ? (jminF + 4) : 0; - const int k6_lo = (kminF + 4 > 0) ? (kminF + 4) : 0; + const int i6_lo = (iminF + 2 > 0) ? (iminF + 2) : 0; + const int j6_lo = (jminF + 2 > 0) ? (jminF + 2) : 0; + const int k6_lo = (kminF + 2 > 0) ? (kminF + 2) : 0; const int i6_hi = ex1 - 4, j6_hi = ex2 - 4, k6_hi = ex3 - 4; - const int i8_lo = (iminF + 5 > 0) ? (iminF + 5) : 0; - const int j8_lo = (jminF + 5 > 0) ? (jminF + 5) : 0; - const int k8_lo = (kminF + 5 > 0) ? (kminF + 5) : 0; + const int i8_lo = (iminF + 3 > 0) ? (iminF + 3) : 0; + const int j8_lo = (jminF + 3 > 0) ? (jminF + 3) : 0; + const int k8_lo = (kminF + 3 > 0) ? (kminF + 3) : 0; const int i8_hi = ex1 - 5, j8_hi = ex2 - 5, k8_hi = ex3 - 5; const int has4 = (i4_lo <= i4_hi && j4_lo <= j4_hi && k4_lo <= k4_hi); @@ -712,13 +712,13 @@ void fdderivs(const int ex[3], /* 2nd-order: skip 4th+6th+8th overlap */ if (i2_lo <= i2_hi && j2_lo <= j2_hi && k2_lo <= k2_hi) { for (int k0 = k2_lo; k0 <= k2_hi; ++k0) { - const int kF = k0 + 3; + const int kF = k0 + 1; for (int j0 = j2_lo; j0 <= j2_hi; ++j0) { - const int jF = j0 + 3; + const int jF = j0 + 1; for (int i0 = i2_lo; i0 <= i2_hi; ++i0) { _Bool in4 = has4 && i0>=i4_lo && i0<=i4_hi && j0>=j4_lo && j0<=j4_hi && k0>=k4_lo && k0<=k4_hi; if (in4) continue; - const int iF = i0 + 3; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fxx[p] = Sdxdx * (fh[idx_fh_F_ord4(iF-1,jF,kF,ex)] - TWO*fh[idx_fh_F_ord4(iF,jF,kF,ex)] + fh[idx_fh_F_ord4(iF+1,jF,kF,ex)]); @@ -736,13 +736,13 @@ void fdderivs(const int ex[3], /* 4th-order: skip 6th+8th overlap */ if (has4) { for (int k0 = k4_lo; k0 <= k4_hi; ++k0) { - const int kF = k0 + 3; + const int kF = k0 + 1; for (int j0 = j4_lo; j0 <= j4_hi; ++j0) { - const int jF = j0 + 3; + const int jF = j0 + 1; for (int i0 = i4_lo; i0 <= i4_hi; ++i0) { _Bool in6 = has6 && i0>=i6_lo && i0<=i6_hi && j0>=j6_lo && j0<=j6_hi && k0>=k6_lo && k0<=k6_hi; if (in6) continue; - const int iF = i0 + 3; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fxx[p] = Fdxdx * (-fh[idx_fh_F_ord4(iF-2,jF,kF,ex)] + F16*fh[idx_fh_F_ord4(iF-1,jF,kF,ex)] - F30*fh[idx_fh_F_ord4(iF,jF,kF,ex)] - fh[idx_fh_F_ord4(iF+2,jF,kF,ex)] + F16*fh[idx_fh_F_ord4(iF+1,jF,kF,ex)]); @@ -778,12 +778,12 @@ void fdderivs(const int ex[3], /* 6th-order: skip 8th overlap */ if (has6) { for (int k0 = k6_lo; k0 <= k6_hi; ++k0) { - const int kF = k0 + 3; + const int kF = k0 + 1; for (int j0 = j6_lo; j0 <= j6_hi; ++j0) { - const int jF = j0 + 3; + const int jF = j0 + 1; for (int i0 = i6_lo; i0 <= i6_hi; ++i0) { if (has8 && i0>=i8_lo && i0<=i8_hi && j0>=j8_lo && j0<=j8_hi && k0>=k8_lo && k0<=k8_hi) continue; - const int iF = i0 + 3; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fxx[p] = Xdxdx * ( @@ -817,11 +817,11 @@ void fdderivs(const int ex[3], /* 8th-order: interior only */ if (has8) { for (int k0 = k8_lo; k0 <= k8_hi; ++k0) { - const int kF = k0 + 3; + const int kF = k0 + 1; for (int j0 = j8_lo; j0 <= j8_hi; ++j0) { - const int jF = j0 + 3; + const int jF = j0 + 1; for (int i0 = i8_lo; i0 <= i8_hi; ++i0) { - const int iF = i0 + 3; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); /* Diagonal: [-9,+128,-1008,+8064,-14350,+8064,-1008,+128,-9] / (5040*dx^2) */ diff --git a/AMSS_NCKU_source/fdderivs_sh_c.C b/AMSS_NCKU_source/fdderivs_sh_c.C index 497ec93..955319d 100644 --- a/AMSS_NCKU_source/fdderivs_sh_c.C +++ b/AMSS_NCKU_source/fdderivs_sh_c.C @@ -48,7 +48,7 @@ extern "C" void fdderivs_sh_(const int ex[3], const size_t all=(size_t)ex1*ex2*ex3; for(size_t p=0;p0)?iminF:0,j2_lo=(jminF>0)?jminF:0,k2_lo=0,i2_hi=ex1-2,j2_hi=ex2-2,k2_hi=ex3-2; + const int i2_lo=(iminF>0)?iminF:0,j2_lo=(jminF>0)?jminF:0,k2_lo=1,i2_hi=ex1-2,j2_hi=ex2-2,k2_hi=ex3-2; #define FH(iF,jF,kF) fh[idx_fh_stbd(iF,jF,kF,ord,ex)] if(i2_lo<=i2_hi&&j2_lo<=j2_hi&&k2_lo<=k2_hi){ for(int k0=k2_lo;k0<=k2_hi;++k0){const int kF=k0+1; @@ -86,8 +86,8 @@ extern "C" void fdderivs_sh_(const int ex[3], const size_t all=(size_t)ex1*ex2*ex3; for(size_t p=0;p0)?iminF:0,j2_lo=(jminF>0)?jminF:0,k2_lo=0,i2_hi=ex1-2,j2_hi=ex2-2,k2_hi=ex3-2; - const int i4_lo=(iminF+1>0)?iminF+1:0,j4_lo=(jminF+1>0)?jminF+1:0,k4_lo=0,i4_hi=ex1-3,j4_hi=ex2-3,k4_hi=ex3-3; + const int i2_lo=(iminF>0)?iminF:0,j2_lo=(jminF>0)?jminF:0,k2_lo=1,i2_hi=ex1-2,j2_hi=ex2-2,k2_hi=ex3-2; + const int i4_lo=(iminF+1>0)?iminF+1:0,j4_lo=(jminF+1>0)?jminF+1:0,k4_lo=2,i4_hi=ex1-3,j4_hi=ex2-3,k4_hi=ex3-3; const int has4=(i4_lo<=i4_hi&&j4_lo<=j4_hi&&k4_lo<=k4_hi); #define FH(iF,jF,kF) fh[idx_fh_stbd(iF,jF,kF,ord,ex)] @@ -155,9 +155,9 @@ extern "C" void fdderivs_sh_(const int ex[3], const size_t all=(size_t)ex1*ex2*ex3; for(size_t p=0;p0)?iminF+1:0,j2_lo=(jminF+1>0)?jminF+1:0,k2_lo=0,i2_hi=ex1-2,j2_hi=ex2-2,k2_hi=ex3-2; - const int i4_lo=(iminF+2>0)?iminF+2:0,j4_lo=(jminF+2>0)?jminF+2:0,k4_lo=0,i4_hi=ex1-3,j4_hi=ex2-3,k4_hi=ex3-3; - const int i6_lo=(iminF+3>0)?iminF+3:0,j6_lo=(jminF+3>0)?jminF+3:0,k6_lo=0,i6_hi=ex1-4,j6_hi=ex2-4,k6_hi=ex3-4; + const int i2_lo=(iminF>0)?iminF:0,j2_lo=(jminF>0)?jminF:0,k2_lo=1,i2_hi=ex1-2,j2_hi=ex2-2,k2_hi=ex3-2; + const int i4_lo=(iminF+1>0)?iminF+1:0,j4_lo=(jminF+1>0)?jminF+1:0,k4_lo=2,i4_hi=ex1-3,j4_hi=ex2-3,k4_hi=ex3-3; + const int i6_lo=(iminF+2>0)?iminF+2:0,j6_lo=(jminF+2>0)?jminF+2:0,k6_lo=3,i6_hi=ex1-4,j6_hi=ex2-4,k6_hi=ex3-4; const int has4=(i4_lo<=i4_hi&&j4_lo<=j4_hi&&k4_lo<=k4_hi),has6=(i6_lo<=i6_hi&&j6_lo<=j6_hi&&k6_lo<=k6_hi); #define FH(iF,jF,kF) fh[idx_fh_stbd(iF,jF,kF,ord,ex)] @@ -238,10 +238,10 @@ extern "C" void fdderivs_sh_(const int ex[3], const size_t all=(size_t)ex1*ex2*ex3; for(size_t p=0;p0)?iminF+2:0,j2_lo=(jminF+2>0)?jminF+2:0,k2_lo=0,i2_hi=ex1-2,j2_hi=ex2-2,k2_hi=ex3-2; - const int i4_lo=(iminF+3>0)?iminF+3:0,j4_lo=(jminF+3>0)?jminF+3:0,k4_lo=0,i4_hi=ex1-3,j4_hi=ex2-3,k4_hi=ex3-3; - const int i6_lo=(iminF+4>0)?iminF+4:0,j6_lo=(jminF+4>0)?jminF+4:0,k6_lo=0,i6_hi=ex1-4,j6_hi=ex2-4,k6_hi=ex3-4; - const int i8_lo=(iminF+5>0)?iminF+5:0,j8_lo=(jminF+5>0)?jminF+5:0,k8_lo=0,i8_hi=ex1-5,j8_hi=ex2-5,k8_hi=ex3-5; + const int i2_lo=(iminF>0)?iminF:0,j2_lo=(jminF>0)?jminF:0,k2_lo=1,i2_hi=ex1-2,j2_hi=ex2-2,k2_hi=ex3-2; + const int i4_lo=(iminF+1>0)?iminF+1:0,j4_lo=(jminF+1>0)?jminF+1:0,k4_lo=2,i4_hi=ex1-3,j4_hi=ex2-3,k4_hi=ex3-3; + const int i6_lo=(iminF+2>0)?iminF+2:0,j6_lo=(jminF+2>0)?jminF+2:0,k6_lo=3,i6_hi=ex1-4,j6_hi=ex2-4,k6_hi=ex3-4; + const int i8_lo=(iminF+3>0)?iminF+3:0,j8_lo=(jminF+3>0)?jminF+3:0,k8_lo=4,i8_hi=ex1-5,j8_hi=ex2-5,k8_hi=ex3-5; const int has4=(i4_lo<=i4_hi&&j4_lo<=j4_hi&&k4_lo<=k4_hi),has6=(i6_lo<=i6_hi&&j6_lo<=j6_hi&&k6_lo<=k6_hi),has8=(i8_lo<=i8_hi&&j8_lo<=j8_hi&&k8_lo<=k8_hi); #define FH(iF,jF,kF) fh[idx_fh_stbd(iF,jF,kF,ord,ex)] diff --git a/AMSS_NCKU_source/fderivs_c.C b/AMSS_NCKU_source/fderivs_c.C index b61c7f8..baa2eb6 100644 --- a/AMSS_NCKU_source/fderivs_c.C +++ b/AMSS_NCKU_source/fderivs_c.C @@ -86,11 +86,11 @@ void fderivs(const int ex[3], if (i2_lo <= i2_hi && j2_lo <= j2_hi && k2_lo <= k2_hi) { for (int k0 = k2_lo; k0 <= k2_hi; ++k0) { - const int kF = k0; // Fortran index for fh: kF = k0 (shift=0) + const int kF = k0 + 1; for (int j0 = j2_lo; j0 <= j2_hi; ++j0) { - const int jF = j0; + const int jF = j0 + 1; for (int i0 = i2_lo; i0 <= i2_hi; ++i0) { - const int iF = i0; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fx[p] = d2dx * ( @@ -277,25 +277,25 @@ void fderivs(const int ex[3], } /* 2nd-order pass: 3pt, widest */ - const int i2_lo = (iminF + 1 > 0) ? (iminF + 1) : 0; - const int j2_lo = (jminF + 1 > 0) ? (jminF + 1) : 0; - const int k2_lo = (kminF + 1 > 0) ? (kminF + 1) : 0; + const int i2_lo = (iminF > 0) ? iminF : 0; + const int j2_lo = (jminF > 0) ? jminF : 0; + const int k2_lo = (kminF > 0) ? kminF : 0; const int i2_hi = ex1 - 2; const int j2_hi = ex2 - 2; const int k2_hi = ex3 - 2; /* 4th-order pass: 5pt */ - const int i4_lo = (iminF + 2 > 0) ? (iminF + 2) : 0; - const int j4_lo = (jminF + 2 > 0) ? (jminF + 2) : 0; - const int k4_lo = (kminF + 2 > 0) ? (kminF + 2) : 0; + const int i4_lo = (iminF + 1 > 0) ? (iminF + 1) : 0; + const int j4_lo = (jminF + 1 > 0) ? (jminF + 1) : 0; + const int k4_lo = (kminF + 1 > 0) ? (kminF + 1) : 0; const int i4_hi = ex1 - 3; const int j4_hi = ex2 - 3; const int k4_hi = ex3 - 3; /* 6th-order pass: 7pt, narrowest interior */ - const int i6_lo = (iminF + 3 > 0) ? (iminF + 3) : 0; - const int j6_lo = (jminF + 3 > 0) ? (jminF + 3) : 0; - const int k6_lo = (kminF + 3 > 0) ? (kminF + 3) : 0; + const int i6_lo = (iminF + 2 > 0) ? (iminF + 2) : 0; + const int j6_lo = (jminF + 2 > 0) ? (jminF + 2) : 0; + const int k6_lo = (kminF + 2 > 0) ? (kminF + 2) : 0; const int i6_hi = ex1 - 4; const int j6_hi = ex2 - 4; const int k6_hi = ex3 - 4; @@ -303,11 +303,11 @@ void fderivs(const int ex[3], /* 2nd-order */ if (i2_lo <= i2_hi && j2_lo <= j2_hi && k2_lo <= k2_hi) { for (int k0 = k2_lo; k0 <= k2_hi; ++k0) { - const int kF = k0 + 2; + const int kF = k0 + 1; for (int j0 = j2_lo; j0 <= j2_hi; ++j0) { - const int jF = j0 + 2; + const int jF = j0 + 1; for (int i0 = i2_lo; i0 <= i2_hi; ++i0) { - const int iF = i0 + 2; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fx[p] = d2dx * ( @@ -327,11 +327,11 @@ void fderivs(const int ex[3], /* 4th-order overwrite */ if (i4_lo <= i4_hi && j4_lo <= j4_hi && k4_lo <= k4_hi) { for (int k0 = k4_lo; k0 <= k4_hi; ++k0) { - const int kF = k0 + 2; + const int kF = k0 + 1; for (int j0 = j4_lo; j0 <= j4_hi; ++j0) { - const int jF = j0 + 2; + const int jF = j0 + 1; for (int i0 = i4_lo; i0 <= i4_hi; ++i0) { - const int iF = i0 + 2; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fx[p] = d12dx * ( @@ -359,11 +359,11 @@ void fderivs(const int ex[3], /* 6th-order overwrite: [-1,+9,-45,0,+45,-9,+1] / (60*dx) */ if (i6_lo <= i6_hi && j6_lo <= j6_hi && k6_lo <= k6_hi) { for (int k0 = k6_lo; k0 <= k6_hi; ++k0) { - const int kF = k0 + 2; + const int kF = k0 + 1; for (int j0 = j6_lo; j0 <= j6_hi; ++j0) { - const int jF = j0 + 2; + const int jF = j0 + 1; for (int i0 = i6_lo; i0 <= i6_hi; ++i0) { - const int iF = i0 + 2; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fx[p] = d60dx * ( @@ -443,44 +443,44 @@ void fderivs(const int ex[3], } /* 2nd: 3pt, widest */ - const int i2_lo = (iminF + 2 > 0) ? (iminF + 2) : 0; - const int j2_lo = (jminF + 2 > 0) ? (jminF + 2) : 0; - const int k2_lo = (kminF + 2 > 0) ? (kminF + 2) : 0; + const int i2_lo = (iminF > 0) ? iminF : 0; + const int j2_lo = (jminF > 0) ? jminF : 0; + const int k2_lo = (kminF > 0) ? kminF : 0; const int i2_hi = ex1 - 2; const int j2_hi = ex2 - 2; const int k2_hi = ex3 - 2; /* 4th: 5pt */ - const int i4_lo = (iminF + 3 > 0) ? (iminF + 3) : 0; - const int j4_lo = (jminF + 3 > 0) ? (jminF + 3) : 0; - const int k4_lo = (kminF + 3 > 0) ? (kminF + 3) : 0; + const int i4_lo = (iminF + 1 > 0) ? (iminF + 1) : 0; + const int j4_lo = (jminF + 1 > 0) ? (jminF + 1) : 0; + const int k4_lo = (kminF + 1 > 0) ? (kminF + 1) : 0; const int i4_hi = ex1 - 3; const int j4_hi = ex2 - 3; const int k4_hi = ex3 - 3; /* 6th: 7pt */ - const int i6_lo = (iminF + 4 > 0) ? (iminF + 4) : 0; - const int j6_lo = (jminF + 4 > 0) ? (jminF + 4) : 0; - const int k6_lo = (kminF + 4 > 0) ? (kminF + 4) : 0; + const int i6_lo = (iminF + 2 > 0) ? (iminF + 2) : 0; + const int j6_lo = (jminF + 2 > 0) ? (jminF + 2) : 0; + const int k6_lo = (kminF + 2 > 0) ? (kminF + 2) : 0; const int i6_hi = ex1 - 4; const int j6_hi = ex2 - 4; const int k6_hi = ex3 - 4; /* 8th: 9pt, narrowest */ - const int i8_lo = (iminF + 5 > 0) ? (iminF + 5) : 0; - const int j8_lo = (jminF + 5 > 0) ? (jminF + 5) : 0; - const int k8_lo = (kminF + 5 > 0) ? (kminF + 5) : 0; + const int i8_lo = (iminF + 3 > 0) ? (iminF + 3) : 0; + const int j8_lo = (jminF + 3 > 0) ? (jminF + 3) : 0; + const int k8_lo = (kminF + 3 > 0) ? (kminF + 3) : 0; const int i8_hi = ex1 - 5; const int j8_hi = ex2 - 5; const int k8_hi = ex3 - 5; if (i2_lo <= i2_hi && j2_lo <= j2_hi && k2_lo <= k2_hi) { for (int k0 = k2_lo; k0 <= k2_hi; ++k0) { - const int kF = k0 + 3; + const int kF = k0 + 1; for (int j0 = j2_lo; j0 <= j2_hi; ++j0) { - const int jF = j0 + 3; + const int jF = j0 + 1; for (int i0 = i2_lo; i0 <= i2_hi; ++i0) { - const int iF = i0 + 3; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fx[p] = d2dx * ( @@ -499,11 +499,11 @@ void fderivs(const int ex[3], if (i4_lo <= i4_hi && j4_lo <= j4_hi && k4_lo <= k4_hi) { for (int k0 = k4_lo; k0 <= k4_hi; ++k0) { - const int kF = k0 + 3; + const int kF = k0 + 1; for (int j0 = j4_lo; j0 <= j4_hi; ++j0) { - const int jF = j0 + 3; + const int jF = j0 + 1; for (int i0 = i4_lo; i0 <= i4_hi; ++i0) { - const int iF = i0 + 3; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fx[p] = d12dx * ( @@ -530,11 +530,11 @@ void fderivs(const int ex[3], if (i6_lo <= i6_hi && j6_lo <= j6_hi && k6_lo <= k6_hi) { for (int k0 = k6_lo; k0 <= k6_hi; ++k0) { - const int kF = k0 + 3; + const int kF = k0 + 1; for (int j0 = j6_lo; j0 <= j6_hi; ++j0) { - const int jF = j0 + 3; + const int jF = j0 + 1; for (int i0 = i6_lo; i0 <= i6_hi; ++i0) { - const int iF = i0 + 3; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fx[p] = d60dx * ( @@ -568,11 +568,11 @@ void fderivs(const int ex[3], /* 8th-order overwrite: [+3,-32,+168,-672,0,+672,-168,+32,-3] / (840*dx) */ if (i8_lo <= i8_hi && j8_lo <= j8_hi && k8_lo <= k8_hi) { for (int k0 = k8_lo; k0 <= k8_hi; ++k0) { - const int kF = k0 + 3; + const int kF = k0 + 1; for (int j0 = j8_lo; j0 <= j8_hi; ++j0) { - const int jF = j0 + 3; + const int jF = j0 + 1; for (int i0 = i8_lo; i0 <= i8_hi; ++i0) { - const int iF = i0 + 3; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); fx[p] = d840dx * ( diff --git a/AMSS_NCKU_source/fderivs_sh_c.C b/AMSS_NCKU_source/fderivs_sh_c.C index 0a28f23..f82051e 100644 --- a/AMSS_NCKU_source/fderivs_sh_c.C +++ b/AMSS_NCKU_source/fderivs_sh_c.C @@ -52,7 +52,7 @@ extern "C" void fderivs_sh_(const int ex[3], const size_t all = (size_t)ex1*ex2*ex3; for (size_t p=0;p0)?iminF:0, j2_lo=(jminF>0)?jminF:0, k2_lo=0; + const int i2_lo=(iminF>0)?iminF:0, j2_lo=(jminF>0)?jminF:0, k2_lo=1; const int i2_hi=ex1-2, j2_hi=ex2-2, k2_hi=ex3-2; if (i2_lo<=i2_hi&&j2_lo<=j2_hi&&k2_lo<=k2_hi) { for (int k0=k2_lo;k0<=k2_hi;++k0) { const int kF=k0+1; @@ -86,9 +86,9 @@ extern "C" void fderivs_sh_(const int ex[3], const size_t all=(size_t)ex1*ex2*ex3; for(size_t p=0;p0)?iminF:0, j2_lo=(jminF>0)?jminF:0, k2_lo=0; + const int i2_lo=(iminF>0)?iminF:0, j2_lo=(jminF>0)?jminF:0, k2_lo=1; const int i2_hi=ex1-2, j2_hi=ex2-2, k2_hi=ex3-2; - const int i4_lo=(iminF+1>0)?iminF+1:0, j4_lo=(jminF+1>0)?jminF+1:0, k4_lo=0; + const int i4_lo=(iminF+1>0)?iminF+1:0, j4_lo=(jminF+1>0)?jminF+1:0, k4_lo=2; const int i4_hi=ex1-3, j4_hi=ex2-3, k4_hi=ex3-3; if (i2_lo<=i2_hi&&j2_lo<=j2_hi&&k2_lo<=k2_hi) { @@ -134,9 +134,9 @@ extern "C" void fderivs_sh_(const int ex[3], const size_t all=(size_t)ex1*ex2*ex3; for(size_t p=0;p0)?iminF+1:0,j2_lo=(jminF+1>0)?jminF+1:0,k2_lo=0,i2_hi=ex1-2,j2_hi=ex2-2,k2_hi=ex3-2; - const int i4_lo=(iminF+2>0)?iminF+2:0,j4_lo=(jminF+2>0)?jminF+2:0,k4_lo=0,i4_hi=ex1-3,j4_hi=ex2-3,k4_hi=ex3-3; - const int i6_lo=(iminF+3>0)?iminF+3:0,j6_lo=(jminF+3>0)?jminF+3:0,k6_lo=0,i6_hi=ex1-4,j6_hi=ex2-4,k6_hi=ex3-4; + const int i2_lo=(iminF>0)?iminF:0,j2_lo=(jminF>0)?jminF:0,k2_lo=1,i2_hi=ex1-2,j2_hi=ex2-2,k2_hi=ex3-2; + const int i4_lo=(iminF+1>0)?iminF+1:0,j4_lo=(jminF+1>0)?jminF+1:0,k4_lo=2,i4_hi=ex1-3,j4_hi=ex2-3,k4_hi=ex3-3; + const int i6_lo=(iminF+2>0)?iminF+2:0,j6_lo=(jminF+2>0)?jminF+2:0,k6_lo=3,i6_hi=ex1-4,j6_hi=ex2-4,k6_hi=ex3-4; if(i2_lo<=i2_hi&&j2_lo<=j2_hi&&k2_lo<=k2_hi){ for(int k0=k2_lo;k0<=k2_hi;++k0){const int kF=k0+1; @@ -192,10 +192,10 @@ extern "C" void fderivs_sh_(const int ex[3], const size_t all=(size_t)ex1*ex2*ex3; for(size_t p=0;p0)?iminF+2:0,j2_lo=(jminF+2>0)?jminF+2:0,k2_lo=0,i2_hi=ex1-2,j2_hi=ex2-2,k2_hi=ex3-2; - const int i4_lo=(iminF+3>0)?iminF+3:0,j4_lo=(jminF+3>0)?jminF+3:0,k4_lo=0,i4_hi=ex1-3,j4_hi=ex2-3,k4_hi=ex3-3; - const int i6_lo=(iminF+4>0)?iminF+4:0,j6_lo=(jminF+4>0)?jminF+4:0,k6_lo=0,i6_hi=ex1-4,j6_hi=ex2-4,k6_hi=ex3-4; - const int i8_lo=(iminF+5>0)?iminF+5:0,j8_lo=(jminF+5>0)?jminF+5:0,k8_lo=0,i8_hi=ex1-5,j8_hi=ex2-5,k8_hi=ex3-5; + const int i2_lo=(iminF>0)?iminF:0,j2_lo=(jminF>0)?jminF:0,k2_lo=1,i2_hi=ex1-2,j2_hi=ex2-2,k2_hi=ex3-2; + const int i4_lo=(iminF+1>0)?iminF+1:0,j4_lo=(jminF+1>0)?jminF+1:0,k4_lo=2,i4_hi=ex1-3,j4_hi=ex2-3,k4_hi=ex3-3; + const int i6_lo=(iminF+2>0)?iminF+2:0,j6_lo=(jminF+2>0)?jminF+2:0,k6_lo=3,i6_hi=ex1-4,j6_hi=ex2-4,k6_hi=ex3-4; + const int i8_lo=(iminF+3>0)?iminF+3:0,j8_lo=(jminF+3>0)?jminF+3:0,k8_lo=4,i8_hi=ex1-5,j8_hi=ex2-5,k8_hi=ex3-5; #define FH_S(iF,jF,kF) fh[idx_fh_stbd(iF,jF,kF,ord,ex)] if(i2_lo<=i2_hi&&j2_lo<=j2_hi&&k2_lo<=k2_hi){for(int k0=k2_lo;k0<=k2_hi;++k0){const int kF=k0+1; diff --git a/AMSS_NCKU_source/kodiss_c.C b/AMSS_NCKU_source/kodiss_c.C index 0176088..09b3b3c 100644 --- a/AMSS_NCKU_source/kodiss_c.C +++ b/AMSS_NCKU_source/kodiss_c.C @@ -197,11 +197,11 @@ void kodis(const int ex[3], if (!(i0_lo > i0_hi || j0_lo > j0_hi || k0_lo > k0_hi)) { for (int k0 = k0_lo; k0 <= k0_hi; ++k0) { - const int kF = k0 + 3; + const int kF = k0 + 1; for (int j0 = j0_lo; j0 <= j0_hi; ++j0) { - const int jF = j0 + 3; + const int jF = j0 + 1; for (int i0 = i0_lo; i0 <= i0_hi; ++i0) { - const int iF = i0 + 3; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); /* Stencil: [1,-8,28,-56,70,-56,28,-8,1] */ @@ -272,11 +272,11 @@ void kodis(const int ex[3], if (!(i0_lo > i0_hi || j0_lo > j0_hi || k0_lo > k0_hi)) { for (int k0 = k0_lo; k0 <= k0_hi; ++k0) { - const int kF = k0 + 4; + const int kF = k0 + 1; for (int j0 = j0_lo; j0 <= j0_hi; ++j0) { - const int jF = j0 + 4; + const int jF = j0 + 1; for (int i0 = i0_lo; i0 <= i0_hi; ++i0) { - const int iF = i0 + 4; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); /* Stencil: [1,-10,45,-120,210,-252,210,-120,45,-10,1] */ diff --git a/AMSS_NCKU_source/kodiss_sh_c.C b/AMSS_NCKU_source/kodiss_sh_c.C index 5dfa69d..4e23768 100644 --- a/AMSS_NCKU_source/kodiss_sh_c.C +++ b/AMSS_NCKU_source/kodiss_sh_c.C @@ -32,7 +32,7 @@ extern "C" void kodis_sh_(const int ex[3], double *fh=(double*)malloc(fh_size*sizeof(double));if(!fh)return; symmetry_stbd(ord,ex,f,fh,SoA); - const int i0_lo=(iminF+1>0)?iminF+1:0,j0_lo=(jminF+1>0)?jminF+1:0,k0_lo=0; + const int i0_lo=(iminF+1>0)?iminF+1:0,j0_lo=(jminF+1>0)?jminF+1:0,k0_lo=2; const int i0_hi=imaxF-3,j0_hi=jmaxF-3,k0_hi=kmaxF-3; if(!(i0_lo>i0_hi||j0_lo>j0_hi||k0_lo>k0_hi)){ for(int k0=k0_lo;k0<=k0_hi;++k0){const int kF=k0+1; @@ -60,7 +60,7 @@ extern "C" void kodis_sh_(const int ex[3], double *fh=(double*)malloc(fh_size*sizeof(double));if(!fh)return; symmetry_stbd(ord,ex,f,fh,SoA); - const int i0_lo=(iminF+2>0)?iminF+2:0,j0_lo=(jminF+2>0)?jminF+2:0,k0_lo=0; + const int i0_lo=(iminF+2>0)?iminF+2:0,j0_lo=(jminF+2>0)?jminF+2:0,k0_lo=3; const int i0_hi=imaxF-4,j0_hi=jmaxF-4,k0_hi=kmaxF-4; if(!(i0_lo>i0_hi||j0_lo>j0_hi||k0_lo>k0_hi)){ for(int k0=k0_lo;k0<=k0_hi;++k0){const int kF=k0+1; @@ -88,7 +88,7 @@ extern "C" void kodis_sh_(const int ex[3], double *fh=(double*)malloc(fh_size*sizeof(double));if(!fh)return; symmetry_stbd(ord,ex,f,fh,SoA); - const int i0_lo=(iminF+3>0)?iminF+3:0,j0_lo=(jminF+3>0)?jminF+3:0,k0_lo=0; + const int i0_lo=(iminF+3>0)?iminF+3:0,j0_lo=(jminF+3>0)?jminF+3:0,k0_lo=4; const int i0_hi=imaxF-5,j0_hi=jmaxF-5,k0_hi=kmaxF-5; if(!(i0_lo>i0_hi||j0_lo>j0_hi||k0_lo>k0_hi)){ for(int k0=k0_lo;k0<=k0_hi;++k0){const int kF=k0+1; @@ -116,7 +116,7 @@ extern "C" void kodis_sh_(const int ex[3], double *fh=(double*)malloc(fh_size*sizeof(double));if(!fh)return; symmetry_stbd(ord,ex,f,fh,SoA); - const int i0_lo=(iminF+4>0)?iminF+4:0,j0_lo=(jminF+4>0)?jminF+4:0,k0_lo=0; + const int i0_lo=(iminF+4>0)?iminF+4:0,j0_lo=(jminF+4>0)?jminF+4:0,k0_lo=5; const int i0_hi=imaxF-6,j0_hi=jmaxF-6,k0_hi=kmaxF-6; if(!(i0_lo>i0_hi||j0_lo>j0_hi||k0_lo>k0_hi)){ for(int k0=k0_lo;k0<=k0_hi;++k0){const int kF=k0+1; diff --git a/AMSS_NCKU_source/lopsided_c.C b/AMSS_NCKU_source/lopsided_c.C index 9df2450..2dd839e 100644 --- a/AMSS_NCKU_source/lopsided_c.C +++ b/AMSS_NCKU_source/lopsided_c.C @@ -311,11 +311,11 @@ void lopsided(const int ex[3], const int imaxF = ex1, jmaxF = ex2, kmaxF = ex3; for (int k0 = 0; k0 <= ex3 - 2; ++k0) { - const int kF = k0 + 3; + const int kF = k0 + 1; for (int j0 = 0; j0 <= ex2 - 2; ++j0) { - const int jF = j0 + 3; + const int jF = j0 + 1; for (int i0 = 0; i0 <= ex1 - 2; ++i0) { - const int iF = i0 + 3; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); /* ---- x-direction ---- */ @@ -466,11 +466,11 @@ void lopsided(const int ex[3], const int imaxF = ex1, jmaxF = ex2, kmaxF = ex3; for (int k0 = 0; k0 <= ex3 - 2; ++k0) { - const int kF = k0 + 4; + const int kF = k0 + 1; for (int j0 = 0; j0 <= ex2 - 2; ++j0) { - const int jF = j0 + 4; + const int jF = j0 + 1; for (int i0 = 0; i0 <= ex1 - 2; ++i0) { - const int iF = i0 + 4; + const int iF = i0 + 1; const size_t p = idx_ex(i0, j0, k0, ex); const double sfx = Sfx[p]; diff --git a/AMSS_NCKU_source/lopsided_kodis_c.C b/AMSS_NCKU_source/lopsided_kodis_c.C index fe66c3e..a8c9f83 100644 --- a/AMSS_NCKU_source/lopsided_kodis_c.C +++ b/AMSS_NCKU_source/lopsided_kodis_c.C @@ -218,9 +218,9 @@ void lopsided_kodis(const int ex[3], const double d2dx=ONE/TWO/dX, d2dy=ONE/TWO/dY, d2dz=ONE/TWO/dZ; /* ---- advection (6th-order lopsided) ---- */ - for (int k0=0;k0<=ex3-2;++k0) { const int kF=k0+3; - for (int j0=0;j0<=ex2-2;++j0) { const int jF=j0+3; - for (int i0=0;i0<=ex1-2;++i0) { const int iF=i0+3; + for (int k0=0;k0<=ex3-2;++k0) { const int kF=k0+1; + for (int j0=0;j0<=ex2-2;++j0) { const int jF=j0+1; + for (int i0=0;i0<=ex1-2;++i0) { const int iF=i0+1; const size_t p=idx_ex(i0,j0,k0,ex); /* x */ const double sfx=Sfx[p]; @@ -276,9 +276,9 @@ void lopsided_kodis(const int ex[3], const int i0_lo=(iminF+3>0)?iminF+3:0, j0_lo=(jminF+3>0)?jminF+3:0, k0_lo=(kminF+3>0)?kminF+3:0; const int i0_hi=imaxF-5, j0_hi=jmaxF-5, k0_hi=kmaxF-5; if (!(i0_lo>i0_hi||j0_lo>j0_hi||k0_lo>k0_hi)) { - for (int k0=k0_lo;k0<=k0_hi;++k0) { const int kF=k0+3; - for (int j0=j0_lo;j0<=j0_hi;++j0) { const int jF=j0+3; - for (int i0=i0_lo;i0<=i0_hi;++i0) { const int iF=i0+3; + for (int k0=k0_lo;k0<=k0_hi;++k0) { const int kF=k0+1; + for (int j0=j0_lo;j0<=j0_hi;++j0) { const int jF=j0+1; + for (int i0=i0_lo;i0<=i0_hi;++i0) { const int iF=i0+1; const size_t p=idx_ex(i0,j0,k0,ex); const double Dx=((fh[idx_fh_F_ord4(iF-4,jF,kF,ex)]+fh[idx_fh_F_ord4(iF+4,jF,kF,ex)])-F8k*(fh[idx_fh_F_ord4(iF-3,jF,kF,ex)]+fh[idx_fh_F_ord4(iF+3,jF,kF,ex)])+F28*(fh[idx_fh_F_ord4(iF-2,jF,kF,ex)]+fh[idx_fh_F_ord4(iF+2,jF,kF,ex)])-F56*(fh[idx_fh_F_ord4(iF-1,jF,kF,ex)]+fh[idx_fh_F_ord4(iF+1,jF,kF,ex)])+F70*fh[idx_fh_F_ord4(iF,jF,kF,ex)])/dX; const double Dy=((fh[idx_fh_F_ord4(iF,jF-4,kF,ex)]+fh[idx_fh_F_ord4(iF,jF+4,kF,ex)])-F8k*(fh[idx_fh_F_ord4(iF,jF-3,kF,ex)]+fh[idx_fh_F_ord4(iF,jF+3,kF,ex)])+F28*(fh[idx_fh_F_ord4(iF,jF-2,kF,ex)]+fh[idx_fh_F_ord4(iF,jF+2,kF,ex)])-F56*(fh[idx_fh_F_ord4(iF,jF-1,kF,ex)]+fh[idx_fh_F_ord4(iF,jF+1,kF,ex)])+F70*fh[idx_fh_F_ord4(iF,jF,kF,ex)])/dY; @@ -311,9 +311,9 @@ void lopsided_kodis(const int ex[3], const double d2dx=ONE/TWO/dX, d2dy=ONE/TWO/dY, d2dz=ONE/TWO/dZ; /* ---- advection (8th-order lopsided) ---- */ - for (int k0=0;k0<=ex3-2;++k0) { const int kF=k0+4; - for (int j0=0;j0<=ex2-2;++j0) { const int jF=j0+4; - for (int i0=0;i0<=ex1-2;++i0) { const int iF=i0+4; + for (int k0=0;k0<=ex3-2;++k0) { const int kF=k0+1; + for (int j0=0;j0<=ex2-2;++j0) { const int jF=j0+1; + for (int i0=0;i0<=ex1-2;++i0) { const int iF=i0+1; const size_t p=idx_ex(i0,j0,k0,ex); const double sfx=Sfx[p]; if (sfx>ZEO) { @@ -366,9 +366,9 @@ void lopsided_kodis(const int ex[3], const int i0_lo=(iminF+4>0)?iminF+4:0, j0_lo=(jminF+4>0)?jminF+4:0, k0_lo=(kminF+4>0)?kminF+4:0; const int i0_hi=imaxF-6, j0_hi=jmaxF-6, k0_hi=kmaxF-6; if (!(i0_lo>i0_hi||j0_lo>j0_hi||k0_lo>k0_hi)) { - for (int k0=k0_lo;k0<=k0_hi;++k0) { const int kF=k0+4; - for (int j0=j0_lo;j0<=j0_hi;++j0) { const int jF=j0+4; - for (int i0=i0_lo;i0<=i0_hi;++i0) { const int iF=i0+4; + for (int k0=k0_lo;k0<=k0_hi;++k0) { const int kF=k0+1; + for (int j0=j0_lo;j0<=j0_hi;++j0) { const int jF=j0+1; + for (int i0=i0_lo;i0<=i0_hi;++i0) { const int iF=i0+1; const size_t p=idx_ex(i0,j0,k0,ex); const double Dx=((fh[idx_fh_F_ord5(iF-5,jF,kF,ex)]+fh[idx_fh_F_ord5(iF+5,jF,kF,ex)])-F10k*(fh[idx_fh_F_ord5(iF-4,jF,kF,ex)]+fh[idx_fh_F_ord5(iF+4,jF,kF,ex)])+F45k*(fh[idx_fh_F_ord5(iF-3,jF,kF,ex)]+fh[idx_fh_F_ord5(iF+3,jF,kF,ex)])-F120*(fh[idx_fh_F_ord5(iF-2,jF,kF,ex)]+fh[idx_fh_F_ord5(iF+2,jF,kF,ex)])+F210*(fh[idx_fh_F_ord5(iF-1,jF,kF,ex)]+fh[idx_fh_F_ord5(iF+1,jF,kF,ex)])-F252*fh[idx_fh_F_ord5(iF,jF,kF,ex)])/dX; const double Dy=((fh[idx_fh_F_ord5(iF,jF-5,kF,ex)]+fh[idx_fh_F_ord5(iF,jF+5,kF,ex)])-F10k*(fh[idx_fh_F_ord5(iF,jF-4,kF,ex)]+fh[idx_fh_F_ord5(iF,jF+4,kF,ex)])+F45k*(fh[idx_fh_F_ord5(iF,jF-3,kF,ex)]+fh[idx_fh_F_ord5(iF,jF+3,kF,ex)])-F120*(fh[idx_fh_F_ord5(iF,jF-2,kF,ex)]+fh[idx_fh_F_ord5(iF,jF+2,kF,ex)])+F210*(fh[idx_fh_F_ord5(iF,jF-1,kF,ex)]+fh[idx_fh_F_ord5(iF,jF+1,kF,ex)])-F252*fh[idx_fh_F_ord5(iF,jF,kF,ex)])/dY; diff --git a/AMSS_NCKU_source/share_func.h b/AMSS_NCKU_source/share_func.h index 659d1fc..125a058 100644 --- a/AMSS_NCKU_source/share_func.h +++ b/AMSS_NCKU_source/share_func.h @@ -317,7 +317,7 @@ static inline void symmetry_stbd(int ord, const double *src = func + (size_t)k0 * (size_t)extc2 * (size_t)extc1; const size_t kbase = (size_t)k0 * splane; for (int j0 = 0; j0 < extc2; ++j0) { - double *dst = funcc + kbase + (size_t)(sh + j0) * snx + (size_t)sh; + double *dst = funcc + kbase + (size_t)(sh + j0 + 1) * snx + (size_t)(sh + 1); const double *s = src + (size_t)j0 * (size_t)extc1; for (int i0 = 0; i0 < extc1; ++i0) dst[i0] = s[i0]; } @@ -328,7 +328,7 @@ static inline void symmetry_stbd(int ord, for (int k0 = 0; k0 < extc3; ++k0) { const size_t kbase = (size_t)k0 * splane; for (int j0 = 0; j0 < extc2; ++j0) { - const size_t off = kbase + (size_t)(sh + j0) * snx; + const size_t off = kbase + (size_t)(sh + j0 + 1) * snx; /* left side: funcc(-i) = funcc(i+2) * s1 */ for (int i = 0; i < ord; ++i) { funcc[off + (size_t)(sh - i)] = funcc[off + (size_t)(sh + i + 2)] * s1;