64-BitBrainstorm_2026/AMSS-NCKU
8
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Fix C derivative ghost-buffer indexing across FD orders

Browse Source
This commit is contained in:
CGH0S7
2026-05-14 16:00:11 +08:00
parent 7e1ee37bf7
commit b489564311
9 changed files with 140 additions and 140 deletions
Download Patch File Download Diff File Expand all files Collapse all files

90
AMSS_NCKU_source/fdderivs_c.C
View File

@@ -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) */

20
AMSS_NCKU_source/fdderivs_sh_c.C
View File

@@ -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;p<all;++p){fxx[p]=fyy[p]=fzz[p]=ZEO;fxy[p]=fxz[p]=fyz[p]=ZEO;}
const int i2_lo=(iminF>0)?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;p<all;++p){fxx[p]=fyy[p]=fzz[p]=fxy[p]=fxz[p]=fyz[p]=ZEO;}
const int i2_lo=(iminF>0)?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;p<all;++p){fxx[p]=fyy[p]=fzz[p]=fxy[p]=fxz[p]=fyz[p]=ZEO;}
const int i2_lo=(iminF+1>0)?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;p<all;++p){fxx[p]=fyy[p]=fzz[p]=fxy[p]=fxz[p]=fyz[p]=ZEO;}
const int i2_lo=(iminF+2>0)?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)]

90
AMSS_NCKU_source/fderivs_c.C
View File

@@ -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 * (

20
AMSS_NCKU_source/fderivs_sh_c.C
View File

@@ -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;p<all;++p) { fx[p]=ZEO; fy[p]=ZEO; fz[p]=ZEO; }
const int i2_lo=(iminF>0)?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;p<all;++p){fx[p]=ZEO;fy[p]=ZEO;fz[p]=ZEO;}
const int i2_lo=(iminF>0)?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;p<all;++p){fx[p]=ZEO;fy[p]=ZEO;fz[p]=ZEO;}
const int i2_lo=(iminF+1>0)?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;p<all;++p){fx[p]=ZEO;fy[p]=ZEO;fz[p]=ZEO;}
const int i2_lo=(iminF+2>0)?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;

12
AMSS_NCKU_source/kodiss_c.C
View File

@@ -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] */

8
AMSS_NCKU_source/kodiss_sh_c.C
View File

@@ -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;

12
AMSS_NCKU_source/lopsided_c.C
View File

@@ -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];

24
AMSS_NCKU_source/lopsided_kodis_c.C
View File

@@ -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;

4
AMSS_NCKU_source/share_func.h
View File

@@ -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;