diff --git a/src/gpu/VirtualFluids_GPU/GPU/GPU_Interface.h b/src/gpu/VirtualFluids_GPU/GPU/GPU_Interface.h
index c2ad940d0c2323ef808dd9b15061e7b4d6a7b722..80f683ea097600da4e29ce7ba4a62298ca9c073e 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/GPU_Interface.h
+++ b/src/gpu/VirtualFluids_GPU/GPU/GPU_Interface.h
@@ -8,6 +8,7 @@
//random numbers
#include <curand.h>
#include <curand_kernel.h>
+#include <cuda_runtime.h>
#include <DataTypes.h>
#include "LBM/LB.h"
@@ -2093,6 +2094,32 @@ extern "C" void ScaleFC_RhoSq_comp_27( real* DC,
unsigned int numberOfThreads,
OffFC offFC);
+extern "C" void ScaleFC_RhoSq_comp_27_Stream( real* DC,
+ real* DF,
+ unsigned int* neighborCX,
+ unsigned int* neighborCY,
+ unsigned int* neighborCZ,
+ unsigned int* neighborFX,
+ unsigned int* neighborFY,
+ unsigned int* neighborFZ,
+ unsigned int size_MatC,
+ unsigned int size_MatF,
+ bool evenOrOdd,
+ unsigned int* posC,
+ unsigned int* posFSWB,
+ unsigned int kFC,
+ real omCoarse,
+ real omFine,
+ real nu,
+ unsigned int nxC,
+ unsigned int nyC,
+ unsigned int nxF,
+ unsigned int nyF,
+ OffFC offFC,
+ unsigned int *fluidNodeIndices,
+ unsigned int numberOfFluidNodes,
+ CUstream_st *stream);
+
extern "C" void ScaleFC_RhoSq_3rdMom_comp_27( real* DC,
real* DF,
unsigned int* neighborCX,
diff --git a/src/gpu/VirtualFluids_GPU/GPU/GPU_Kernels.cuh b/src/gpu/VirtualFluids_GPU/GPU/GPU_Kernels.cuh
index 288db43e7bcd36dc4d187982b86178d345601094..fc35b1c661515f08c231be28f7041cc91a849373 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/GPU_Kernels.cuh
+++ b/src/gpu/VirtualFluids_GPU/GPU/GPU_Kernels.cuh
@@ -2049,6 +2049,14 @@ extern "C" __global__ void scaleFC_RhoSq_comp_27( real* DC,
unsigned int nyF,
OffFC offFC);
+extern "C" __global__ void
+scaleFC_RhoSq_comp_27_Stream(real *DC, real *DF, unsigned int *neighborCX, unsigned int *neighborCY,
+ unsigned int *neighborCZ, unsigned int *neighborFX, unsigned int *neighborFY,
+ unsigned int *neighborFZ, unsigned int size_MatC, unsigned int size_MatF, bool evenOrOdd,
+ unsigned int *posC, unsigned int *posFSWB, unsigned int kFC, real omCoarse, real omFine,
+ real nu, unsigned int nxC, unsigned int nyC, unsigned int nxF, unsigned int nyF,
+ OffFC offFC, const unsigned int *fluidNodeIndices, unsigned int numberOfFluidNodes);
+
extern "C" __global__ void scaleFC_RhoSq_3rdMom_comp_27(real* DC,
real* DF,
unsigned int* neighborCX,
diff --git a/src/gpu/VirtualFluids_GPU/GPU/LBMKernel.cu b/src/gpu/VirtualFluids_GPU/GPU/LBMKernel.cu
index 27a27bbdfc6e5649e07ec09f552cb29dcc2b1b34..b2a848e98b1b9145d282dde29cdb7058180d9338 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/LBMKernel.cu
+++ b/src/gpu/VirtualFluids_GPU/GPU/LBMKernel.cu
@@ -6242,6 +6242,48 @@ extern "C" void ScaleFC_RhoSq_comp_27(real* DC,
offFC);
getLastCudaError("scaleFC_RhoSq_27 execution failed");
}
+extern "C" void ScaleFC_RhoSq_comp_27_Stream(real *DC,
+ real * DF,
+ unsigned int * neighborCX,
+ unsigned int * neighborCY,
+ unsigned int * neighborCZ,
+ unsigned int * neighborFX,
+ unsigned int * neighborFY,
+ unsigned int * neighborFZ,
+ unsigned int size_MatC,
+ unsigned int size_MatF,
+ bool evenOrOdd,
+ unsigned int * posC,
+ unsigned int * posFSWB,
+ unsigned int kFC,
+ real omCoarse,
+ real omFine,
+ real nu,
+ unsigned int nxC,
+ unsigned int nyC,
+ unsigned int nxF,
+ unsigned int nyF,
+ OffFC offFC,
+ unsigned int *fluidNodeIndices,
+ unsigned int numberOfFluidNodes,
+ CUstream_st* stream)
+{
+ int Grid = (kFC / numberOfFluidNodes) + 1;
+ int Grid1, Grid2;
+ if (Grid > 512) {
+ Grid1 = 512;
+ Grid2 = (Grid / Grid1) + 1;
+ } else {
+ Grid1 = 1;
+ Grid2 = Grid;
+ }
+ dim3 gridINT_FC(Grid1, Grid2);
+ dim3 threads(numberOfFluidNodes, 1, 1);
+
+ scaleFC_RhoSq_comp_27_Stream<<<gridINT_FC, threads, 0, stream>>>(DC, DF, neighborCX, neighborCY, neighborCZ, neighborFX, neighborFY, neighborFZ, size_MatC, size_MatF, evenOrOdd,
+ posC, posFSWB, kFC, omCoarse, omFine, nu, nxC, nyC, nxF, nyF, offFC, fluidNodeIndices, numberOfFluidNodes);
+ getLastCudaError("scaleFC_RhoSq_27_Stream execution failed");
+}
//////////////////////////////////////////////////////////////////////////
extern "C" void ScaleFC_RhoSq_3rdMom_comp_27( real* DC,
real* DF,
diff --git a/src/gpu/VirtualFluids_GPU/GPU/ScaleFC27.cu b/src/gpu/VirtualFluids_GPU/GPU/ScaleFC27.cu
index 06d5086cc7f57c8245517f922662e8eea4571d05..a56c75e5859d8a722e148d07643763fa31086718 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/ScaleFC27.cu
+++ b/src/gpu/VirtualFluids_GPU/GPU/ScaleFC27.cu
@@ -9586,6 +9586,1475 @@ extern "C" __global__ void scaleFC_RhoSq_3rdMom_comp_27(real* DC,
//////////////////////////////////////////////////////////////////////////
+__device__ void scaleFC_RhoSq_comp_27_Calculation(real *DC, real *DF, unsigned int *neighborCX, unsigned int *neighborCY,
+ unsigned int *neighborCZ, unsigned int *neighborFX, unsigned int *neighborFY,
+ unsigned int *neighborFZ, unsigned int size_MatC, unsigned int size_MatF,
+ bool evenOrOdd, unsigned int *posC, unsigned int *posFSWB, unsigned int kFC,
+ real omCoarse, real omFine, real nu, unsigned int nxC, unsigned int nyC,
+ unsigned int nxF, unsigned int nyF, OffFC offFC, const unsigned k)
+{
+ real *feF, *fwF, *fnF, *fsF, *ftF, *fbF, *fneF, *fswF, *fseF, *fnwF, *fteF, *fbwF, *fbeF, *ftwF, *ftnF, *fbsF,
+ *fbnF, *ftsF, *fzeroF, *ftneF, *ftswF, *ftseF, *ftnwF, *fbneF, *fbswF, *fbseF, *fbnwF;
+
+ feF = &DF[dirE * size_MatF];
+ fwF = &DF[dirW * size_MatF];
+ fnF = &DF[dirN * size_MatF];
+ fsF = &DF[dirS * size_MatF];
+ ftF = &DF[dirT * size_MatF];
+ fbF = &DF[dirB * size_MatF];
+ fneF = &DF[dirNE * size_MatF];
+ fswF = &DF[dirSW * size_MatF];
+ fseF = &DF[dirSE * size_MatF];
+ fnwF = &DF[dirNW * size_MatF];
+ fteF = &DF[dirTE * size_MatF];
+ fbwF = &DF[dirBW * size_MatF];
+ fbeF = &DF[dirBE * size_MatF];
+ ftwF = &DF[dirTW * size_MatF];
+ ftnF = &DF[dirTN * size_MatF];
+ fbsF = &DF[dirBS * size_MatF];
+ fbnF = &DF[dirBN * size_MatF];
+ ftsF = &DF[dirTS * size_MatF];
+ fzeroF = &DF[dirZERO * size_MatF];
+ ftneF = &DF[dirTNE * size_MatF];
+ ftswF = &DF[dirTSW * size_MatF];
+ ftseF = &DF[dirTSE * size_MatF];
+ ftnwF = &DF[dirTNW * size_MatF];
+ fbneF = &DF[dirBNE * size_MatF];
+ fbswF = &DF[dirBSW * size_MatF];
+ fbseF = &DF[dirBSE * size_MatF];
+ fbnwF = &DF[dirBNW * size_MatF];
+
+ real *feC, *fwC, *fnC, *fsC, *ftC, *fbC, *fneC, *fswC, *fseC, *fnwC, *fteC, *fbwC, *fbeC, *ftwC, *ftnC, *fbsC,
+ *fbnC, *ftsC, *fzeroC, *ftneC, *ftswC, *ftseC, *ftnwC, *fbneC, *fbswC, *fbseC, *fbnwC;
+
+ if (evenOrOdd == true) {
+ feC = &DC[dirE * size_MatC];
+ fwC = &DC[dirW * size_MatC];
+ fnC = &DC[dirN * size_MatC];
+ fsC = &DC[dirS * size_MatC];
+ ftC = &DC[dirT * size_MatC];
+ fbC = &DC[dirB * size_MatC];
+ fneC = &DC[dirNE * size_MatC];
+ fswC = &DC[dirSW * size_MatC];
+ fseC = &DC[dirSE * size_MatC];
+ fnwC = &DC[dirNW * size_MatC];
+ fteC = &DC[dirTE * size_MatC];
+ fbwC = &DC[dirBW * size_MatC];
+ fbeC = &DC[dirBE * size_MatC];
+ ftwC = &DC[dirTW * size_MatC];
+ ftnC = &DC[dirTN * size_MatC];
+ fbsC = &DC[dirBS * size_MatC];
+ fbnC = &DC[dirBN * size_MatC];
+ ftsC = &DC[dirTS * size_MatC];
+ fzeroC = &DC[dirZERO * size_MatC];
+ ftneC = &DC[dirTNE * size_MatC];
+ ftswC = &DC[dirTSW * size_MatC];
+ ftseC = &DC[dirTSE * size_MatC];
+ ftnwC = &DC[dirTNW * size_MatC];
+ fbneC = &DC[dirBNE * size_MatC];
+ fbswC = &DC[dirBSW * size_MatC];
+ fbseC = &DC[dirBSE * size_MatC];
+ fbnwC = &DC[dirBNW * size_MatC];
+ } else {
+ fwC = &DC[dirE * size_MatC];
+ feC = &DC[dirW * size_MatC];
+ fsC = &DC[dirN * size_MatC];
+ fnC = &DC[dirS * size_MatC];
+ fbC = &DC[dirT * size_MatC];
+ ftC = &DC[dirB * size_MatC];
+ fswC = &DC[dirNE * size_MatC];
+ fneC = &DC[dirSW * size_MatC];
+ fnwC = &DC[dirSE * size_MatC];
+ fseC = &DC[dirNW * size_MatC];
+ fbwC = &DC[dirTE * size_MatC];
+ fteC = &DC[dirBW * size_MatC];
+ ftwC = &DC[dirBE * size_MatC];
+ fbeC = &DC[dirTW * size_MatC];
+ fbsC = &DC[dirTN * size_MatC];
+ ftnC = &DC[dirBS * size_MatC];
+ ftsC = &DC[dirBN * size_MatC];
+ fbnC = &DC[dirTS * size_MatC];
+ fzeroC = &DC[dirZERO * size_MatC];
+ fbswC = &DC[dirTNE * size_MatC];
+ fbneC = &DC[dirTSW * size_MatC];
+ fbnwC = &DC[dirTSE * size_MatC];
+ fbseC = &DC[dirTNW * size_MatC];
+ ftswC = &DC[dirBNE * size_MatC];
+ ftneC = &DC[dirBSW * size_MatC];
+ ftnwC = &DC[dirBSE * size_MatC];
+ ftseC = &DC[dirBNW * size_MatC];
+ }
+
+ ////////////////////////////////////////////////////////////////////////////////
+ real eps_new = c2o1;
+ real omegaS = omFine; //-omFine;
+ real o = omCoarse; //-omCoarse;
+ // real op = one;
+ // real cu_sq;
+
+ real xoff, yoff, zoff;
+ real xoff_sq, yoff_sq, zoff_sq;
+
+ real press; //,drho,vx1,vx2,vx3;
+ real /*press_SWT,*/ drho_SWT, vx1_SWT, vx2_SWT, vx3_SWT;
+ real /*press_NWT,*/ drho_NWT, vx1_NWT, vx2_NWT, vx3_NWT;
+ real /*press_NET,*/ drho_NET, vx1_NET, vx2_NET, vx3_NET;
+ real /*press_SET,*/ drho_SET, vx1_SET, vx2_SET, vx3_SET;
+ real /*press_SWB,*/ drho_SWB, vx1_SWB, vx2_SWB, vx3_SWB;
+ real /*press_NWB,*/ drho_NWB, vx1_NWB, vx2_NWB, vx3_NWB;
+ real /*press_NEB,*/ drho_NEB, vx1_NEB, vx2_NEB, vx3_NEB;
+ real /*press_SEB,*/ drho_SEB, vx1_SEB, vx2_SEB, vx3_SEB;
+ real f_E, f_W, f_N, f_S, f_T, f_B, f_NE, f_SW, f_SE, f_NW, f_TE, f_BW, f_BE, f_TW, f_TN, f_BS, f_BN, f_TS, f_ZERO,
+ f_TNE, f_TSW, f_TSE, f_TNW, f_BNE, f_BSW, f_BSE, f_BNW;
+ // real
+ // feq_E,feq_W,feq_N,feq_S,feq_T,feq_B,feq_NE,feq_SW,feq_SE,feq_NW,feq_TE,feq_BW,feq_BE,feq_TW,feq_TN,feq_BS,feq_BN,feq_TS,feq_ZERO,feq_TNE,
+ // feq_TSW, feq_TSE, feq_TNW, feq_BNE, feq_BSW, feq_BSE, feq_BNW;
+ real kxyFromfcNEQ_SWT, kyzFromfcNEQ_SWT, kxzFromfcNEQ_SWT, kxxMyyFromfcNEQ_SWT, kxxMzzFromfcNEQ_SWT;
+ real kxyFromfcNEQ_NWT, kyzFromfcNEQ_NWT, kxzFromfcNEQ_NWT, kxxMyyFromfcNEQ_NWT, kxxMzzFromfcNEQ_NWT;
+ real kxyFromfcNEQ_NET, kyzFromfcNEQ_NET, kxzFromfcNEQ_NET, kxxMyyFromfcNEQ_NET, kxxMzzFromfcNEQ_NET;
+ real kxyFromfcNEQ_SET, kyzFromfcNEQ_SET, kxzFromfcNEQ_SET, kxxMyyFromfcNEQ_SET, kxxMzzFromfcNEQ_SET;
+ real kxyFromfcNEQ_SWB, kyzFromfcNEQ_SWB, kxzFromfcNEQ_SWB, kxxMyyFromfcNEQ_SWB, kxxMzzFromfcNEQ_SWB;
+ real kxyFromfcNEQ_NWB, kyzFromfcNEQ_NWB, kxzFromfcNEQ_NWB, kxxMyyFromfcNEQ_NWB, kxxMzzFromfcNEQ_NWB;
+ real kxyFromfcNEQ_NEB, kyzFromfcNEQ_NEB, kxzFromfcNEQ_NEB, kxxMyyFromfcNEQ_NEB, kxxMzzFromfcNEQ_NEB;
+ real kxyFromfcNEQ_SEB, kyzFromfcNEQ_SEB, kxzFromfcNEQ_SEB, kxxMyyFromfcNEQ_SEB, kxxMzzFromfcNEQ_SEB;
+ real a0, ax, ay, az, axx, ayy, azz, axy, axz, ayz, b0, bx, by, bz, bxx, byy, bzz, bxy, bxz, byz, c0, cx, cy, cz,
+ cxx, cyy, czz, cxy, cxz, cyz /*, axyz, bxyz, cxyz*/;
+ real d0, dx, dy, dz, dxy, dxz, dyz /*, dxyz*/;
+
+ if (k < kFC) {
+ //////////////////////////////////////////////////////////////////////////
+ xoff = offFC.xOffFC[k];
+ yoff = offFC.yOffFC[k];
+ zoff = offFC.zOffFC[k];
+ xoff_sq = xoff * xoff;
+ yoff_sq = yoff * yoff;
+ zoff_sq = zoff * zoff;
+ //////////////////////////////////////////////////////////////////////////
+ // SWB//
+ //////////////////////////////////////////////////////////////////////////
+ // index 0
+ unsigned int k0zero = posFSWB[k];
+ unsigned int k0w = neighborFX[k0zero];
+ unsigned int k0s = neighborFY[k0zero];
+ unsigned int k0b = neighborFZ[k0zero];
+ unsigned int k0sw = neighborFY[k0w];
+ unsigned int k0bw = neighborFZ[k0w];
+ unsigned int k0bs = neighborFZ[k0s];
+ unsigned int k0bsw = neighborFZ[k0sw];
+ //////////////////////////////////////////////////////////////////////////
+ // index
+ unsigned int kzero = k0zero;
+ unsigned int kw = k0w;
+ unsigned int ks = k0s;
+ unsigned int kb = k0b;
+ unsigned int ksw = k0sw;
+ unsigned int kbw = k0bw;
+ unsigned int kbs = k0bs;
+ unsigned int kbsw = k0bsw;
+ ////////////////////////////////////////////////////////////////////////////////
+ f_E = feF[kzero];
+ f_W = fwF[kw];
+ f_N = fnF[kzero];
+ f_S = fsF[ks];
+ f_T = ftF[kzero];
+ f_B = fbF[kb];
+ f_NE = fneF[kzero];
+ f_SW = fswF[ksw];
+ f_SE = fseF[ks];
+ f_NW = fnwF[kw];
+ f_TE = fteF[kzero];
+ f_BW = fbwF[kbw];
+ f_BE = fbeF[kb];
+ f_TW = ftwF[kw];
+ f_TN = ftnF[kzero];
+ f_BS = fbsF[kbs];
+ f_BN = fbnF[kb];
+ f_TS = ftsF[ks];
+ f_ZERO = fzeroF[kzero];
+ f_TNE = ftneF[kzero];
+ f_TSW = ftswF[ksw];
+ f_TSE = ftseF[ks];
+ f_TNW = ftnwF[kw];
+ f_BNE = fbneF[kb];
+ f_BSW = fbswF[kbsw];
+ f_BSE = fbseF[kbs];
+ f_BNW = fbnwF[kbw];
+
+ drho_SWB = f_E + f_W + f_N + f_S + f_T + f_B + f_NE + f_SW + f_SE + f_NW + f_TE + f_BW + f_BE + f_TW + f_TN +
+ f_BS + f_BN + f_TS + f_ZERO + f_TNE + f_TSW + f_TSE + f_TNW + f_BNE + f_BSW + f_BSE + f_BNW;
+ vx1_SWB = (((f_TNE - f_BSW) + (f_TSE - f_BNW) + (f_BNE - f_TSW) + (f_BSE - f_TNW)) +
+ (((f_NE - f_SW) + (f_TE - f_BW)) + ((f_SE - f_NW) + (f_BE - f_TW))) + (f_E - f_W)) /
+ (c1o1 + drho_SWB);
+ vx2_SWB = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_BNE - f_TSW) + (f_BNW - f_TSE)) +
+ (((f_NE - f_SW) + (f_TN - f_BS)) + ((f_BN - f_TS) + (f_NW - f_SE))) + (f_N - f_S)) /
+ (c1o1 + drho_SWB);
+ vx3_SWB = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_TSE - f_BNW) + (f_TSW - f_BNE)) +
+ (((f_TE - f_BW) + (f_TN - f_BS)) + ((f_TW - f_BE) + (f_TS - f_BN))) + (f_T - f_B)) /
+ (c1o1 + drho_SWB);
+
+ kxyFromfcNEQ_SWB =
+ -c3o1 * omegaS *
+ ((f_SW + f_BSW + f_TSW - f_NW - f_BNW - f_TNW - f_SE - f_BSE - f_TSE + f_NE + f_BNE + f_TNE) /
+ (c1o1 + drho_SWB) -
+ ((vx1_SWB * vx2_SWB)));
+ kyzFromfcNEQ_SWB =
+ -c3o1 * omegaS *
+ ((f_BS + f_BSE + f_BSW - f_TS - f_TSE - f_TSW - f_BN - f_BNE - f_BNW + f_TN + f_TNE + f_TNW) /
+ (c1o1 + drho_SWB) -
+ ((vx2_SWB * vx3_SWB)));
+ kxzFromfcNEQ_SWB =
+ -c3o1 * omegaS *
+ ((f_BW + f_BSW + f_BNW - f_TW - f_TSW - f_TNW - f_BE - f_BSE - f_BNE + f_TE + f_TSE + f_TNE) /
+ (c1o1 + drho_SWB) -
+ ((vx1_SWB * vx3_SWB)));
+ kxxMyyFromfcNEQ_SWB =
+ -c3o2 * omegaS *
+ ((f_BW + f_W + f_TW - f_BS - f_S - f_TS - f_BN - f_N - f_TN + f_BE + f_E + f_TE) / (c1o1 + drho_SWB) -
+ ((vx1_SWB * vx1_SWB - vx2_SWB * vx2_SWB)));
+ kxxMzzFromfcNEQ_SWB =
+ -c3o2 * omegaS *
+ ((f_SW + f_W + f_NW - f_BS - f_TS - f_B - f_T - f_BN - f_TN + f_SE + f_E + f_NE) / (c1o1 + drho_SWB) -
+ ((vx1_SWB * vx1_SWB - vx3_SWB * vx3_SWB)));
+
+ //////////////////////////////////////////////////////////////////////////
+ // SWT//
+ //////////////////////////////////////////////////////////////////////////
+ // index
+ kzero = kb;
+ kw = kbw;
+ ks = kbs;
+ kb = neighborFZ[kb];
+ ksw = kbsw;
+ kbw = neighborFZ[kbw];
+ kbs = neighborFZ[kbs];
+ kbsw = neighborFZ[kbsw];
+ ////////////////////////////////////////////////////////////////////////////////
+ f_E = feF[kzero];
+ f_W = fwF[kw];
+ f_N = fnF[kzero];
+ f_S = fsF[ks];
+ f_T = ftF[kzero];
+ f_B = fbF[kb];
+ f_NE = fneF[kzero];
+ f_SW = fswF[ksw];
+ f_SE = fseF[ks];
+ f_NW = fnwF[kw];
+ f_TE = fteF[kzero];
+ f_BW = fbwF[kbw];
+ f_BE = fbeF[kb];
+ f_TW = ftwF[kw];
+ f_TN = ftnF[kzero];
+ f_BS = fbsF[kbs];
+ f_BN = fbnF[kb];
+ f_TS = ftsF[ks];
+ f_ZERO = fzeroF[kzero];
+ f_TNE = ftneF[kzero];
+ f_TSW = ftswF[ksw];
+ f_TSE = ftseF[ks];
+ f_TNW = ftnwF[kw];
+ f_BNE = fbneF[kb];
+ f_BSW = fbswF[kbsw];
+ f_BSE = fbseF[kbs];
+ f_BNW = fbnwF[kbw];
+
+ drho_SWT = f_E + f_W + f_N + f_S + f_T + f_B + f_NE + f_SW + f_SE + f_NW + f_TE + f_BW + f_BE + f_TW + f_TN +
+ f_BS + f_BN + f_TS + f_ZERO + f_TNE + f_TSW + f_TSE + f_TNW + f_BNE + f_BSW + f_BSE + f_BNW;
+ vx1_SWT = (((f_TNE - f_BSW) + (f_TSE - f_BNW) + (f_BNE - f_TSW) + (f_BSE - f_TNW)) +
+ (((f_NE - f_SW) + (f_TE - f_BW)) + ((f_SE - f_NW) + (f_BE - f_TW))) + (f_E - f_W)) /
+ (c1o1 + drho_SWT);
+ vx2_SWT = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_BNE - f_TSW) + (f_BNW - f_TSE)) +
+ (((f_NE - f_SW) + (f_TN - f_BS)) + ((f_BN - f_TS) + (f_NW - f_SE))) + (f_N - f_S)) /
+ (c1o1 + drho_SWT);
+ vx3_SWT = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_TSE - f_BNW) + (f_TSW - f_BNE)) +
+ (((f_TE - f_BW) + (f_TN - f_BS)) + ((f_TW - f_BE) + (f_TS - f_BN))) + (f_T - f_B)) /
+ (c1o1 + drho_SWT);
+
+ kxyFromfcNEQ_SWT =
+ -c3o1 * omegaS *
+ ((f_SW + f_BSW + f_TSW - f_NW - f_BNW - f_TNW - f_SE - f_BSE - f_TSE + f_NE + f_BNE + f_TNE) /
+ (c1o1 + drho_SWT) -
+ ((vx1_SWT * vx2_SWT)));
+ kyzFromfcNEQ_SWT =
+ -c3o1 * omegaS *
+ ((f_BS + f_BSE + f_BSW - f_TS - f_TSE - f_TSW - f_BN - f_BNE - f_BNW + f_TN + f_TNE + f_TNW) /
+ (c1o1 + drho_SWT) -
+ ((vx2_SWT * vx3_SWT)));
+ kxzFromfcNEQ_SWT =
+ -c3o1 * omegaS *
+ ((f_BW + f_BSW + f_BNW - f_TW - f_TSW - f_TNW - f_BE - f_BSE - f_BNE + f_TE + f_TSE + f_TNE) /
+ (c1o1 + drho_SWT) -
+ ((vx1_SWT * vx3_SWT)));
+ kxxMyyFromfcNEQ_SWT =
+ -c3o2 * omegaS *
+ ((f_BW + f_W + f_TW - f_BS - f_S - f_TS - f_BN - f_N - f_TN + f_BE + f_E + f_TE) / (c1o1 + drho_SWT) -
+ ((vx1_SWT * vx1_SWT - vx2_SWT * vx2_SWT)));
+ kxxMzzFromfcNEQ_SWT =
+ -c3o2 * omegaS *
+ ((f_SW + f_W + f_NW - f_BS - f_TS - f_B - f_T - f_BN - f_TN + f_SE + f_E + f_NE) / (c1o1 + drho_SWT) -
+ ((vx1_SWT * vx1_SWT - vx3_SWT * vx3_SWT)));
+
+ //////////////////////////////////////////////////////////////////////////
+ // SET//
+ //////////////////////////////////////////////////////////////////////////
+ // index
+ kzero = kw;
+ kw = neighborFX[kw];
+ ks = ksw;
+ kb = kbw;
+ ksw = neighborFX[ksw];
+ kbw = neighborFX[kbw];
+ kbs = kbsw;
+ kbsw = neighborFX[kbsw];
+ ////////////////////////////////////////////////////////////////////////////////
+ f_E = feF[kzero];
+ f_W = fwF[kw];
+ f_N = fnF[kzero];
+ f_S = fsF[ks];
+ f_T = ftF[kzero];
+ f_B = fbF[kb];
+ f_NE = fneF[kzero];
+ f_SW = fswF[ksw];
+ f_SE = fseF[ks];
+ f_NW = fnwF[kw];
+ f_TE = fteF[kzero];
+ f_BW = fbwF[kbw];
+ f_BE = fbeF[kb];
+ f_TW = ftwF[kw];
+ f_TN = ftnF[kzero];
+ f_BS = fbsF[kbs];
+ f_BN = fbnF[kb];
+ f_TS = ftsF[ks];
+ f_ZERO = fzeroF[kzero];
+ f_TNE = ftneF[kzero];
+ f_TSW = ftswF[ksw];
+ f_TSE = ftseF[ks];
+ f_TNW = ftnwF[kw];
+ f_BNE = fbneF[kb];
+ f_BSW = fbswF[kbsw];
+ f_BSE = fbseF[kbs];
+ f_BNW = fbnwF[kbw];
+
+ drho_SET = f_E + f_W + f_N + f_S + f_T + f_B + f_NE + f_SW + f_SE + f_NW + f_TE + f_BW + f_BE + f_TW + f_TN +
+ f_BS + f_BN + f_TS + f_ZERO + f_TNE + f_TSW + f_TSE + f_TNW + f_BNE + f_BSW + f_BSE + f_BNW;
+ vx1_SET = (((f_TNE - f_BSW) + (f_TSE - f_BNW) + (f_BNE - f_TSW) + (f_BSE - f_TNW)) +
+ (((f_NE - f_SW) + (f_TE - f_BW)) + ((f_SE - f_NW) + (f_BE - f_TW))) + (f_E - f_W)) /
+ (c1o1 + drho_SET);
+ vx2_SET = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_BNE - f_TSW) + (f_BNW - f_TSE)) +
+ (((f_NE - f_SW) + (f_TN - f_BS)) + ((f_BN - f_TS) + (f_NW - f_SE))) + (f_N - f_S)) /
+ (c1o1 + drho_SET);
+ vx3_SET = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_TSE - f_BNW) + (f_TSW - f_BNE)) +
+ (((f_TE - f_BW) + (f_TN - f_BS)) + ((f_TW - f_BE) + (f_TS - f_BN))) + (f_T - f_B)) /
+ (c1o1 + drho_SET);
+
+ kxyFromfcNEQ_SET =
+ -c3o1 * omegaS *
+ ((f_SW + f_BSW + f_TSW - f_NW - f_BNW - f_TNW - f_SE - f_BSE - f_TSE + f_NE + f_BNE + f_TNE) /
+ (c1o1 + drho_SET) -
+ ((vx1_SET * vx2_SET)));
+ kyzFromfcNEQ_SET =
+ -c3o1 * omegaS *
+ ((f_BS + f_BSE + f_BSW - f_TS - f_TSE - f_TSW - f_BN - f_BNE - f_BNW + f_TN + f_TNE + f_TNW) /
+ (c1o1 + drho_SET) -
+ ((vx2_SET * vx3_SET)));
+ kxzFromfcNEQ_SET =
+ -c3o1 * omegaS *
+ ((f_BW + f_BSW + f_BNW - f_TW - f_TSW - f_TNW - f_BE - f_BSE - f_BNE + f_TE + f_TSE + f_TNE) /
+ (c1o1 + drho_SET) -
+ ((vx1_SET * vx3_SET)));
+ kxxMyyFromfcNEQ_SET =
+ -c3o2 * omegaS *
+ ((f_BW + f_W + f_TW - f_BS - f_S - f_TS - f_BN - f_N - f_TN + f_BE + f_E + f_TE) / (c1o1 + drho_SET) -
+ ((vx1_SET * vx1_SET - vx2_SET * vx2_SET)));
+ kxxMzzFromfcNEQ_SET =
+ -c3o2 * omegaS *
+ ((f_SW + f_W + f_NW - f_BS - f_TS - f_B - f_T - f_BN - f_TN + f_SE + f_E + f_NE) / (c1o1 + drho_SET) -
+ ((vx1_SET * vx1_SET - vx3_SET * vx3_SET)));
+
+ //////////////////////////////////////////////////////////////////////////
+ // SEB//
+ //////////////////////////////////////////////////////////////////////////
+ // index
+ kb = kzero;
+ kbw = kw;
+ kbs = ks;
+ kbsw = ksw;
+ kzero = k0w;
+ kw = neighborFX[k0w];
+ ks = k0sw;
+ ksw = neighborFX[k0sw];
+ ////////////////////////////////////////////////////////////////////////////////
+ f_E = feF[kzero];
+ f_W = fwF[kw];
+ f_N = fnF[kzero];
+ f_S = fsF[ks];
+ f_T = ftF[kzero];
+ f_B = fbF[kb];
+ f_NE = fneF[kzero];
+ f_SW = fswF[ksw];
+ f_SE = fseF[ks];
+ f_NW = fnwF[kw];
+ f_TE = fteF[kzero];
+ f_BW = fbwF[kbw];
+ f_BE = fbeF[kb];
+ f_TW = ftwF[kw];
+ f_TN = ftnF[kzero];
+ f_BS = fbsF[kbs];
+ f_BN = fbnF[kb];
+ f_TS = ftsF[ks];
+ f_ZERO = fzeroF[kzero];
+ f_TNE = ftneF[kzero];
+ f_TSW = ftswF[ksw];
+ f_TSE = ftseF[ks];
+ f_TNW = ftnwF[kw];
+ f_BNE = fbneF[kb];
+ f_BSW = fbswF[kbsw];
+ f_BSE = fbseF[kbs];
+ f_BNW = fbnwF[kbw];
+
+ drho_SEB = f_E + f_W + f_N + f_S + f_T + f_B + f_NE + f_SW + f_SE + f_NW + f_TE + f_BW + f_BE + f_TW + f_TN +
+ f_BS + f_BN + f_TS + f_ZERO + f_TNE + f_TSW + f_TSE + f_TNW + f_BNE + f_BSW + f_BSE + f_BNW;
+ vx1_SEB = (((f_TNE - f_BSW) + (f_TSE - f_BNW) + (f_BNE - f_TSW) + (f_BSE - f_TNW)) +
+ (((f_NE - f_SW) + (f_TE - f_BW)) + ((f_SE - f_NW) + (f_BE - f_TW))) + (f_E - f_W)) /
+ (c1o1 + drho_SEB);
+ vx2_SEB = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_BNE - f_TSW) + (f_BNW - f_TSE)) +
+ (((f_NE - f_SW) + (f_TN - f_BS)) + ((f_BN - f_TS) + (f_NW - f_SE))) + (f_N - f_S)) /
+ (c1o1 + drho_SEB);
+ vx3_SEB = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_TSE - f_BNW) + (f_TSW - f_BNE)) +
+ (((f_TE - f_BW) + (f_TN - f_BS)) + ((f_TW - f_BE) + (f_TS - f_BN))) + (f_T - f_B)) /
+ (c1o1 + drho_SEB);
+
+ kxyFromfcNEQ_SEB =
+ -c3o1 * omegaS *
+ ((f_SW + f_BSW + f_TSW - f_NW - f_BNW - f_TNW - f_SE - f_BSE - f_TSE + f_NE + f_BNE + f_TNE) /
+ (c1o1 + drho_SEB) -
+ ((vx1_SEB * vx2_SEB)));
+ kyzFromfcNEQ_SEB =
+ -c3o1 * omegaS *
+ ((f_BS + f_BSE + f_BSW - f_TS - f_TSE - f_TSW - f_BN - f_BNE - f_BNW + f_TN + f_TNE + f_TNW) /
+ (c1o1 + drho_SEB) -
+ ((vx2_SEB * vx3_SEB)));
+ kxzFromfcNEQ_SEB =
+ -c3o1 * omegaS *
+ ((f_BW + f_BSW + f_BNW - f_TW - f_TSW - f_TNW - f_BE - f_BSE - f_BNE + f_TE + f_TSE + f_TNE) /
+ (c1o1 + drho_SEB) -
+ ((vx1_SEB * vx3_SEB)));
+ kxxMyyFromfcNEQ_SEB =
+ -c3o2 * omegaS *
+ ((f_BW + f_W + f_TW - f_BS - f_S - f_TS - f_BN - f_N - f_TN + f_BE + f_E + f_TE) / (c1o1 + drho_SEB) -
+ ((vx1_SEB * vx1_SEB - vx2_SEB * vx2_SEB)));
+ kxxMzzFromfcNEQ_SEB =
+ -c3o2 * omegaS *
+ ((f_SW + f_W + f_NW - f_BS - f_TS - f_B - f_T - f_BN - f_TN + f_SE + f_E + f_NE) / (c1o1 + drho_SEB) -
+ ((vx1_SEB * vx1_SEB - vx3_SEB * vx3_SEB)));
+
+ //////////////////////////////////////////////////////////////////////////
+ // NWB//
+ //////////////////////////////////////////////////////////////////////////
+ // index 0
+ k0zero = k0s;
+ k0w = k0sw;
+ k0s = neighborFY[k0s];
+ k0b = k0bs;
+ k0sw = neighborFY[k0sw];
+ k0bw = k0bsw;
+ k0bs = neighborFY[k0bs];
+ k0bsw = neighborFY[k0bsw];
+ //////////////////////////////////////////////////////////////////////////
+ // index
+ kzero = k0zero;
+ kw = k0w;
+ ks = k0s;
+ kb = k0b;
+ ksw = k0sw;
+ kbw = k0bw;
+ kbs = k0bs;
+ kbsw = k0bsw;
+ ////////////////////////////////////////////////////////////////////////////////
+ f_E = feF[kzero];
+ f_W = fwF[kw];
+ f_N = fnF[kzero];
+ f_S = fsF[ks];
+ f_T = ftF[kzero];
+ f_B = fbF[kb];
+ f_NE = fneF[kzero];
+ f_SW = fswF[ksw];
+ f_SE = fseF[ks];
+ f_NW = fnwF[kw];
+ f_TE = fteF[kzero];
+ f_BW = fbwF[kbw];
+ f_BE = fbeF[kb];
+ f_TW = ftwF[kw];
+ f_TN = ftnF[kzero];
+ f_BS = fbsF[kbs];
+ f_BN = fbnF[kb];
+ f_TS = ftsF[ks];
+ f_ZERO = fzeroF[kzero];
+ f_TNE = ftneF[kzero];
+ f_TSW = ftswF[ksw];
+ f_TSE = ftseF[ks];
+ f_TNW = ftnwF[kw];
+ f_BNE = fbneF[kb];
+ f_BSW = fbswF[kbsw];
+ f_BSE = fbseF[kbs];
+ f_BNW = fbnwF[kbw];
+
+ drho_NWB = f_E + f_W + f_N + f_S + f_T + f_B + f_NE + f_SW + f_SE + f_NW + f_TE + f_BW + f_BE + f_TW + f_TN +
+ f_BS + f_BN + f_TS + f_ZERO + f_TNE + f_TSW + f_TSE + f_TNW + f_BNE + f_BSW + f_BSE + f_BNW;
+ vx1_NWB = (((f_TNE - f_BSW) + (f_TSE - f_BNW) + (f_BNE - f_TSW) + (f_BSE - f_TNW)) +
+ (((f_NE - f_SW) + (f_TE - f_BW)) + ((f_SE - f_NW) + (f_BE - f_TW))) + (f_E - f_W)) /
+ (c1o1 + drho_NWB);
+ vx2_NWB = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_BNE - f_TSW) + (f_BNW - f_TSE)) +
+ (((f_NE - f_SW) + (f_TN - f_BS)) + ((f_BN - f_TS) + (f_NW - f_SE))) + (f_N - f_S)) /
+ (c1o1 + drho_NWB);
+ vx3_NWB = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_TSE - f_BNW) + (f_TSW - f_BNE)) +
+ (((f_TE - f_BW) + (f_TN - f_BS)) + ((f_TW - f_BE) + (f_TS - f_BN))) + (f_T - f_B)) /
+ (c1o1 + drho_NWB);
+
+ kxyFromfcNEQ_NWB =
+ -c3o1 * omegaS *
+ ((f_SW + f_BSW + f_TSW - f_NW - f_BNW - f_TNW - f_SE - f_BSE - f_TSE + f_NE + f_BNE + f_TNE) /
+ (c1o1 + drho_NWB) -
+ ((vx1_NWB * vx2_NWB)));
+ kyzFromfcNEQ_NWB =
+ -c3o1 * omegaS *
+ ((f_BS + f_BSE + f_BSW - f_TS - f_TSE - f_TSW - f_BN - f_BNE - f_BNW + f_TN + f_TNE + f_TNW) /
+ (c1o1 + drho_NWB) -
+ ((vx2_NWB * vx3_NWB)));
+ kxzFromfcNEQ_NWB =
+ -c3o1 * omegaS *
+ ((f_BW + f_BSW + f_BNW - f_TW - f_TSW - f_TNW - f_BE - f_BSE - f_BNE + f_TE + f_TSE + f_TNE) /
+ (c1o1 + drho_NWB) -
+ ((vx1_NWB * vx3_NWB)));
+ kxxMyyFromfcNEQ_NWB =
+ -c3o2 * omegaS *
+ ((f_BW + f_W + f_TW - f_BS - f_S - f_TS - f_BN - f_N - f_TN + f_BE + f_E + f_TE) / (c1o1 + drho_NWB) -
+ ((vx1_NWB * vx1_NWB - vx2_NWB * vx2_NWB)));
+ kxxMzzFromfcNEQ_NWB =
+ -c3o2 * omegaS *
+ ((f_SW + f_W + f_NW - f_BS - f_TS - f_B - f_T - f_BN - f_TN + f_SE + f_E + f_NE) / (c1o1 + drho_NWB) -
+ ((vx1_NWB * vx1_NWB - vx3_NWB * vx3_NWB)));
+
+ //////////////////////////////////////////////////////////////////////////
+ // NWT//
+ //////////////////////////////////////////////////////////////////////////
+ // index
+ kzero = kb;
+ kw = kbw;
+ ks = kbs;
+ kb = neighborFZ[kb];
+ ksw = kbsw;
+ kbw = neighborFZ[kbw];
+ kbs = neighborFZ[kbs];
+ kbsw = neighborFZ[kbsw];
+ ////////////////////////////////////////////////////////////////////////////////
+ f_E = feF[kzero];
+ f_W = fwF[kw];
+ f_N = fnF[kzero];
+ f_S = fsF[ks];
+ f_T = ftF[kzero];
+ f_B = fbF[kb];
+ f_NE = fneF[kzero];
+ f_SW = fswF[ksw];
+ f_SE = fseF[ks];
+ f_NW = fnwF[kw];
+ f_TE = fteF[kzero];
+ f_BW = fbwF[kbw];
+ f_BE = fbeF[kb];
+ f_TW = ftwF[kw];
+ f_TN = ftnF[kzero];
+ f_BS = fbsF[kbs];
+ f_BN = fbnF[kb];
+ f_TS = ftsF[ks];
+ f_ZERO = fzeroF[kzero];
+ f_TNE = ftneF[kzero];
+ f_TSW = ftswF[ksw];
+ f_TSE = ftseF[ks];
+ f_TNW = ftnwF[kw];
+ f_BNE = fbneF[kb];
+ f_BSW = fbswF[kbsw];
+ f_BSE = fbseF[kbs];
+ f_BNW = fbnwF[kbw];
+
+ drho_NWT = f_E + f_W + f_N + f_S + f_T + f_B + f_NE + f_SW + f_SE + f_NW + f_TE + f_BW + f_BE + f_TW + f_TN +
+ f_BS + f_BN + f_TS + f_ZERO + f_TNE + f_TSW + f_TSE + f_TNW + f_BNE + f_BSW + f_BSE + f_BNW;
+ vx1_NWT = (((f_TNE - f_BSW) + (f_TSE - f_BNW) + (f_BNE - f_TSW) + (f_BSE - f_TNW)) +
+ (((f_NE - f_SW) + (f_TE - f_BW)) + ((f_SE - f_NW) + (f_BE - f_TW))) + (f_E - f_W)) /
+ (c1o1 + drho_NWT);
+ vx2_NWT = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_BNE - f_TSW) + (f_BNW - f_TSE)) +
+ (((f_NE - f_SW) + (f_TN - f_BS)) + ((f_BN - f_TS) + (f_NW - f_SE))) + (f_N - f_S)) /
+ (c1o1 + drho_NWT);
+ vx3_NWT = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_TSE - f_BNW) + (f_TSW - f_BNE)) +
+ (((f_TE - f_BW) + (f_TN - f_BS)) + ((f_TW - f_BE) + (f_TS - f_BN))) + (f_T - f_B)) /
+ (c1o1 + drho_NWT);
+
+ kxyFromfcNEQ_NWT =
+ -c3o1 * omegaS *
+ ((f_SW + f_BSW + f_TSW - f_NW - f_BNW - f_TNW - f_SE - f_BSE - f_TSE + f_NE + f_BNE + f_TNE) /
+ (c1o1 + drho_NWT) -
+ ((vx1_NWT * vx2_NWT)));
+ kyzFromfcNEQ_NWT =
+ -c3o1 * omegaS *
+ ((f_BS + f_BSE + f_BSW - f_TS - f_TSE - f_TSW - f_BN - f_BNE - f_BNW + f_TN + f_TNE + f_TNW) /
+ (c1o1 + drho_NWT) -
+ ((vx2_NWT * vx3_NWT)));
+ kxzFromfcNEQ_NWT =
+ -c3o1 * omegaS *
+ ((f_BW + f_BSW + f_BNW - f_TW - f_TSW - f_TNW - f_BE - f_BSE - f_BNE + f_TE + f_TSE + f_TNE) /
+ (c1o1 + drho_NWT) -
+ ((vx1_NWT * vx3_NWT)));
+ kxxMyyFromfcNEQ_NWT =
+ -c3o2 * omegaS *
+ ((f_BW + f_W + f_TW - f_BS - f_S - f_TS - f_BN - f_N - f_TN + f_BE + f_E + f_TE) / (c1o1 + drho_NWT) -
+ ((vx1_NWT * vx1_NWT - vx2_NWT * vx2_NWT)));
+ kxxMzzFromfcNEQ_NWT =
+ -c3o2 * omegaS *
+ ((f_SW + f_W + f_NW - f_BS - f_TS - f_B - f_T - f_BN - f_TN + f_SE + f_E + f_NE) / (c1o1 + drho_NWT) -
+ ((vx1_NWT * vx1_NWT - vx3_NWT * vx3_NWT)));
+
+ //////////////////////////////////////////////////////////////////////////
+ // NET//
+ //////////////////////////////////////////////////////////////////////////
+ // index
+ kzero = kw;
+ kw = neighborFX[kw];
+ ks = ksw;
+ kb = kbw;
+ ksw = neighborFX[ksw];
+ kbw = neighborFX[kbw];
+ kbs = kbsw;
+ kbsw = neighborFX[kbsw];
+ ////////////////////////////////////////////////////////////////////////////////
+ f_E = feF[kzero];
+ f_W = fwF[kw];
+ f_N = fnF[kzero];
+ f_S = fsF[ks];
+ f_T = ftF[kzero];
+ f_B = fbF[kb];
+ f_NE = fneF[kzero];
+ f_SW = fswF[ksw];
+ f_SE = fseF[ks];
+ f_NW = fnwF[kw];
+ f_TE = fteF[kzero];
+ f_BW = fbwF[kbw];
+ f_BE = fbeF[kb];
+ f_TW = ftwF[kw];
+ f_TN = ftnF[kzero];
+ f_BS = fbsF[kbs];
+ f_BN = fbnF[kb];
+ f_TS = ftsF[ks];
+ f_ZERO = fzeroF[kzero];
+ f_TNE = ftneF[kzero];
+ f_TSW = ftswF[ksw];
+ f_TSE = ftseF[ks];
+ f_TNW = ftnwF[kw];
+ f_BNE = fbneF[kb];
+ f_BSW = fbswF[kbsw];
+ f_BSE = fbseF[kbs];
+ f_BNW = fbnwF[kbw];
+
+ drho_NET = f_E + f_W + f_N + f_S + f_T + f_B + f_NE + f_SW + f_SE + f_NW + f_TE + f_BW + f_BE + f_TW + f_TN +
+ f_BS + f_BN + f_TS + f_ZERO + f_TNE + f_TSW + f_TSE + f_TNW + f_BNE + f_BSW + f_BSE + f_BNW;
+ vx1_NET = (((f_TNE - f_BSW) + (f_TSE - f_BNW) + (f_BNE - f_TSW) + (f_BSE - f_TNW)) +
+ (((f_NE - f_SW) + (f_TE - f_BW)) + ((f_SE - f_NW) + (f_BE - f_TW))) + (f_E - f_W)) /
+ (c1o1 + drho_NET);
+ vx2_NET = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_BNE - f_TSW) + (f_BNW - f_TSE)) +
+ (((f_NE - f_SW) + (f_TN - f_BS)) + ((f_BN - f_TS) + (f_NW - f_SE))) + (f_N - f_S)) /
+ (c1o1 + drho_NET);
+ vx3_NET = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_TSE - f_BNW) + (f_TSW - f_BNE)) +
+ (((f_TE - f_BW) + (f_TN - f_BS)) + ((f_TW - f_BE) + (f_TS - f_BN))) + (f_T - f_B)) /
+ (c1o1 + drho_NET);
+
+ kxyFromfcNEQ_NET =
+ -c3o1 * omegaS *
+ ((f_SW + f_BSW + f_TSW - f_NW - f_BNW - f_TNW - f_SE - f_BSE - f_TSE + f_NE + f_BNE + f_TNE) /
+ (c1o1 + drho_NET) -
+ ((vx1_NET * vx2_NET)));
+ kyzFromfcNEQ_NET =
+ -c3o1 * omegaS *
+ ((f_BS + f_BSE + f_BSW - f_TS - f_TSE - f_TSW - f_BN - f_BNE - f_BNW + f_TN + f_TNE + f_TNW) /
+ (c1o1 + drho_NET) -
+ ((vx2_NET * vx3_NET)));
+ kxzFromfcNEQ_NET =
+ -c3o1 * omegaS *
+ ((f_BW + f_BSW + f_BNW - f_TW - f_TSW - f_TNW - f_BE - f_BSE - f_BNE + f_TE + f_TSE + f_TNE) /
+ (c1o1 + drho_NET) -
+ ((vx1_NET * vx3_NET)));
+ kxxMyyFromfcNEQ_NET =
+ -c3o2 * omegaS *
+ ((f_BW + f_W + f_TW - f_BS - f_S - f_TS - f_BN - f_N - f_TN + f_BE + f_E + f_TE) / (c1o1 + drho_NET) -
+ ((vx1_NET * vx1_NET - vx2_NET * vx2_NET)));
+ kxxMzzFromfcNEQ_NET =
+ -c3o2 * omegaS *
+ ((f_SW + f_W + f_NW - f_BS - f_TS - f_B - f_T - f_BN - f_TN + f_SE + f_E + f_NE) / (c1o1 + drho_NET) -
+ ((vx1_NET * vx1_NET - vx3_NET * vx3_NET)));
+
+ //////////////////////////////////////////////////////////////////////////
+ // NEB//
+ //////////////////////////////////////////////////////////////////////////
+ // index
+ kb = kzero;
+ kbw = kw;
+ kbs = ks;
+ kbsw = ksw;
+ kzero = k0w;
+ kw = neighborFX[k0w];
+ ks = k0sw;
+ ksw = neighborFX[k0sw];
+ ////////////////////////////////////////////////////////////////////////////////
+ f_E = feF[kzero];
+ f_W = fwF[kw];
+ f_N = fnF[kzero];
+ f_S = fsF[ks];
+ f_T = ftF[kzero];
+ f_B = fbF[kb];
+ f_NE = fneF[kzero];
+ f_SW = fswF[ksw];
+ f_SE = fseF[ks];
+ f_NW = fnwF[kw];
+ f_TE = fteF[kzero];
+ f_BW = fbwF[kbw];
+ f_BE = fbeF[kb];
+ f_TW = ftwF[kw];
+ f_TN = ftnF[kzero];
+ f_BS = fbsF[kbs];
+ f_BN = fbnF[kb];
+ f_TS = ftsF[ks];
+ f_ZERO = fzeroF[kzero];
+ f_TNE = ftneF[kzero];
+ f_TSW = ftswF[ksw];
+ f_TSE = ftseF[ks];
+ f_TNW = ftnwF[kw];
+ f_BNE = fbneF[kb];
+ f_BSW = fbswF[kbsw];
+ f_BSE = fbseF[kbs];
+ f_BNW = fbnwF[kbw];
+
+ drho_NEB = f_E + f_W + f_N + f_S + f_T + f_B + f_NE + f_SW + f_SE + f_NW + f_TE + f_BW + f_BE + f_TW + f_TN +
+ f_BS + f_BN + f_TS + f_ZERO + f_TNE + f_TSW + f_TSE + f_TNW + f_BNE + f_BSW + f_BSE + f_BNW;
+ vx1_NEB = (((f_TNE - f_BSW) + (f_TSE - f_BNW) + (f_BNE - f_TSW) + (f_BSE - f_TNW)) +
+ (((f_NE - f_SW) + (f_TE - f_BW)) + ((f_SE - f_NW) + (f_BE - f_TW))) + (f_E - f_W)) /
+ (c1o1 + drho_NEB);
+ vx2_NEB = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_BNE - f_TSW) + (f_BNW - f_TSE)) +
+ (((f_NE - f_SW) + (f_TN - f_BS)) + ((f_BN - f_TS) + (f_NW - f_SE))) + (f_N - f_S)) /
+ (c1o1 + drho_NEB);
+ vx3_NEB = (((f_TNE - f_BSW) + (f_TNW - f_BSE) + (f_TSE - f_BNW) + (f_TSW - f_BNE)) +
+ (((f_TE - f_BW) + (f_TN - f_BS)) + ((f_TW - f_BE) + (f_TS - f_BN))) + (f_T - f_B)) /
+ (c1o1 + drho_NEB);
+
+ kxyFromfcNEQ_NEB =
+ -c3o1 * omegaS *
+ ((f_SW + f_BSW + f_TSW - f_NW - f_BNW - f_TNW - f_SE - f_BSE - f_TSE + f_NE + f_BNE + f_TNE) /
+ (c1o1 + drho_NEB) -
+ ((vx1_NEB * vx2_NEB)));
+ kyzFromfcNEQ_NEB =
+ -c3o1 * omegaS *
+ ((f_BS + f_BSE + f_BSW - f_TS - f_TSE - f_TSW - f_BN - f_BNE - f_BNW + f_TN + f_TNE + f_TNW) /
+ (c1o1 + drho_NEB) -
+ ((vx2_NEB * vx3_NEB)));
+ kxzFromfcNEQ_NEB =
+ -c3o1 * omegaS *
+ ((f_BW + f_BSW + f_BNW - f_TW - f_TSW - f_TNW - f_BE - f_BSE - f_BNE + f_TE + f_TSE + f_TNE) /
+ (c1o1 + drho_NEB) -
+ ((vx1_NEB * vx3_NEB)));
+ kxxMyyFromfcNEQ_NEB =
+ -c3o2 * omegaS *
+ ((f_BW + f_W + f_TW - f_BS - f_S - f_TS - f_BN - f_N - f_TN + f_BE + f_E + f_TE) / (c1o1 + drho_NEB) -
+ ((vx1_NEB * vx1_NEB - vx2_NEB * vx2_NEB)));
+ kxxMzzFromfcNEQ_NEB =
+ -c3o2 * omegaS *
+ ((f_SW + f_W + f_NW - f_BS - f_TS - f_B - f_T - f_BN - f_TN + f_SE + f_E + f_NE) / (c1o1 + drho_NEB) -
+ ((vx1_NEB * vx1_NEB - vx3_NEB * vx3_NEB)));
+
+ //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ // kxyFromfcNEQ_SWB = zero;
+ // kyzFromfcNEQ_SWB = zero;
+ // kxzFromfcNEQ_SWB = zero;
+ // kxxMyyFromfcNEQ_SWB = zero;
+ // kxxMzzFromfcNEQ_SWB = zero;
+ // kxyFromfcNEQ_SWT = zero;
+ // kyzFromfcNEQ_SWT = zero;
+ // kxzFromfcNEQ_SWT = zero;
+ // kxxMyyFromfcNEQ_SWT = zero;
+ // kxxMzzFromfcNEQ_SWT = zero;
+ // kxyFromfcNEQ_SET = zero;
+ // kyzFromfcNEQ_SET = zero;
+ // kxzFromfcNEQ_SET = zero;
+ // kxxMyyFromfcNEQ_SET = zero;
+ // kxxMzzFromfcNEQ_SET = zero;
+ // kxyFromfcNEQ_SEB = zero;
+ // kyzFromfcNEQ_SEB = zero;
+ // kxzFromfcNEQ_SEB = zero;
+ // kxxMyyFromfcNEQ_SEB = zero;
+ // kxxMzzFromfcNEQ_SEB = zero;
+ // kxyFromfcNEQ_NWB = zero;
+ // kyzFromfcNEQ_NWB = zero;
+ // kxzFromfcNEQ_NWB = zero;
+ // kxxMyyFromfcNEQ_NWB = zero;
+ // kxxMzzFromfcNEQ_NWB = zero;
+ // kxyFromfcNEQ_NWT = zero;
+ // kyzFromfcNEQ_NWT = zero;
+ // kxzFromfcNEQ_NWT = zero;
+ // kxxMyyFromfcNEQ_NWT = zero;
+ // kxxMzzFromfcNEQ_NWT = zero;
+ // kxyFromfcNEQ_NET = zero;
+ // kyzFromfcNEQ_NET = zero;
+ // kxzFromfcNEQ_NET = zero;
+ // kxxMyyFromfcNEQ_NET = zero;
+ // kxxMzzFromfcNEQ_NET = zero;
+ // kxyFromfcNEQ_NEB = zero;
+ // kyzFromfcNEQ_NEB = zero;
+ // kxzFromfcNEQ_NEB = zero;
+ // kxxMyyFromfcNEQ_NEB = zero;
+ // kxxMzzFromfcNEQ_NEB = zero;
+ //////////////////////////////////////////////////////////////////////////
+ // 3
+ //////////////////////////////////////////////////////////////////////////
+ a0 = (-kxxMyyFromfcNEQ_NEB - kxxMyyFromfcNEQ_NET + kxxMyyFromfcNEQ_NWB + kxxMyyFromfcNEQ_NWT -
+ kxxMyyFromfcNEQ_SEB - kxxMyyFromfcNEQ_SET + kxxMyyFromfcNEQ_SWB + kxxMyyFromfcNEQ_SWT -
+ kxxMzzFromfcNEQ_NEB - kxxMzzFromfcNEQ_NET + kxxMzzFromfcNEQ_NWB + kxxMzzFromfcNEQ_NWT -
+ kxxMzzFromfcNEQ_SEB - kxxMzzFromfcNEQ_SET + kxxMzzFromfcNEQ_SWB + kxxMzzFromfcNEQ_SWT -
+ c2o1 * kxyFromfcNEQ_NEB - c2o1 * kxyFromfcNEQ_NET - c2o1 * kxyFromfcNEQ_NWB - c2o1 * kxyFromfcNEQ_NWT +
+ c2o1 * kxyFromfcNEQ_SEB + c2o1 * kxyFromfcNEQ_SET + c2o1 * kxyFromfcNEQ_SWB + c2o1 * kxyFromfcNEQ_SWT +
+ c2o1 * kxzFromfcNEQ_NEB - c2o1 * kxzFromfcNEQ_NET + c2o1 * kxzFromfcNEQ_NWB - c2o1 * kxzFromfcNEQ_NWT +
+ c2o1 * kxzFromfcNEQ_SEB - c2o1 * kxzFromfcNEQ_SET + c2o1 * kxzFromfcNEQ_SWB - c2o1 * kxzFromfcNEQ_SWT +
+ c8o1 * vx1_NEB + c8o1 * vx1_NET + c8o1 * vx1_NWB + c8o1 * vx1_NWT + c8o1 * vx1_SEB + c8o1 * vx1_SET +
+ c8o1 * vx1_SWB + c8o1 * vx1_SWT + c2o1 * vx2_NEB + c2o1 * vx2_NET - c2o1 * vx2_NWB - c2o1 * vx2_NWT -
+ c2o1 * vx2_SEB - c2o1 * vx2_SET + c2o1 * vx2_SWB + c2o1 * vx2_SWT - c2o1 * vx3_NEB + c2o1 * vx3_NET +
+ c2o1 * vx3_NWB - c2o1 * vx3_NWT - c2o1 * vx3_SEB + c2o1 * vx3_SET + c2o1 * vx3_SWB - c2o1 * vx3_SWT) /
+ c64o1;
+ b0 = (c2o1 * kxxMyyFromfcNEQ_NEB + c2o1 * kxxMyyFromfcNEQ_NET + c2o1 * kxxMyyFromfcNEQ_NWB +
+ c2o1 * kxxMyyFromfcNEQ_NWT - c2o1 * kxxMyyFromfcNEQ_SEB - c2o1 * kxxMyyFromfcNEQ_SET -
+ c2o1 * kxxMyyFromfcNEQ_SWB - c2o1 * kxxMyyFromfcNEQ_SWT - kxxMzzFromfcNEQ_NEB - kxxMzzFromfcNEQ_NET -
+ kxxMzzFromfcNEQ_NWB - kxxMzzFromfcNEQ_NWT + kxxMzzFromfcNEQ_SEB + kxxMzzFromfcNEQ_SET +
+ kxxMzzFromfcNEQ_SWB + kxxMzzFromfcNEQ_SWT - c2o1 * kxyFromfcNEQ_NEB - c2o1 * kxyFromfcNEQ_NET +
+ c2o1 * kxyFromfcNEQ_NWB + c2o1 * kxyFromfcNEQ_NWT - c2o1 * kxyFromfcNEQ_SEB - c2o1 * kxyFromfcNEQ_SET +
+ c2o1 * kxyFromfcNEQ_SWB + c2o1 * kxyFromfcNEQ_SWT + c2o1 * kyzFromfcNEQ_NEB - c2o1 * kyzFromfcNEQ_NET +
+ c2o1 * kyzFromfcNEQ_NWB - c2o1 * kyzFromfcNEQ_NWT + c2o1 * kyzFromfcNEQ_SEB - c2o1 * kyzFromfcNEQ_SET +
+ c2o1 * kyzFromfcNEQ_SWB - c2o1 * kyzFromfcNEQ_SWT + c2o1 * vx1_NEB + c2o1 * vx1_NET - c2o1 * vx1_NWB -
+ c2o1 * vx1_NWT - c2o1 * vx1_SEB - c2o1 * vx1_SET + c2o1 * vx1_SWB + c2o1 * vx1_SWT + c8o1 * vx2_NEB +
+ c8o1 * vx2_NET + c8o1 * vx2_NWB + c8o1 * vx2_NWT + c8o1 * vx2_SEB + c8o1 * vx2_SET + c8o1 * vx2_SWB +
+ c8o1 * vx2_SWT - c2o1 * vx3_NEB + c2o1 * vx3_NET - c2o1 * vx3_NWB + c2o1 * vx3_NWT + c2o1 * vx3_SEB -
+ c2o1 * vx3_SET + c2o1 * vx3_SWB - c2o1 * vx3_SWT) /
+ c64o1;
+ c0 = (kxxMyyFromfcNEQ_NEB - kxxMyyFromfcNEQ_NET + kxxMyyFromfcNEQ_NWB - kxxMyyFromfcNEQ_NWT +
+ kxxMyyFromfcNEQ_SEB - kxxMyyFromfcNEQ_SET + kxxMyyFromfcNEQ_SWB - kxxMyyFromfcNEQ_SWT -
+ c2o1 * kxxMzzFromfcNEQ_NEB + c2o1 * kxxMzzFromfcNEQ_NET - c2o1 * kxxMzzFromfcNEQ_NWB +
+ c2o1 * kxxMzzFromfcNEQ_NWT - c2o1 * kxxMzzFromfcNEQ_SEB + c2o1 * kxxMzzFromfcNEQ_SET -
+ c2o1 * kxxMzzFromfcNEQ_SWB + c2o1 * kxxMzzFromfcNEQ_SWT - c2o1 * kxzFromfcNEQ_NEB -
+ c2o1 * kxzFromfcNEQ_NET + c2o1 * kxzFromfcNEQ_NWB + c2o1 * kxzFromfcNEQ_NWT - c2o1 * kxzFromfcNEQ_SEB -
+ c2o1 * kxzFromfcNEQ_SET + c2o1 * kxzFromfcNEQ_SWB + c2o1 * kxzFromfcNEQ_SWT - c2o1 * kyzFromfcNEQ_NEB -
+ c2o1 * kyzFromfcNEQ_NET - c2o1 * kyzFromfcNEQ_NWB - c2o1 * kyzFromfcNEQ_NWT + c2o1 * kyzFromfcNEQ_SEB +
+ c2o1 * kyzFromfcNEQ_SET + c2o1 * kyzFromfcNEQ_SWB + c2o1 * kyzFromfcNEQ_SWT - c2o1 * vx1_NEB +
+ c2o1 * vx1_NET + c2o1 * vx1_NWB - c2o1 * vx1_NWT - c2o1 * vx1_SEB + c2o1 * vx1_SET + c2o1 * vx1_SWB -
+ c2o1 * vx1_SWT - c2o1 * vx2_NEB + c2o1 * vx2_NET - c2o1 * vx2_NWB + c2o1 * vx2_NWT + c2o1 * vx2_SEB -
+ c2o1 * vx2_SET + c2o1 * vx2_SWB - c2o1 * vx2_SWT + c8o1 * vx3_NEB + c8o1 * vx3_NET + c8o1 * vx3_NWB +
+ c8o1 * vx3_NWT + c8o1 * vx3_SEB + c8o1 * vx3_SET + c8o1 * vx3_SWB + c8o1 * vx3_SWT) /
+ c64o1;
+ ax = (vx1_NEB + vx1_NET - vx1_NWB - vx1_NWT + vx1_SEB + vx1_SET - vx1_SWB - vx1_SWT) / c4o1;
+ bx = (vx2_NEB + vx2_NET - vx2_NWB - vx2_NWT + vx2_SEB + vx2_SET - vx2_SWB - vx2_SWT) / c4o1;
+ cx = (vx3_NEB + vx3_NET - vx3_NWB - vx3_NWT + vx3_SEB + vx3_SET - vx3_SWB - vx3_SWT) / c4o1;
+ axx = (kxxMyyFromfcNEQ_NEB + kxxMyyFromfcNEQ_NET - kxxMyyFromfcNEQ_NWB - kxxMyyFromfcNEQ_NWT +
+ kxxMyyFromfcNEQ_SEB + kxxMyyFromfcNEQ_SET - kxxMyyFromfcNEQ_SWB - kxxMyyFromfcNEQ_SWT +
+ kxxMzzFromfcNEQ_NEB + kxxMzzFromfcNEQ_NET - kxxMzzFromfcNEQ_NWB - kxxMzzFromfcNEQ_NWT +
+ kxxMzzFromfcNEQ_SEB + kxxMzzFromfcNEQ_SET - kxxMzzFromfcNEQ_SWB - kxxMzzFromfcNEQ_SWT + c2o1 * vx2_NEB +
+ c2o1 * vx2_NET - c2o1 * vx2_NWB - c2o1 * vx2_NWT - c2o1 * vx2_SEB - c2o1 * vx2_SET + c2o1 * vx2_SWB +
+ c2o1 * vx2_SWT - c2o1 * vx3_NEB + c2o1 * vx3_NET + c2o1 * vx3_NWB - c2o1 * vx3_NWT - c2o1 * vx3_SEB +
+ c2o1 * vx3_SET + c2o1 * vx3_SWB - c2o1 * vx3_SWT) /
+ c16o1;
+ bxx = (kxyFromfcNEQ_NEB + kxyFromfcNEQ_NET - kxyFromfcNEQ_NWB - kxyFromfcNEQ_NWT + kxyFromfcNEQ_SEB +
+ kxyFromfcNEQ_SET - kxyFromfcNEQ_SWB - kxyFromfcNEQ_SWT - c2o1 * vx1_NEB - c2o1 * vx1_NET +
+ c2o1 * vx1_NWB + c2o1 * vx1_NWT + c2o1 * vx1_SEB + c2o1 * vx1_SET - c2o1 * vx1_SWB - c2o1 * vx1_SWT) /
+ c8o1;
+ cxx = (kxzFromfcNEQ_NEB + kxzFromfcNEQ_NET - kxzFromfcNEQ_NWB - kxzFromfcNEQ_NWT + kxzFromfcNEQ_SEB +
+ kxzFromfcNEQ_SET - kxzFromfcNEQ_SWB - kxzFromfcNEQ_SWT + c2o1 * vx1_NEB - c2o1 * vx1_NET -
+ c2o1 * vx1_NWB + c2o1 * vx1_NWT + c2o1 * vx1_SEB - c2o1 * vx1_SET - c2o1 * vx1_SWB + c2o1 * vx1_SWT) /
+ c8o1;
+ ay = (vx1_NEB + vx1_NET + vx1_NWB + vx1_NWT - vx1_SEB - vx1_SET - vx1_SWB - vx1_SWT) / c4o1;
+ by = (vx2_NEB + vx2_NET + vx2_NWB + vx2_NWT - vx2_SEB - vx2_SET - vx2_SWB - vx2_SWT) / c4o1;
+ cy = (vx3_NEB + vx3_NET + vx3_NWB + vx3_NWT - vx3_SEB - vx3_SET - vx3_SWB - vx3_SWT) / c4o1;
+ ayy = (kxyFromfcNEQ_NEB + kxyFromfcNEQ_NET + kxyFromfcNEQ_NWB + kxyFromfcNEQ_NWT - kxyFromfcNEQ_SEB -
+ kxyFromfcNEQ_SET - kxyFromfcNEQ_SWB - kxyFromfcNEQ_SWT - c2o1 * vx2_NEB - c2o1 * vx2_NET +
+ c2o1 * vx2_NWB + c2o1 * vx2_NWT + c2o1 * vx2_SEB + c2o1 * vx2_SET - c2o1 * vx2_SWB - c2o1 * vx2_SWT) /
+ c8o1;
+ byy = (-c2o1 * kxxMyyFromfcNEQ_NEB - c2o1 * kxxMyyFromfcNEQ_NET - c2o1 * kxxMyyFromfcNEQ_NWB -
+ c2o1 * kxxMyyFromfcNEQ_NWT + c2o1 * kxxMyyFromfcNEQ_SEB + c2o1 * kxxMyyFromfcNEQ_SET +
+ c2o1 * kxxMyyFromfcNEQ_SWB + c2o1 * kxxMyyFromfcNEQ_SWT + kxxMzzFromfcNEQ_NEB + kxxMzzFromfcNEQ_NET +
+ kxxMzzFromfcNEQ_NWB + kxxMzzFromfcNEQ_NWT - kxxMzzFromfcNEQ_SEB - kxxMzzFromfcNEQ_SET -
+ kxxMzzFromfcNEQ_SWB - kxxMzzFromfcNEQ_SWT + c2o1 * vx1_NEB + c2o1 * vx1_NET - c2o1 * vx1_NWB -
+ c2o1 * vx1_NWT - c2o1 * vx1_SEB - c2o1 * vx1_SET + c2o1 * vx1_SWB + c2o1 * vx1_SWT - c2o1 * vx3_NEB +
+ c2o1 * vx3_NET - c2o1 * vx3_NWB + c2o1 * vx3_NWT + c2o1 * vx3_SEB - c2o1 * vx3_SET + c2o1 * vx3_SWB -
+ c2o1 * vx3_SWT) /
+ c16o1;
+ cyy = (kyzFromfcNEQ_NEB + kyzFromfcNEQ_NET + kyzFromfcNEQ_NWB + kyzFromfcNEQ_NWT - kyzFromfcNEQ_SEB -
+ kyzFromfcNEQ_SET - kyzFromfcNEQ_SWB - kyzFromfcNEQ_SWT + c2o1 * vx2_NEB - c2o1 * vx2_NET +
+ c2o1 * vx2_NWB - c2o1 * vx2_NWT - c2o1 * vx2_SEB + c2o1 * vx2_SET - c2o1 * vx2_SWB + c2o1 * vx2_SWT) /
+ c8o1;
+ az = (-vx1_NEB + vx1_NET - vx1_NWB + vx1_NWT - vx1_SEB + vx1_SET - vx1_SWB + vx1_SWT) / c4o1;
+ bz = (-vx2_NEB + vx2_NET - vx2_NWB + vx2_NWT - vx2_SEB + vx2_SET - vx2_SWB + vx2_SWT) / c4o1;
+ cz = (-vx3_NEB + vx3_NET - vx3_NWB + vx3_NWT - vx3_SEB + vx3_SET - vx3_SWB + vx3_SWT) / c4o1;
+ azz = (-kxzFromfcNEQ_NEB + kxzFromfcNEQ_NET - kxzFromfcNEQ_NWB + kxzFromfcNEQ_NWT - kxzFromfcNEQ_SEB +
+ kxzFromfcNEQ_SET - kxzFromfcNEQ_SWB + kxzFromfcNEQ_SWT + c2o1 * vx3_NEB - c2o1 * vx3_NET -
+ c2o1 * vx3_NWB + c2o1 * vx3_NWT + c2o1 * vx3_SEB - c2o1 * vx3_SET - c2o1 * vx3_SWB + c2o1 * vx3_SWT) /
+ c8o1;
+ bzz = (-kyzFromfcNEQ_NEB + kyzFromfcNEQ_NET - kyzFromfcNEQ_NWB + kyzFromfcNEQ_NWT - kyzFromfcNEQ_SEB +
+ kyzFromfcNEQ_SET - kyzFromfcNEQ_SWB + kyzFromfcNEQ_SWT + c2o1 * vx3_NEB - c2o1 * vx3_NET +
+ c2o1 * vx3_NWB - c2o1 * vx3_NWT - c2o1 * vx3_SEB + c2o1 * vx3_SET - c2o1 * vx3_SWB + c2o1 * vx3_SWT) /
+ c8o1;
+ czz = (-kxxMyyFromfcNEQ_NEB + kxxMyyFromfcNEQ_NET - kxxMyyFromfcNEQ_NWB + kxxMyyFromfcNEQ_NWT -
+ kxxMyyFromfcNEQ_SEB + kxxMyyFromfcNEQ_SET - kxxMyyFromfcNEQ_SWB + kxxMyyFromfcNEQ_SWT +
+ c2o1 * kxxMzzFromfcNEQ_NEB - c2o1 * kxxMzzFromfcNEQ_NET + c2o1 * kxxMzzFromfcNEQ_NWB -
+ c2o1 * kxxMzzFromfcNEQ_NWT + c2o1 * kxxMzzFromfcNEQ_SEB - c2o1 * kxxMzzFromfcNEQ_SET +
+ c2o1 * kxxMzzFromfcNEQ_SWB - c2o1 * kxxMzzFromfcNEQ_SWT - c2o1 * vx1_NEB + c2o1 * vx1_NET +
+ c2o1 * vx1_NWB - c2o1 * vx1_NWT - c2o1 * vx1_SEB + c2o1 * vx1_SET + c2o1 * vx1_SWB - c2o1 * vx1_SWT -
+ c2o1 * vx2_NEB + c2o1 * vx2_NET - c2o1 * vx2_NWB + c2o1 * vx2_NWT + c2o1 * vx2_SEB - c2o1 * vx2_SET +
+ c2o1 * vx2_SWB - c2o1 * vx2_SWT) /
+ c16o1;
+ axy = (vx1_NEB + vx1_NET - vx1_NWB - vx1_NWT - vx1_SEB - vx1_SET + vx1_SWB + vx1_SWT) / c2o1;
+ bxy = (vx2_NEB + vx2_NET - vx2_NWB - vx2_NWT - vx2_SEB - vx2_SET + vx2_SWB + vx2_SWT) / c2o1;
+ cxy = (vx3_NEB + vx3_NET - vx3_NWB - vx3_NWT - vx3_SEB - vx3_SET + vx3_SWB + vx3_SWT) / c2o1;
+ axz = (-vx1_NEB + vx1_NET + vx1_NWB - vx1_NWT - vx1_SEB + vx1_SET + vx1_SWB - vx1_SWT) / c2o1;
+ bxz = (-vx2_NEB + vx2_NET + vx2_NWB - vx2_NWT - vx2_SEB + vx2_SET + vx2_SWB - vx2_SWT) / c2o1;
+ cxz = (-vx3_NEB + vx3_NET + vx3_NWB - vx3_NWT - vx3_SEB + vx3_SET + vx3_SWB - vx3_SWT) / c2o1;
+ ayz = (-vx1_NEB + vx1_NET - vx1_NWB + vx1_NWT + vx1_SEB - vx1_SET + vx1_SWB - vx1_SWT) / c2o1;
+ byz = (-vx2_NEB + vx2_NET - vx2_NWB + vx2_NWT + vx2_SEB - vx2_SET + vx2_SWB - vx2_SWT) / c2o1;
+ cyz = (-vx3_NEB + vx3_NET - vx3_NWB + vx3_NWT + vx3_SEB - vx3_SET + vx3_SWB - vx3_SWT) / c2o1;
+ // axyz=-vx1_NEB + vx1_NET + vx1_NWB - vx1_NWT + vx1_SEB - vx1_SET - vx1_SWB + vx1_SWT;
+ // bxyz=-vx2_NEB + vx2_NET + vx2_NWB - vx2_NWT + vx2_SEB - vx2_SET - vx2_SWB + vx2_SWT;
+ // cxyz=-vx3_NEB + vx3_NET + vx3_NWB - vx3_NWT + vx3_SEB - vx3_SET - vx3_SWB + vx3_SWT;
+ //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ real kxyAverage = c0o1;
+ real kyzAverage = c0o1;
+ real kxzAverage = c0o1;
+ real kxxMyyAverage = c0o1;
+ real kxxMzzAverage = c0o1;
+ // real kxyAverage =(kxyFromfcNEQ_SWB+
+ // kxyFromfcNEQ_SWT+
+ // kxyFromfcNEQ_SET+
+ // kxyFromfcNEQ_SEB+
+ // kxyFromfcNEQ_NWB+
+ // kxyFromfcNEQ_NWT+
+ // kxyFromfcNEQ_NET+
+ // kxyFromfcNEQ_NEB)*c1o8-(ay+bx);
+ // real kyzAverage =(kyzFromfcNEQ_SWB+
+ // kyzFromfcNEQ_SWT+
+ // kyzFromfcNEQ_SET+
+ // kyzFromfcNEQ_SEB+
+ // kyzFromfcNEQ_NWB+
+ // kyzFromfcNEQ_NWT+
+ // kyzFromfcNEQ_NET+
+ // kyzFromfcNEQ_NEB)*c1o8-(bz+cy);
+ // real kxzAverage =(kxzFromfcNEQ_SWB+
+ // kxzFromfcNEQ_SWT+
+ // kxzFromfcNEQ_SET+
+ // kxzFromfcNEQ_SEB+
+ // kxzFromfcNEQ_NWB+
+ // kxzFromfcNEQ_NWT+
+ // kxzFromfcNEQ_NET+
+ // kxzFromfcNEQ_NEB)*c1o8-(az+cx);
+ // real kxxMyyAverage =(kxxMyyFromfcNEQ_SWB+
+ // kxxMyyFromfcNEQ_SWT+
+ // kxxMyyFromfcNEQ_SET+
+ // kxxMyyFromfcNEQ_SEB+
+ // kxxMyyFromfcNEQ_NWB+
+ // kxxMyyFromfcNEQ_NWT+
+ // kxxMyyFromfcNEQ_NET+
+ // kxxMyyFromfcNEQ_NEB)*c1o8-(ax-by);
+ // real kxxMzzAverage =(kxxMzzFromfcNEQ_SWB+
+ // kxxMzzFromfcNEQ_SWT+
+ // kxxMzzFromfcNEQ_SET+
+ // kxxMzzFromfcNEQ_SEB+
+ // kxxMzzFromfcNEQ_NWB+
+ // kxxMzzFromfcNEQ_NWT+
+ // kxxMzzFromfcNEQ_NET+
+ // kxxMzzFromfcNEQ_NEB)*c1o8-(ax-cz);
+
+ //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ ////Press
+ // d0 = ( press_NEB + press_NET + press_NWB + press_NWT + press_SEB + press_SET + press_SWB + press_SWT) *
+ // c1o8; dx = ( press_NEB + press_NET - press_NWB - press_NWT + press_SEB + press_SET - press_SWB - press_SWT)
+ // * c1o4; dy = ( press_NEB + press_NET + press_NWB + press_NWT - press_SEB - press_SET - press_SWB -
+ // press_SWT) * c1o4; dz = (-press_NEB + press_NET - press_NWB + press_NWT - press_SEB + press_SET - press_SWB
+ // + press_SWT) * c1o4; dxy = ( press_NEB + press_NET - press_NWB - press_NWT - press_SEB - press_SET +
+ // press_SWB + press_SWT) * c1o2; dxz = (-press_NEB + press_NET + press_NWB - press_NWT - press_SEB + press_SET
+ // + press_SWB - press_SWT) * c1o2; dyz = (-press_NEB + press_NET - press_NWB + press_NWT + press_SEB -
+ // press_SET + press_SWB - press_SWT) * c1o2; dxyz = -press_NEB + press_NET + press_NWB - press_NWT + press_SEB
+ // - press_SET - press_SWB + press_SWT;
+ //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ // drho
+ real LapRho = ((xoff != c0o1) || (yoff != c0o1) || (zoff != c0o1))
+ ? c0o1
+ : -c3o1 * (ax * ax + by * by + cz * cz) - c6o1 * (bx * ay + cx * az + cy * bz);
+ d0 = (drho_NEB + drho_NET + drho_NWB + drho_NWT + drho_SEB + drho_SET + drho_SWB + drho_SWT - c2o1 * LapRho) *
+ c1o8;
+ dx = (drho_NEB + drho_NET - drho_NWB - drho_NWT + drho_SEB + drho_SET - drho_SWB - drho_SWT) * c1o4;
+ dy = (drho_NEB + drho_NET + drho_NWB + drho_NWT - drho_SEB - drho_SET - drho_SWB - drho_SWT) * c1o4;
+ dz = (-drho_NEB + drho_NET - drho_NWB + drho_NWT - drho_SEB + drho_SET - drho_SWB + drho_SWT) * c1o4;
+ dxy = (drho_NEB + drho_NET - drho_NWB - drho_NWT - drho_SEB - drho_SET + drho_SWB + drho_SWT) * c1o2;
+ dxz = (-drho_NEB + drho_NET + drho_NWB - drho_NWT - drho_SEB + drho_SET + drho_SWB - drho_SWT) * c1o2;
+ dyz = (-drho_NEB + drho_NET - drho_NWB + drho_NWT + drho_SEB - drho_SET + drho_SWB - drho_SWT) * c1o2;
+ // dxyz = -drho_NEB + drho_NET + drho_NWB - drho_NWT + drho_SEB - drho_SET - drho_SWB + drho_SWT;
+ // d0 = zero;
+ // dx = zero;
+ // dy = zero;
+ // dz = zero;
+ // dxy = zero;
+ // dxz = zero;
+ // dyz = zero;
+ // dxyz = zero;
+ ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ //
+ // Bernd das Brot
+ //
+ //
+ // x------x
+ // | |
+ // | ---+--->X
+ // | | \
+ // x------x \
+ // off-vector
+ ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ a0 = a0 + xoff * ax + yoff * ay + zoff * az + xoff_sq * axx + yoff_sq * ayy + zoff_sq * azz +
+ xoff * yoff * axy + xoff * zoff * axz + yoff * zoff * ayz;
+ ax = ax + c2o1 * xoff * axx + yoff * axy + zoff * axz;
+ ay = ay + c2o1 * yoff * ayy + xoff * axy + zoff * ayz;
+ az = az + c2o1 * zoff * azz + xoff * axz + yoff * ayz;
+ b0 = b0 + xoff * bx + yoff * by + zoff * bz + xoff_sq * bxx + yoff_sq * byy + zoff_sq * bzz +
+ xoff * yoff * bxy + xoff * zoff * bxz + yoff * zoff * byz;
+ bx = bx + c2o1 * xoff * bxx + yoff * bxy + zoff * bxz;
+ by = by + c2o1 * yoff * byy + xoff * bxy + zoff * byz;
+ bz = bz + c2o1 * zoff * bzz + xoff * bxz + yoff * byz;
+ c0 = c0 + xoff * cx + yoff * cy + zoff * cz + xoff_sq * cxx + yoff_sq * cyy + zoff_sq * czz +
+ xoff * yoff * cxy + xoff * zoff * cxz + yoff * zoff * cyz;
+ cx = cx + c2o1 * xoff * cxx + yoff * cxy + zoff * cxz;
+ cy = cy + c2o1 * yoff * cyy + xoff * cxy + zoff * cyz;
+ cz = cz + c2o1 * zoff * czz + xoff * cxz + yoff * cyz;
+ d0 = d0 + xoff * dx + yoff * dy + zoff * dz + xoff * yoff * dxy + xoff * zoff * dxz + yoff * zoff * dyz;
+ ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ // FIX
+ // ///////////////////////////////////////////////////////////////////////////////////////////////////////////
+ ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ // AAAAAAAAAAAAHHHHHHHHHHHH!!!!! Mieser Test!!!
+ // b0= bx= by= bz= bxx= byy= bzz= bxy= bxz= byz= c0= cx= cy= cz= cxx= cyy= czz= cxy= cxz= cyz= axyz= bxyz=
+ // cxyz=zero; b0=zero; bx=zero; by=zero; bz=zero; bxx=zero; byy=zero; bzz=zero; bxy=zero; bxz=zero; byz=zero;
+ // c0=zero;
+ // cx=zero;
+ // cy=zero;
+ // cz=zero;
+ // cxx=zero;
+ // cyy=zero;
+ // czz=zero;
+ // cxy=zero;
+ // cxz=zero;
+ // cyz=zero;
+ ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ real mfcbb = c0o1;
+ real mfabb = c0o1;
+ real mfbcb = c0o1;
+ real mfbab = c0o1;
+ real mfbbc = c0o1;
+ real mfbba = c0o1;
+ real mfccb = c0o1;
+ real mfaab = c0o1;
+ real mfcab = c0o1;
+ real mfacb = c0o1;
+ real mfcbc = c0o1;
+ real mfaba = c0o1;
+ real mfcba = c0o1;
+ real mfabc = c0o1;
+ real mfbcc = c0o1;
+ real mfbaa = c0o1;
+ real mfbca = c0o1;
+ real mfbac = c0o1;
+ real mfbbb = c0o1;
+ real mfccc = c0o1;
+ real mfaac = c0o1;
+ real mfcac = c0o1;
+ real mfacc = c0o1;
+ real mfcca = c0o1;
+ real mfaaa = c0o1;
+ real mfcaa = c0o1;
+ real mfaca = c0o1;
+ ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ real m0, m1, m2, vvx, vvy, vvz, vx2, vy2, vz2, oMdrho;
+ real mxxPyyPzz, mxxMyy, mxxMzz, mxxyPyzz, mxxyMyzz, mxxzPyyz, mxxzMyyz, mxyyPxzz, mxyyMxzz;
+ // real qudricLimit = c1o100;//ganz schlechte Idee -> muss global sein
+ // real O3 = c2o1 - o;
+ // real residu, residutmp;
+ // residutmp = c0o1;///*-*/ c2o9 * (1./o - c1o2) * eps_new * eps_new;
+ real NeqOn = c1o1; // zero;//one; //.... one = on ..... zero = off
+ ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+
+ ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ //
+ // Position C 0., 0., 0.
+ //
+ ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ // x = 0.;
+ // y = 0.;
+ // z = 0.;
+ ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ // real mxoff = -xoff;
+ // real myoff = -yoff;
+ // real mzoff = -zoff;
+ ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ // press = press_NET * (c1o8 - c1o4 * mxoff - c1o4 * myoff - c1o4 * mzoff) +
+ // press_NWT * (c1o8 + c1o4 * mxoff - c1o4 * myoff - c1o4 * mzoff) +
+ // press_SET * (c1o8 - c1o4 * mxoff + c1o4 * myoff - c1o4 * mzoff) +
+ // press_SWT * (c1o8 + c1o4 * mxoff + c1o4 * myoff - c1o4 * mzoff) +
+ // press_NEB * (c1o8 - c1o4 * mxoff - c1o4 * myoff + c1o4 * mzoff) +
+ // press_NWB * (c1o8 + c1o4 * mxoff - c1o4 * myoff + c1o4 * mzoff) +
+ // press_SEB * (c1o8 - c1o4 * mxoff + c1o4 * myoff + c1o4 * mzoff) +
+ // press_SWB * (c1o8 + c1o4 * mxoff + c1o4 * myoff + c1o4 * mzoff);
+ // drho = drho_NET * (c1o8 - c1o4 * xoff - c1o4 * yoff - c1o4 * zoff) +
+ // drho_NWT * (c1o8 + c1o4 * xoff - c1o4 * yoff - c1o4 * zoff) +
+ // drho_SET * (c1o8 - c1o4 * xoff + c1o4 * yoff - c1o4 * zoff) +
+ // drho_SWT * (c1o8 + c1o4 * xoff + c1o4 * yoff - c1o4 * zoff) +
+ // drho_NEB * (c1o8 - c1o4 * xoff - c1o4 * yoff + c1o4 * zoff) +
+ // drho_NWB * (c1o8 + c1o4 * xoff - c1o4 * yoff + c1o4 * zoff) +
+ // drho_SEB * (c1o8 - c1o4 * xoff + c1o4 * yoff + c1o4 * zoff) +
+ // drho_SWB * (c1o8 + c1o4 * xoff + c1o4 * yoff + c1o4 * zoff);
+ press = d0;
+ vvx = a0;
+ vvy = b0;
+ vvz = c0;
+
+ // mfaaa = drho;
+ // mfaaa = press + (ax+by+cz)/three; // 1/3 = 2/3*(1/op-1/2)
+ mfaaa = press; // if drho is interpolated directly
+
+ vx2 = vvx * vvx;
+ vy2 = vvy * vvy;
+ vz2 = vvz * vvz;
+ oMdrho = c1o1;
+ // oMdrho = one - mfaaa;
+
+ // two
+ // linear combinations
+ mxxPyyPzz = mfaaa;
+ // mxxMyy = -c2o3*(ax - by)*eps_new/o;
+ // mxxMzz = -c2o3*(ax - cz)*eps_new/o;
+
+ // mfabb = -c1o3 * (bz + cy)*eps_new/o;
+ // mfbab = -c1o3 * (az + cx)*eps_new/o;
+ // mfbba = -c1o3 * (ay + bx)*eps_new/o;
+ mxxMyy = -c2o3 * ((ax - by) + kxxMyyAverage) * eps_new / o * (c1o1 + press);
+ mxxMzz = -c2o3 * ((ax - cz) + kxxMzzAverage) * eps_new / o * (c1o1 + press);
+
+ mfabb = -c1o3 * ((bz + cy) + kyzAverage) * eps_new / o * (c1o1 + press);
+ mfbab = -c1o3 * ((az + cx) + kxzAverage) * eps_new / o * (c1o1 + press);
+ mfbba = -c1o3 * ((ay + bx) + kxyAverage) * eps_new / o * (c1o1 + press);
+
+ // linear combinations back
+ mfcaa = c1o3 * (mxxMyy + mxxMzz + mxxPyyPzz) * NeqOn;
+ mfaca = c1o3 * (-c2o1 * mxxMyy + mxxMzz + mxxPyyPzz) * NeqOn;
+ mfaac = c1o3 * (mxxMyy - c2o1 * mxxMzz + mxxPyyPzz) * NeqOn;
+
+ // 3.
+ // linear combinations
+ // residu = residutmp * (ayz + bxz + cxy );
+ // mfbbb = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
+ mfbbb = c0o1;
+
+ // residu = residutmp * (axy + two*bxx + two*bzz + cyz );
+ // residu = -(c1o9*(axy - 2*bxx - 2*bzz + cyz ));
+ // mxxyPyzz = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
+ mxxyPyzz = c0o1;
+
+ // residu = residutmp * (axy + two*bxx - two*bzz - cyz );
+ // residu = c1o9*(axy - 2*bxx + 2*bzz - cyz );
+ // mxxyMyzz = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
+ mxxyMyzz = c0o1;
+
+ // residu = residutmp * (axz + byz + two*cxx + two*cyy );
+ // residu = -(c1o9*(axz + byz - 2*cxx - 2*cyy ));
+ // mxxzPyyz = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
+ mxxzPyyz = c0o1;
+
+ // residu = residutmp * (axz - byz + two*cxx - two*cyy );
+ // residu = c1o9*(axz - byz - 2*cxx + 2*cyy );
+ // mxxzMyyz = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
+ mxxzMyyz = c0o1;
+
+ // residu = residutmp * (two*ayy + two*azz + bxy + cxz );
+ // residu = c1o9*(2*ayy + 2*azz - bxy - cxz );
+ // mxyyPxzz = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
+ mxyyPxzz = c0o1;
+
+ // residu = residutmp * (two*ayy - two*azz + bxy - cxz );
+ // residu = c1o9*(-2*ayy + 2*azz + bxy - cxz );
+ // mxyyMxzz = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
+ mxyyMxzz = c0o1;
+
+ // linear combinations back
+ mfcba = (mxxyMyzz + mxxyPyzz) * c1o2;
+ mfabc = (-mxxyMyzz + mxxyPyzz) * c1o2;
+ mfcab = (mxxzMyyz + mxxzPyyz) * c1o2;
+ mfacb = (-mxxzMyyz + mxxzPyyz) * c1o2;
+ mfbca = (mxyyMxzz + mxyyPxzz) * c1o2;
+ mfbac = (-mxyyMxzz + mxyyPxzz) * c1o2;
+
+ // 4.
+ mfacc = mfaaa * c1o9;
+ mfcac = mfacc;
+ mfcca = mfacc;
+ // 5.
+
+ // 6.
+ mfccc = mfaaa * c1o27;
+ ////////////////////////////////////////////////////////////////////////////////////
+ // back
+ ////////////////////////////////////////////////////////////////////////////////////
+ // mit 1, 0, 1/3, 0, 0, 0, 1/3, 0, 1/9 Konditionieren
+ ////////////////////////////////////////////////////////////////////////////////////
+ // Z - Dir
+ m0 = mfaac * c1o2 + mfaab * (vvz - c1o2) + (mfaaa + c1o1 * oMdrho) * (vz2 - vvz) * c1o2;
+ m1 = -mfaac - c2o1 * mfaab * vvz + mfaaa * (c1o1 - vz2) - c1o1 * oMdrho * vz2;
+ m2 = mfaac * c1o2 + mfaab * (vvz + c1o2) + (mfaaa + c1o1 * oMdrho) * (vz2 + vvz) * c1o2;
+ mfaaa = m0;
+ mfaab = m1;
+ mfaac = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ m0 = mfabc * c1o2 + mfabb * (vvz - c1o2) + mfaba * (vz2 - vvz) * c1o2;
+ m1 = -mfabc - c2o1 * mfabb * vvz + mfaba * (c1o1 - vz2);
+ m2 = mfabc * c1o2 + mfabb * (vvz + c1o2) + mfaba * (vz2 + vvz) * c1o2;
+ mfaba = m0;
+ mfabb = m1;
+ mfabc = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ m0 = mfacc * c1o2 + mfacb * (vvz - c1o2) + (mfaca + c1o3 * oMdrho) * (vz2 - vvz) * c1o2;
+ m1 = -mfacc - c2o1 * mfacb * vvz + mfaca * (c1o1 - vz2) - c1o3 * oMdrho * vz2;
+ m2 = mfacc * c1o2 + mfacb * (vvz + c1o2) + (mfaca + c1o3 * oMdrho) * (vz2 + vvz) * c1o2;
+ mfaca = m0;
+ mfacb = m1;
+ mfacc = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ ////////////////////////////////////////////////////////////////////////////////////
+ m0 = mfbac * c1o2 + mfbab * (vvz - c1o2) + mfbaa * (vz2 - vvz) * c1o2;
+ m1 = -mfbac - c2o1 * mfbab * vvz + mfbaa * (c1o1 - vz2);
+ m2 = mfbac * c1o2 + mfbab * (vvz + c1o2) + mfbaa * (vz2 + vvz) * c1o2;
+ mfbaa = m0;
+ mfbab = m1;
+ mfbac = m2;
+ /////////b//////////////////////////////////////////////////////////////////////////
+ m0 = mfbbc * c1o2 + mfbbb * (vvz - c1o2) + mfbba * (vz2 - vvz) * c1o2;
+ m1 = -mfbbc - c2o1 * mfbbb * vvz + mfbba * (c1o1 - vz2);
+ m2 = mfbbc * c1o2 + mfbbb * (vvz + c1o2) + mfbba * (vz2 + vvz) * c1o2;
+ mfbba = m0;
+ mfbbb = m1;
+ mfbbc = m2;
+ /////////b//////////////////////////////////////////////////////////////////////////
+ m0 = mfbcc * c1o2 + mfbcb * (vvz - c1o2) + mfbca * (vz2 - vvz) * c1o2;
+ m1 = -mfbcc - c2o1 * mfbcb * vvz + mfbca * (c1o1 - vz2);
+ m2 = mfbcc * c1o2 + mfbcb * (vvz + c1o2) + mfbca * (vz2 + vvz) * c1o2;
+ mfbca = m0;
+ mfbcb = m1;
+ mfbcc = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ ////////////////////////////////////////////////////////////////////////////////////
+ m0 = mfcac * c1o2 + mfcab * (vvz - c1o2) + (mfcaa + c1o3 * oMdrho) * (vz2 - vvz) * c1o2;
+ m1 = -mfcac - c2o1 * mfcab * vvz + mfcaa * (c1o1 - vz2) - c1o3 * oMdrho * vz2;
+ m2 = mfcac * c1o2 + mfcab * (vvz + c1o2) + (mfcaa + c1o3 * oMdrho) * (vz2 + vvz) * c1o2;
+ mfcaa = m0;
+ mfcab = m1;
+ mfcac = m2;
+ /////////c//////////////////////////////////////////////////////////////////////////
+ m0 = mfcbc * c1o2 + mfcbb * (vvz - c1o2) + mfcba * (vz2 - vvz) * c1o2;
+ m1 = -mfcbc - c2o1 * mfcbb * vvz + mfcba * (c1o1 - vz2);
+ m2 = mfcbc * c1o2 + mfcbb * (vvz + c1o2) + mfcba * (vz2 + vvz) * c1o2;
+ mfcba = m0;
+ mfcbb = m1;
+ mfcbc = m2;
+ /////////c//////////////////////////////////////////////////////////////////////////
+ m0 = mfccc * c1o2 + mfccb * (vvz - c1o2) + (mfcca + c1o9 * oMdrho) * (vz2 - vvz) * c1o2;
+ m1 = -mfccc - c2o1 * mfccb * vvz + mfcca * (c1o1 - vz2) - c1o9 * oMdrho * vz2;
+ m2 = mfccc * c1o2 + mfccb * (vvz + c1o2) + (mfcca + c1o9 * oMdrho) * (vz2 + vvz) * c1o2;
+ mfcca = m0;
+ mfccb = m1;
+ mfccc = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ ////////////////////////////////////////////////////////////////////////////////////
+ // mit 1/6, 2/3, 1/6, 0, 0, 0, 1/18, 2/9, 1/18 Konditionieren
+ ////////////////////////////////////////////////////////////////////////////////////
+ // Y - Dir
+ m0 = mfaca * c1o2 + mfaba * (vvy - c1o2) + (mfaaa + c1o6 * oMdrho) * (vy2 - vvy) * c1o2;
+ m1 = -mfaca - c2o1 * mfaba * vvy + mfaaa * (c1o1 - vy2) - c1o6 * oMdrho * vy2;
+ m2 = mfaca * c1o2 + mfaba * (vvy + c1o2) + (mfaaa + c1o6 * oMdrho) * (vy2 + vvy) * c1o2;
+ mfaaa = m0;
+ mfaba = m1;
+ mfaca = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ m0 = mfacb * c1o2 + mfabb * (vvy - c1o2) + (mfaab + c2o3 * oMdrho) * (vy2 - vvy) * c1o2;
+ m1 = -mfacb - c2o1 * mfabb * vvy + mfaab * (c1o1 - vy2) - c2o3 * oMdrho * vy2;
+ m2 = mfacb * c1o2 + mfabb * (vvy + c1o2) + (mfaab + c2o3 * oMdrho) * (vy2 + vvy) * c1o2;
+ mfaab = m0;
+ mfabb = m1;
+ mfacb = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ m0 = mfacc * c1o2 + mfabc * (vvy - c1o2) + (mfaac + c1o6 * oMdrho) * (vy2 - vvy) * c1o2;
+ m1 = -mfacc - c2o1 * mfabc * vvy + mfaac * (c1o1 - vy2) - c1o6 * oMdrho * vy2;
+ m2 = mfacc * c1o2 + mfabc * (vvy + c1o2) + (mfaac + c1o6 * oMdrho) * (vy2 + vvy) * c1o2;
+ mfaac = m0;
+ mfabc = m1;
+ mfacc = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ ////////////////////////////////////////////////////////////////////////////////////
+ m0 = mfbca * c1o2 + mfbba * (vvy - c1o2) + mfbaa * (vy2 - vvy) * c1o2;
+ m1 = -mfbca - c2o1 * mfbba * vvy + mfbaa * (c1o1 - vy2);
+ m2 = mfbca * c1o2 + mfbba * (vvy + c1o2) + mfbaa * (vy2 + vvy) * c1o2;
+ mfbaa = m0;
+ mfbba = m1;
+ mfbca = m2;
+ /////////b//////////////////////////////////////////////////////////////////////////
+ m0 = mfbcb * c1o2 + mfbbb * (vvy - c1o2) + mfbab * (vy2 - vvy) * c1o2;
+ m1 = -mfbcb - c2o1 * mfbbb * vvy + mfbab * (c1o1 - vy2);
+ m2 = mfbcb * c1o2 + mfbbb * (vvy + c1o2) + mfbab * (vy2 + vvy) * c1o2;
+ mfbab = m0;
+ mfbbb = m1;
+ mfbcb = m2;
+ /////////b//////////////////////////////////////////////////////////////////////////
+ m0 = mfbcc * c1o2 + mfbbc * (vvy - c1o2) + mfbac * (vy2 - vvy) * c1o2;
+ m1 = -mfbcc - c2o1 * mfbbc * vvy + mfbac * (c1o1 - vy2);
+ m2 = mfbcc * c1o2 + mfbbc * (vvy + c1o2) + mfbac * (vy2 + vvy) * c1o2;
+ mfbac = m0;
+ mfbbc = m1;
+ mfbcc = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ ////////////////////////////////////////////////////////////////////////////////////
+ m0 = mfcca * c1o2 + mfcba * (vvy - c1o2) + (mfcaa + c1o18 * oMdrho) * (vy2 - vvy) * c1o2;
+ m1 = -mfcca - c2o1 * mfcba * vvy + mfcaa * (c1o1 - vy2) - c1o18 * oMdrho * vy2;
+ m2 = mfcca * c1o2 + mfcba * (vvy + c1o2) + (mfcaa + c1o18 * oMdrho) * (vy2 + vvy) * c1o2;
+ mfcaa = m0;
+ mfcba = m1;
+ mfcca = m2;
+ /////////c//////////////////////////////////////////////////////////////////////////
+ m0 = mfccb * c1o2 + mfcbb * (vvy - c1o2) + (mfcab + c2o9 * oMdrho) * (vy2 - vvy) * c1o2;
+ m1 = -mfccb - c2o1 * mfcbb * vvy + mfcab * (c1o1 - vy2) - c2o9 * oMdrho * vy2;
+ m2 = mfccb * c1o2 + mfcbb * (vvy + c1o2) + (mfcab + c2o9 * oMdrho) * (vy2 + vvy) * c1o2;
+ mfcab = m0;
+ mfcbb = m1;
+ mfccb = m2;
+ /////////c//////////////////////////////////////////////////////////////////////////
+ m0 = mfccc * c1o2 + mfcbc * (vvy - c1o2) + (mfcac + c1o18 * oMdrho) * (vy2 - vvy) * c1o2;
+ m1 = -mfccc - c2o1 * mfcbc * vvy + mfcac * (c1o1 - vy2) - c1o18 * oMdrho * vy2;
+ m2 = mfccc * c1o2 + mfcbc * (vvy + c1o2) + (mfcac + c1o18 * oMdrho) * (vy2 + vvy) * c1o2;
+ mfcac = m0;
+ mfcbc = m1;
+ mfccc = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ ////////////////////////////////////////////////////////////////////////////////////
+ // mit 1/36, 1/9, 1/36, 1/9, 4/9, 1/9, 1/36, 1/9, 1/36 Konditionieren
+ ////////////////////////////////////////////////////////////////////////////////////
+ // X - Dir
+ m0 = mfcaa * c1o2 + mfbaa * (vvx - c1o2) + (mfaaa + c1o36 * oMdrho) * (vx2 - vvx) * c1o2;
+ m1 = -mfcaa - c2o1 * mfbaa * vvx + mfaaa * (c1o1 - vx2) - c1o36 * oMdrho * vx2;
+ m2 = mfcaa * c1o2 + mfbaa * (vvx + c1o2) + (mfaaa + c1o36 * oMdrho) * (vx2 + vvx) * c1o2;
+ mfaaa = m0;
+ mfbaa = m1;
+ mfcaa = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ m0 = mfcba * c1o2 + mfbba * (vvx - c1o2) + (mfaba + c1o9 * oMdrho) * (vx2 - vvx) * c1o2;
+ m1 = -mfcba - c2o1 * mfbba * vvx + mfaba * (c1o1 - vx2) - c1o9 * oMdrho * vx2;
+ m2 = mfcba * c1o2 + mfbba * (vvx + c1o2) + (mfaba + c1o9 * oMdrho) * (vx2 + vvx) * c1o2;
+ mfaba = m0;
+ mfbba = m1;
+ mfcba = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ m0 = mfcca * c1o2 + mfbca * (vvx - c1o2) + (mfaca + c1o36 * oMdrho) * (vx2 - vvx) * c1o2;
+ m1 = -mfcca - c2o1 * mfbca * vvx + mfaca * (c1o1 - vx2) - c1o36 * oMdrho * vx2;
+ m2 = mfcca * c1o2 + mfbca * (vvx + c1o2) + (mfaca + c1o36 * oMdrho) * (vx2 + vvx) * c1o2;
+ mfaca = m0;
+ mfbca = m1;
+ mfcca = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ ////////////////////////////////////////////////////////////////////////////////////
+ m0 = mfcab * c1o2 + mfbab * (vvx - c1o2) + (mfaab + c1o9 * oMdrho) * (vx2 - vvx) * c1o2;
+ m1 = -mfcab - c2o1 * mfbab * vvx + mfaab * (c1o1 - vx2) - c1o9 * oMdrho * vx2;
+ m2 = mfcab * c1o2 + mfbab * (vvx + c1o2) + (mfaab + c1o9 * oMdrho) * (vx2 + vvx) * c1o2;
+ mfaab = m0;
+ mfbab = m1;
+ mfcab = m2;
+ ///////////b////////////////////////////////////////////////////////////////////////
+ m0 = mfcbb * c1o2 + mfbbb * (vvx - c1o2) + (mfabb + c4o9 * oMdrho) * (vx2 - vvx) * c1o2;
+ m1 = -mfcbb - c2o1 * mfbbb * vvx + mfabb * (c1o1 - vx2) - c4o9 * oMdrho * vx2;
+ m2 = mfcbb * c1o2 + mfbbb * (vvx + c1o2) + (mfabb + c4o9 * oMdrho) * (vx2 + vvx) * c1o2;
+ mfabb = m0;
+ mfbbb = m1;
+ mfcbb = m2;
+ ///////////b////////////////////////////////////////////////////////////////////////
+ m0 = mfccb * c1o2 + mfbcb * (vvx - c1o2) + (mfacb + c1o9 * oMdrho) * (vx2 - vvx) * c1o2;
+ m1 = -mfccb - c2o1 * mfbcb * vvx + mfacb * (c1o1 - vx2) - c1o9 * oMdrho * vx2;
+ m2 = mfccb * c1o2 + mfbcb * (vvx + c1o2) + (mfacb + c1o9 * oMdrho) * (vx2 + vvx) * c1o2;
+ mfacb = m0;
+ mfbcb = m1;
+ mfccb = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+ ////////////////////////////////////////////////////////////////////////////////////
+ m0 = mfcac * c1o2 + mfbac * (vvx - c1o2) + (mfaac + c1o36 * oMdrho) * (vx2 - vvx) * c1o2;
+ m1 = -mfcac - c2o1 * mfbac * vvx + mfaac * (c1o1 - vx2) - c1o36 * oMdrho * vx2;
+ m2 = mfcac * c1o2 + mfbac * (vvx + c1o2) + (mfaac + c1o36 * oMdrho) * (vx2 + vvx) * c1o2;
+ mfaac = m0;
+ mfbac = m1;
+ mfcac = m2;
+ ///////////c////////////////////////////////////////////////////////////////////////
+ m0 = mfcbc * c1o2 + mfbbc * (vvx - c1o2) + (mfabc + c1o9 * oMdrho) * (vx2 - vvx) * c1o2;
+ m1 = -mfcbc - c2o1 * mfbbc * vvx + mfabc * (c1o1 - vx2) - c1o9 * oMdrho * vx2;
+ m2 = mfcbc * c1o2 + mfbbc * (vvx + c1o2) + (mfabc + c1o9 * oMdrho) * (vx2 + vvx) * c1o2;
+ mfabc = m0;
+ mfbbc = m1;
+ mfcbc = m2;
+ ///////////c////////////////////////////////////////////////////////////////////////
+ m0 = mfccc * c1o2 + mfbcc * (vvx - c1o2) + (mfacc + c1o36 * oMdrho) * (vx2 - vvx) * c1o2;
+ m1 = -mfccc - c2o1 * mfbcc * vvx + mfacc * (c1o1 - vx2) - c1o36 * oMdrho * vx2;
+ m2 = mfccc * c1o2 + mfbcc * (vvx + c1o2) + (mfacc + c1o36 * oMdrho) * (vx2 + vvx) * c1o2;
+ mfacc = m0;
+ mfbcc = m1;
+ mfccc = m2;
+ ////////////////////////////////////////////////////////////////////////////////////
+
+ ////////////////////////////////////////////////////////////////////////////////////
+ // index 0
+ kzero = posC[k];
+ kw = neighborCX[kzero];
+ ks = neighborCY[kzero];
+ kb = neighborCZ[kzero];
+ ksw = neighborCY[kw];
+ kbw = neighborCZ[kw];
+ kbs = neighborCZ[ks];
+ kbsw = neighborCZ[ksw];
+ ////////////////////////////////////////////////////////////////////////////////////
+
+ ////////////////////////////////////////////////////////////////////////////////////
+ feC[kzero] = mfcbb;
+ fwC[kw] = mfabb;
+ fnC[kzero] = mfbcb;
+ fsC[ks] = mfbab;
+ ftC[kzero] = mfbbc;
+ fbC[kb] = mfbba;
+ fneC[kzero] = mfccb;
+ fswC[ksw] = mfaab;
+ fseC[ks] = mfcab;
+ fnwC[kw] = mfacb;
+ fteC[kzero] = mfcbc;
+ fbwC[kbw] = mfaba;
+ fbeC[kb] = mfcba;
+ ftwC[kw] = mfabc;
+ ftnC[kzero] = mfbcc;
+ fbsC[kbs] = mfbaa;
+ fbnC[kb] = mfbca;
+ ftsC[ks] = mfbac;
+ fzeroC[kzero] = mfbbb;
+ ftneC[kzero] = mfccc;
+ ftseC[ks] = mfcac;
+ fbneC[kb] = mfcca;
+ fbseC[kbs] = mfcaa;
+ ftnwC[kw] = mfacc;
+ ftswC[ksw] = mfaac;
+ fbnwC[kbw] = mfaca;
+ fbswC[kbsw] = mfaaa;
+ ////////////////////////////////////////////////////////////////////////////////////
+ }
+}
+
extern "C" __global__ void scaleFC_RhoSq_comp_27(real* DC,
real* DF,
unsigned int* neighborCX,
@@ -9609,100 +11078,6 @@ extern "C" __global__ void scaleFC_RhoSq_comp_27(real* DC,
unsigned int nyF,
OffFC offFC)
{
- real *feF, *fwF, *fnF, *fsF, *ftF, *fbF, *fneF, *fswF, *fseF, *fnwF, *fteF, *fbwF, *fbeF, *ftwF, *ftnF, *fbsF, *fbnF, *ftsF, *fzeroF,
- *ftneF, *ftswF, *ftseF, *ftnwF, *fbneF, *fbswF, *fbseF, *fbnwF;
-
- feF = &DF[dirE *size_MatF];
- fwF = &DF[dirW *size_MatF];
- fnF = &DF[dirN *size_MatF];
- fsF = &DF[dirS *size_MatF];
- ftF = &DF[dirT *size_MatF];
- fbF = &DF[dirB *size_MatF];
- fneF = &DF[dirNE *size_MatF];
- fswF = &DF[dirSW *size_MatF];
- fseF = &DF[dirSE *size_MatF];
- fnwF = &DF[dirNW *size_MatF];
- fteF = &DF[dirTE *size_MatF];
- fbwF = &DF[dirBW *size_MatF];
- fbeF = &DF[dirBE *size_MatF];
- ftwF = &DF[dirTW *size_MatF];
- ftnF = &DF[dirTN *size_MatF];
- fbsF = &DF[dirBS *size_MatF];
- fbnF = &DF[dirBN *size_MatF];
- ftsF = &DF[dirTS *size_MatF];
- fzeroF = &DF[dirZERO*size_MatF];
- ftneF = &DF[dirTNE *size_MatF];
- ftswF = &DF[dirTSW *size_MatF];
- ftseF = &DF[dirTSE *size_MatF];
- ftnwF = &DF[dirTNW *size_MatF];
- fbneF = &DF[dirBNE *size_MatF];
- fbswF = &DF[dirBSW *size_MatF];
- fbseF = &DF[dirBSE *size_MatF];
- fbnwF = &DF[dirBNW *size_MatF];
-
- real *feC, *fwC, *fnC, *fsC, *ftC, *fbC, *fneC, *fswC, *fseC, *fnwC, *fteC, *fbwC, *fbeC, *ftwC, *ftnC, *fbsC, *fbnC, *ftsC, *fzeroC,
- *ftneC, *ftswC, *ftseC, *ftnwC, *fbneC, *fbswC, *fbseC, *fbnwC;
-
- if (evenOrOdd==true)
- {
- feC = &DC[dirE *size_MatC];
- fwC = &DC[dirW *size_MatC];
- fnC = &DC[dirN *size_MatC];
- fsC = &DC[dirS *size_MatC];
- ftC = &DC[dirT *size_MatC];
- fbC = &DC[dirB *size_MatC];
- fneC = &DC[dirNE *size_MatC];
- fswC = &DC[dirSW *size_MatC];
- fseC = &DC[dirSE *size_MatC];
- fnwC = &DC[dirNW *size_MatC];
- fteC = &DC[dirTE *size_MatC];
- fbwC = &DC[dirBW *size_MatC];
- fbeC = &DC[dirBE *size_MatC];
- ftwC = &DC[dirTW *size_MatC];
- ftnC = &DC[dirTN *size_MatC];
- fbsC = &DC[dirBS *size_MatC];
- fbnC = &DC[dirBN *size_MatC];
- ftsC = &DC[dirTS *size_MatC];
- fzeroC = &DC[dirZERO*size_MatC];
- ftneC = &DC[dirTNE *size_MatC];
- ftswC = &DC[dirTSW *size_MatC];
- ftseC = &DC[dirTSE *size_MatC];
- ftnwC = &DC[dirTNW *size_MatC];
- fbneC = &DC[dirBNE *size_MatC];
- fbswC = &DC[dirBSW *size_MatC];
- fbseC = &DC[dirBSE *size_MatC];
- fbnwC = &DC[dirBNW *size_MatC];
- }
- else
- {
- fwC = &DC[dirE *size_MatC];
- feC = &DC[dirW *size_MatC];
- fsC = &DC[dirN *size_MatC];
- fnC = &DC[dirS *size_MatC];
- fbC = &DC[dirT *size_MatC];
- ftC = &DC[dirB *size_MatC];
- fswC = &DC[dirNE *size_MatC];
- fneC = &DC[dirSW *size_MatC];
- fnwC = &DC[dirSE *size_MatC];
- fseC = &DC[dirNW *size_MatC];
- fbwC = &DC[dirTE *size_MatC];
- fteC = &DC[dirBW *size_MatC];
- ftwC = &DC[dirBE *size_MatC];
- fbeC = &DC[dirTW *size_MatC];
- fbsC = &DC[dirTN *size_MatC];
- ftnC = &DC[dirBS *size_MatC];
- ftsC = &DC[dirBN *size_MatC];
- fbnC = &DC[dirTS *size_MatC];
- fzeroC = &DC[dirZERO*size_MatC];
- fbswC = &DC[dirTNE *size_MatC];
- fbneC = &DC[dirTSW *size_MatC];
- fbnwC = &DC[dirTSE *size_MatC];
- fbseC = &DC[dirTNW *size_MatC];
- ftswC = &DC[dirBNE *size_MatC];
- ftneC = &DC[dirBSW *size_MatC];
- ftnwC = &DC[dirBSE *size_MatC];
- ftseC = &DC[dirBNW *size_MatC];
- }
////////////////////////////////////////////////////////////////////////////////
const unsigned ix = threadIdx.x; // Globaler x-Index
const unsigned iy = blockIdx.x; // Globaler y-Index
@@ -9714,1178 +11089,36 @@ extern "C" __global__ void scaleFC_RhoSq_comp_27(real* DC,
const unsigned k = nx*(ny*iz + iy) + ix;
//////////////////////////////////////////////////////////////////////////
- ////////////////////////////////////////////////////////////////////////////////
- real eps_new = c2o1;
- real omegaS = omFine;//-omFine;
- real o = omCoarse;//-omCoarse;
- //real op = one;
- //real cu_sq;
-
- real xoff, yoff, zoff;
- real xoff_sq, yoff_sq, zoff_sq;
-
- real press;//,drho,vx1,vx2,vx3;
- real /*press_SWT,*/drho_SWT,vx1_SWT,vx2_SWT,vx3_SWT;
- real /*press_NWT,*/drho_NWT,vx1_NWT,vx2_NWT,vx3_NWT;
- real /*press_NET,*/drho_NET,vx1_NET,vx2_NET,vx3_NET;
- real /*press_SET,*/drho_SET,vx1_SET,vx2_SET,vx3_SET;
- real /*press_SWB,*/drho_SWB,vx1_SWB,vx2_SWB,vx3_SWB;
- real /*press_NWB,*/drho_NWB,vx1_NWB,vx2_NWB,vx3_NWB;
- real /*press_NEB,*/drho_NEB,vx1_NEB,vx2_NEB,vx3_NEB;
- real /*press_SEB,*/drho_SEB,vx1_SEB,vx2_SEB,vx3_SEB;
- real f_E,f_W,f_N,f_S,f_T,f_B,f_NE,f_SW,f_SE,f_NW,f_TE,f_BW,f_BE,f_TW,f_TN,f_BS,f_BN,f_TS,f_ZERO,f_TNE, f_TSW, f_TSE, f_TNW, f_BNE, f_BSW, f_BSE, f_BNW;
- //real feq_E,feq_W,feq_N,feq_S,feq_T,feq_B,feq_NE,feq_SW,feq_SE,feq_NW,feq_TE,feq_BW,feq_BE,feq_TW,feq_TN,feq_BS,feq_BN,feq_TS,feq_ZERO,feq_TNE, feq_TSW, feq_TSE, feq_TNW, feq_BNE, feq_BSW, feq_BSE, feq_BNW;
- real kxyFromfcNEQ_SWT, kyzFromfcNEQ_SWT, kxzFromfcNEQ_SWT, kxxMyyFromfcNEQ_SWT, kxxMzzFromfcNEQ_SWT;
- real kxyFromfcNEQ_NWT, kyzFromfcNEQ_NWT, kxzFromfcNEQ_NWT, kxxMyyFromfcNEQ_NWT, kxxMzzFromfcNEQ_NWT;
- real kxyFromfcNEQ_NET, kyzFromfcNEQ_NET, kxzFromfcNEQ_NET, kxxMyyFromfcNEQ_NET, kxxMzzFromfcNEQ_NET;
- real kxyFromfcNEQ_SET, kyzFromfcNEQ_SET, kxzFromfcNEQ_SET, kxxMyyFromfcNEQ_SET, kxxMzzFromfcNEQ_SET;
- real kxyFromfcNEQ_SWB, kyzFromfcNEQ_SWB, kxzFromfcNEQ_SWB, kxxMyyFromfcNEQ_SWB, kxxMzzFromfcNEQ_SWB;
- real kxyFromfcNEQ_NWB, kyzFromfcNEQ_NWB, kxzFromfcNEQ_NWB, kxxMyyFromfcNEQ_NWB, kxxMzzFromfcNEQ_NWB;
- real kxyFromfcNEQ_NEB, kyzFromfcNEQ_NEB, kxzFromfcNEQ_NEB, kxxMyyFromfcNEQ_NEB, kxxMzzFromfcNEQ_NEB;
- real kxyFromfcNEQ_SEB, kyzFromfcNEQ_SEB, kxzFromfcNEQ_SEB, kxxMyyFromfcNEQ_SEB, kxxMzzFromfcNEQ_SEB;
- real a0, ax, ay, az, axx, ayy, azz, axy, axz, ayz, b0, bx, by, bz, bxx, byy, bzz, bxy, bxz, byz, c0, cx, cy, cz, cxx, cyy, czz, cxy, cxz, cyz/*, axyz, bxyz, cxyz*/;
- real d0, dx, dy, dz, dxy, dxz, dyz/*, dxyz*/;
-
- if(k<kFC)
- {
- //////////////////////////////////////////////////////////////////////////
- xoff = offFC.xOffFC[k];
- yoff = offFC.yOffFC[k];
- zoff = offFC.zOffFC[k];
- xoff_sq = xoff * xoff;
- yoff_sq = yoff * yoff;
- zoff_sq = zoff * zoff;
- //////////////////////////////////////////////////////////////////////////
- //SWB//
- //////////////////////////////////////////////////////////////////////////
- //index 0
- unsigned int k0zero= posFSWB[k];
- unsigned int k0w = neighborFX[k0zero];
- unsigned int k0s = neighborFY[k0zero];
- unsigned int k0b = neighborFZ[k0zero];
- unsigned int k0sw = neighborFY[k0w];
- unsigned int k0bw = neighborFZ[k0w];
- unsigned int k0bs = neighborFZ[k0s];
- unsigned int k0bsw = neighborFZ[k0sw];
- //////////////////////////////////////////////////////////////////////////
- //index
- unsigned int kzero= k0zero;
- unsigned int kw = k0w;
- unsigned int ks = k0s;
- unsigned int kb = k0b;
- unsigned int ksw = k0sw;
- unsigned int kbw = k0bw;
- unsigned int kbs = k0bs;
- unsigned int kbsw = k0bsw;
- ////////////////////////////////////////////////////////////////////////////////
- f_E = feF[kzero];
- f_W = fwF[kw];
- f_N = fnF[kzero];
- f_S = fsF[ks];
- f_T = ftF[kzero];
- f_B = fbF[kb];
- f_NE = fneF[kzero];
- f_SW = fswF[ksw];
- f_SE = fseF[ks];
- f_NW = fnwF[kw];
- f_TE = fteF[kzero];
- f_BW = fbwF[kbw];
- f_BE = fbeF[kb];
- f_TW = ftwF[kw];
- f_TN = ftnF[kzero];
- f_BS = fbsF[kbs];
- f_BN = fbnF[kb];
- f_TS = ftsF[ks];
- f_ZERO = fzeroF[kzero];
- f_TNE = ftneF[kzero];
- f_TSW = ftswF[ksw];
- f_TSE = ftseF[ks];
- f_TNW = ftnwF[kw];
- f_BNE = fbneF[kb];
- f_BSW = fbswF[kbsw];
- f_BSE = fbseF[kbs];
- f_BNW = fbnwF[kbw];
-
- drho_SWB = f_E+f_W+f_N+f_S+f_T+f_B+f_NE+f_SW+f_SE+f_NW+f_TE+f_BW+f_BE+f_TW+f_TN+f_BS+f_BN+f_TS+f_ZERO+f_TNE+f_TSW+f_TSE+f_TNW+f_BNE+f_BSW+f_BSE+f_BNW;
- vx1_SWB = (((f_TNE-f_BSW)+(f_TSE-f_BNW)+(f_BNE-f_TSW)+(f_BSE-f_TNW)) + (((f_NE-f_SW)+(f_TE-f_BW))+((f_SE-f_NW)+(f_BE-f_TW))) + (f_E-f_W))/(c1o1 + drho_SWB);
- vx2_SWB = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_BNE-f_TSW)+(f_BNW-f_TSE)) + (((f_NE-f_SW)+(f_TN-f_BS))+((f_BN-f_TS)+(f_NW-f_SE))) + (f_N-f_S))/(c1o1 + drho_SWB);
- vx3_SWB = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_TSE-f_BNW)+(f_TSW-f_BNE)) + (((f_TE-f_BW)+(f_TN-f_BS))+((f_TW-f_BE)+(f_TS-f_BN))) + (f_T-f_B))/(c1o1 + drho_SWB);
-
- kxyFromfcNEQ_SWB = -c3o1*omegaS*((f_SW+f_BSW+f_TSW-f_NW-f_BNW-f_TNW-f_SE-f_BSE-f_TSE+f_NE+f_BNE+f_TNE ) / (c1o1 + drho_SWB) - ((vx1_SWB*vx2_SWB)));
- kyzFromfcNEQ_SWB = -c3o1*omegaS*((f_BS+f_BSE+f_BSW-f_TS-f_TSE-f_TSW-f_BN-f_BNE-f_BNW+f_TN+f_TNE+f_TNW ) / (c1o1 + drho_SWB) - ((vx2_SWB*vx3_SWB)));
- kxzFromfcNEQ_SWB = -c3o1*omegaS*((f_BW+f_BSW+f_BNW-f_TW-f_TSW-f_TNW-f_BE-f_BSE-f_BNE+f_TE+f_TSE+f_TNE ) / (c1o1 + drho_SWB) - ((vx1_SWB*vx3_SWB)));
- kxxMyyFromfcNEQ_SWB = -c3o2*omegaS *((f_BW+f_W+f_TW-f_BS-f_S-f_TS-f_BN-f_N-f_TN+f_BE+f_E+f_TE ) / (c1o1 + drho_SWB) - ((vx1_SWB*vx1_SWB-vx2_SWB*vx2_SWB)));
- kxxMzzFromfcNEQ_SWB = -c3o2*omegaS *((f_SW+f_W+f_NW-f_BS-f_TS-f_B-f_T-f_BN-f_TN+f_SE+f_E+f_NE ) / (c1o1 + drho_SWB) - ((vx1_SWB*vx1_SWB-vx3_SWB*vx3_SWB)));
-
-
- //////////////////////////////////////////////////////////////////////////
- //SWT//
- //////////////////////////////////////////////////////////////////////////
- //index
- kzero= kb;
- kw = kbw;
- ks = kbs;
- kb = neighborFZ[kb];
- ksw = kbsw;
- kbw = neighborFZ[kbw];
- kbs = neighborFZ[kbs];
- kbsw = neighborFZ[kbsw];
- ////////////////////////////////////////////////////////////////////////////////
- f_E = feF[kzero];
- f_W = fwF[kw];
- f_N = fnF[kzero];
- f_S = fsF[ks];
- f_T = ftF[kzero];
- f_B = fbF[kb];
- f_NE = fneF[kzero];
- f_SW = fswF[ksw];
- f_SE = fseF[ks];
- f_NW = fnwF[kw];
- f_TE = fteF[kzero];
- f_BW = fbwF[kbw];
- f_BE = fbeF[kb];
- f_TW = ftwF[kw];
- f_TN = ftnF[kzero];
- f_BS = fbsF[kbs];
- f_BN = fbnF[kb];
- f_TS = ftsF[ks];
- f_ZERO = fzeroF[kzero];
- f_TNE = ftneF[kzero];
- f_TSW = ftswF[ksw];
- f_TSE = ftseF[ks];
- f_TNW = ftnwF[kw];
- f_BNE = fbneF[kb];
- f_BSW = fbswF[kbsw];
- f_BSE = fbseF[kbs];
- f_BNW = fbnwF[kbw];
-
- drho_SWT = f_E+f_W+f_N+f_S+f_T+f_B+f_NE+f_SW+f_SE+f_NW+f_TE+f_BW+f_BE+f_TW+f_TN+f_BS+f_BN+f_TS+f_ZERO+f_TNE+f_TSW+f_TSE+f_TNW+f_BNE+f_BSW+f_BSE+f_BNW;
- vx1_SWT = (((f_TNE-f_BSW)+(f_TSE-f_BNW)+(f_BNE-f_TSW)+(f_BSE-f_TNW)) + (((f_NE-f_SW)+(f_TE-f_BW))+((f_SE-f_NW)+(f_BE-f_TW))) + (f_E-f_W))/(c1o1 + drho_SWT);
- vx2_SWT = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_BNE-f_TSW)+(f_BNW-f_TSE)) + (((f_NE-f_SW)+(f_TN-f_BS))+((f_BN-f_TS)+(f_NW-f_SE))) + (f_N-f_S))/(c1o1 + drho_SWT);
- vx3_SWT = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_TSE-f_BNW)+(f_TSW-f_BNE)) + (((f_TE-f_BW)+(f_TN-f_BS))+((f_TW-f_BE)+(f_TS-f_BN))) + (f_T-f_B))/(c1o1 + drho_SWT);
-
- kxyFromfcNEQ_SWT = -c3o1*omegaS*((f_SW+f_BSW+f_TSW-f_NW-f_BNW-f_TNW-f_SE-f_BSE-f_TSE+f_NE+f_BNE+f_TNE ) / (c1o1 + drho_SWT) - ((vx1_SWT*vx2_SWT)));
- kyzFromfcNEQ_SWT = -c3o1*omegaS*((f_BS+f_BSE+f_BSW-f_TS-f_TSE-f_TSW-f_BN-f_BNE-f_BNW+f_TN+f_TNE+f_TNW ) / (c1o1 + drho_SWT) - ((vx2_SWT*vx3_SWT)));
- kxzFromfcNEQ_SWT = -c3o1*omegaS*((f_BW+f_BSW+f_BNW-f_TW-f_TSW-f_TNW-f_BE-f_BSE-f_BNE+f_TE+f_TSE+f_TNE ) / (c1o1 + drho_SWT) - ((vx1_SWT*vx3_SWT)));
- kxxMyyFromfcNEQ_SWT = -c3o2*omegaS *((f_BW+f_W+f_TW-f_BS-f_S-f_TS-f_BN-f_N-f_TN+f_BE+f_E+f_TE ) / (c1o1 + drho_SWT) - ((vx1_SWT*vx1_SWT-vx2_SWT*vx2_SWT)));
- kxxMzzFromfcNEQ_SWT = -c3o2*omegaS *((f_SW+f_W+f_NW-f_BS-f_TS-f_B-f_T-f_BN-f_TN+f_SE+f_E+f_NE ) / (c1o1 + drho_SWT) - ((vx1_SWT*vx1_SWT-vx3_SWT*vx3_SWT)));
+ scaleFC_RhoSq_comp_27_Calculation(DC, DF, neighborCX, neighborCY, neighborCZ, neighborFX, neighborFY, neighborFZ,
+ size_MatC, size_MatF, evenOrOdd, posC, posFSWB, kFC, omCoarse, omFine, nu, nxC,
+ nyC, nxF, nyF, offFC, k);
+}
- //////////////////////////////////////////////////////////////////////////
- //SET//
- //////////////////////////////////////////////////////////////////////////
- //index
- kzero= kw;
- kw = neighborFX[kw];
- ks = ksw;
- kb = kbw;
- ksw = neighborFX[ksw];
- kbw = neighborFX[kbw];
- kbs = kbsw;
- kbsw = neighborFX[kbsw];
- ////////////////////////////////////////////////////////////////////////////////
- f_E = feF[kzero];
- f_W = fwF[kw];
- f_N = fnF[kzero];
- f_S = fsF[ks];
- f_T = ftF[kzero];
- f_B = fbF[kb];
- f_NE = fneF[kzero];
- f_SW = fswF[ksw];
- f_SE = fseF[ks];
- f_NW = fnwF[kw];
- f_TE = fteF[kzero];
- f_BW = fbwF[kbw];
- f_BE = fbeF[kb];
- f_TW = ftwF[kw];
- f_TN = ftnF[kzero];
- f_BS = fbsF[kbs];
- f_BN = fbnF[kb];
- f_TS = ftsF[ks];
- f_ZERO = fzeroF[kzero];
- f_TNE = ftneF[kzero];
- f_TSW = ftswF[ksw];
- f_TSE = ftseF[ks];
- f_TNW = ftnwF[kw];
- f_BNE = fbneF[kb];
- f_BSW = fbswF[kbsw];
- f_BSE = fbseF[kbs];
- f_BNW = fbnwF[kbw];
-
- drho_SET = f_E+f_W+f_N+f_S+f_T+f_B+f_NE+f_SW+f_SE+f_NW+f_TE+f_BW+f_BE+f_TW+f_TN+f_BS+f_BN+f_TS+f_ZERO+f_TNE+f_TSW+f_TSE+f_TNW+f_BNE+f_BSW+f_BSE+f_BNW;
- vx1_SET = (((f_TNE-f_BSW)+(f_TSE-f_BNW)+(f_BNE-f_TSW)+(f_BSE-f_TNW)) + (((f_NE-f_SW)+(f_TE-f_BW))+((f_SE-f_NW)+(f_BE-f_TW))) + (f_E-f_W))/(c1o1 + drho_SET);
- vx2_SET = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_BNE-f_TSW)+(f_BNW-f_TSE)) + (((f_NE-f_SW)+(f_TN-f_BS))+((f_BN-f_TS)+(f_NW-f_SE))) + (f_N-f_S))/(c1o1 + drho_SET);
- vx3_SET = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_TSE-f_BNW)+(f_TSW-f_BNE)) + (((f_TE-f_BW)+(f_TN-f_BS))+((f_TW-f_BE)+(f_TS-f_BN))) + (f_T-f_B))/(c1o1 + drho_SET);
-
- kxyFromfcNEQ_SET = -c3o1*omegaS*((f_SW+f_BSW+f_TSW-f_NW-f_BNW-f_TNW-f_SE-f_BSE-f_TSE+f_NE+f_BNE+f_TNE ) / (c1o1 + drho_SET) - ((vx1_SET*vx2_SET)));
- kyzFromfcNEQ_SET = -c3o1*omegaS*((f_BS+f_BSE+f_BSW-f_TS-f_TSE-f_TSW-f_BN-f_BNE-f_BNW+f_TN+f_TNE+f_TNW ) / (c1o1 + drho_SET) - ((vx2_SET*vx3_SET)));
- kxzFromfcNEQ_SET = -c3o1*omegaS*((f_BW+f_BSW+f_BNW-f_TW-f_TSW-f_TNW-f_BE-f_BSE-f_BNE+f_TE+f_TSE+f_TNE ) / (c1o1 + drho_SET) - ((vx1_SET*vx3_SET)));
- kxxMyyFromfcNEQ_SET = -c3o2*omegaS *((f_BW+f_W+f_TW-f_BS-f_S-f_TS-f_BN-f_N-f_TN+f_BE+f_E+f_TE ) / (c1o1 + drho_SET) - ((vx1_SET*vx1_SET-vx2_SET*vx2_SET)));
- kxxMzzFromfcNEQ_SET = -c3o2*omegaS *((f_SW+f_W+f_NW-f_BS-f_TS-f_B-f_T-f_BN-f_TN+f_SE+f_E+f_NE ) / (c1o1 + drho_SET) - ((vx1_SET*vx1_SET-vx3_SET*vx3_SET)));
-
- //////////////////////////////////////////////////////////////////////////
- //SEB//
- //////////////////////////////////////////////////////////////////////////
- //index
- kb = kzero;
- kbw = kw;
- kbs = ks;
- kbsw = ksw;
- kzero= k0w;
- kw = neighborFX[k0w];
- ks = k0sw;
- ksw = neighborFX[k0sw];
- ////////////////////////////////////////////////////////////////////////////////
- f_E = feF[kzero];
- f_W = fwF[kw];
- f_N = fnF[kzero];
- f_S = fsF[ks];
- f_T = ftF[kzero];
- f_B = fbF[kb];
- f_NE = fneF[kzero];
- f_SW = fswF[ksw];
- f_SE = fseF[ks];
- f_NW = fnwF[kw];
- f_TE = fteF[kzero];
- f_BW = fbwF[kbw];
- f_BE = fbeF[kb];
- f_TW = ftwF[kw];
- f_TN = ftnF[kzero];
- f_BS = fbsF[kbs];
- f_BN = fbnF[kb];
- f_TS = ftsF[ks];
- f_ZERO = fzeroF[kzero];
- f_TNE = ftneF[kzero];
- f_TSW = ftswF[ksw];
- f_TSE = ftseF[ks];
- f_TNW = ftnwF[kw];
- f_BNE = fbneF[kb];
- f_BSW = fbswF[kbsw];
- f_BSE = fbseF[kbs];
- f_BNW = fbnwF[kbw];
-
- drho_SEB = f_E+f_W+f_N+f_S+f_T+f_B+f_NE+f_SW+f_SE+f_NW+f_TE+f_BW+f_BE+f_TW+f_TN+f_BS+f_BN+f_TS+f_ZERO+f_TNE+f_TSW+f_TSE+f_TNW+f_BNE+f_BSW+f_BSE+f_BNW;
- vx1_SEB = (((f_TNE-f_BSW)+(f_TSE-f_BNW)+(f_BNE-f_TSW)+(f_BSE-f_TNW)) + (((f_NE-f_SW)+(f_TE-f_BW))+((f_SE-f_NW)+(f_BE-f_TW))) + (f_E-f_W))/(c1o1 + drho_SEB);
- vx2_SEB = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_BNE-f_TSW)+(f_BNW-f_TSE)) + (((f_NE-f_SW)+(f_TN-f_BS))+((f_BN-f_TS)+(f_NW-f_SE))) + (f_N-f_S))/(c1o1 + drho_SEB);
- vx3_SEB = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_TSE-f_BNW)+(f_TSW-f_BNE)) + (((f_TE-f_BW)+(f_TN-f_BS))+((f_TW-f_BE)+(f_TS-f_BN))) + (f_T-f_B))/(c1o1 + drho_SEB);
-
- kxyFromfcNEQ_SEB = -c3o1*omegaS*((f_SW+f_BSW+f_TSW-f_NW-f_BNW-f_TNW-f_SE-f_BSE-f_TSE+f_NE+f_BNE+f_TNE ) / (c1o1 + drho_SEB) - ((vx1_SEB*vx2_SEB)));
- kyzFromfcNEQ_SEB = -c3o1*omegaS*((f_BS+f_BSE+f_BSW-f_TS-f_TSE-f_TSW-f_BN-f_BNE-f_BNW+f_TN+f_TNE+f_TNW ) / (c1o1 + drho_SEB) - ((vx2_SEB*vx3_SEB)));
- kxzFromfcNEQ_SEB = -c3o1*omegaS*((f_BW+f_BSW+f_BNW-f_TW-f_TSW-f_TNW-f_BE-f_BSE-f_BNE+f_TE+f_TSE+f_TNE ) / (c1o1 + drho_SEB) - ((vx1_SEB*vx3_SEB)));
- kxxMyyFromfcNEQ_SEB = -c3o2*omegaS *((f_BW+f_W+f_TW-f_BS-f_S-f_TS-f_BN-f_N-f_TN+f_BE+f_E+f_TE ) / (c1o1 + drho_SEB) - ((vx1_SEB*vx1_SEB-vx2_SEB*vx2_SEB)));
- kxxMzzFromfcNEQ_SEB = -c3o2*omegaS *((f_SW+f_W+f_NW-f_BS-f_TS-f_B-f_T-f_BN-f_TN+f_SE+f_E+f_NE ) / (c1o1 + drho_SEB) - ((vx1_SEB*vx1_SEB-vx3_SEB*vx3_SEB)));
-
- //////////////////////////////////////////////////////////////////////////
- //NWB//
- //////////////////////////////////////////////////////////////////////////
- //index 0
- k0zero= k0s;
- k0w = k0sw;
- k0s = neighborFY[k0s];
- k0b = k0bs;
- k0sw = neighborFY[k0sw];
- k0bw = k0bsw;
- k0bs = neighborFY[k0bs];
- k0bsw = neighborFY[k0bsw];
- //////////////////////////////////////////////////////////////////////////
- //index
- kzero= k0zero;
- kw = k0w;
- ks = k0s;
- kb = k0b;
- ksw = k0sw;
- kbw = k0bw;
- kbs = k0bs;
- kbsw = k0bsw;
- ////////////////////////////////////////////////////////////////////////////////
- f_E = feF[kzero];
- f_W = fwF[kw];
- f_N = fnF[kzero];
- f_S = fsF[ks];
- f_T = ftF[kzero];
- f_B = fbF[kb];
- f_NE = fneF[kzero];
- f_SW = fswF[ksw];
- f_SE = fseF[ks];
- f_NW = fnwF[kw];
- f_TE = fteF[kzero];
- f_BW = fbwF[kbw];
- f_BE = fbeF[kb];
- f_TW = ftwF[kw];
- f_TN = ftnF[kzero];
- f_BS = fbsF[kbs];
- f_BN = fbnF[kb];
- f_TS = ftsF[ks];
- f_ZERO = fzeroF[kzero];
- f_TNE = ftneF[kzero];
- f_TSW = ftswF[ksw];
- f_TSE = ftseF[ks];
- f_TNW = ftnwF[kw];
- f_BNE = fbneF[kb];
- f_BSW = fbswF[kbsw];
- f_BSE = fbseF[kbs];
- f_BNW = fbnwF[kbw];
-
- drho_NWB = f_E+f_W+f_N+f_S+f_T+f_B+f_NE+f_SW+f_SE+f_NW+f_TE+f_BW+f_BE+f_TW+f_TN+f_BS+f_BN+f_TS+f_ZERO+f_TNE+f_TSW+f_TSE+f_TNW+f_BNE+f_BSW+f_BSE+f_BNW;
- vx1_NWB = (((f_TNE-f_BSW)+(f_TSE-f_BNW)+(f_BNE-f_TSW)+(f_BSE-f_TNW)) + (((f_NE-f_SW)+(f_TE-f_BW))+((f_SE-f_NW)+(f_BE-f_TW))) + (f_E-f_W))/(c1o1 + drho_NWB);
- vx2_NWB = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_BNE-f_TSW)+(f_BNW-f_TSE)) + (((f_NE-f_SW)+(f_TN-f_BS))+((f_BN-f_TS)+(f_NW-f_SE))) + (f_N-f_S))/(c1o1 + drho_NWB);
- vx3_NWB = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_TSE-f_BNW)+(f_TSW-f_BNE)) + (((f_TE-f_BW)+(f_TN-f_BS))+((f_TW-f_BE)+(f_TS-f_BN))) + (f_T-f_B))/(c1o1 + drho_NWB);
-
- kxyFromfcNEQ_NWB = -c3o1*omegaS*((f_SW+f_BSW+f_TSW-f_NW-f_BNW-f_TNW-f_SE-f_BSE-f_TSE+f_NE+f_BNE+f_TNE ) / (c1o1 + drho_NWB) - ((vx1_NWB*vx2_NWB)));
- kyzFromfcNEQ_NWB = -c3o1*omegaS*((f_BS+f_BSE+f_BSW-f_TS-f_TSE-f_TSW-f_BN-f_BNE-f_BNW+f_TN+f_TNE+f_TNW ) / (c1o1 + drho_NWB) - ((vx2_NWB*vx3_NWB)));
- kxzFromfcNEQ_NWB = -c3o1*omegaS*((f_BW+f_BSW+f_BNW-f_TW-f_TSW-f_TNW-f_BE-f_BSE-f_BNE+f_TE+f_TSE+f_TNE ) / (c1o1 + drho_NWB) - ((vx1_NWB*vx3_NWB)));
- kxxMyyFromfcNEQ_NWB = -c3o2*omegaS *((f_BW+f_W+f_TW-f_BS-f_S-f_TS-f_BN-f_N-f_TN+f_BE+f_E+f_TE ) / (c1o1 + drho_NWB) - ((vx1_NWB*vx1_NWB-vx2_NWB*vx2_NWB)));
- kxxMzzFromfcNEQ_NWB = -c3o2*omegaS *((f_SW+f_W+f_NW-f_BS-f_TS-f_B-f_T-f_BN-f_TN+f_SE+f_E+f_NE ) / (c1o1 + drho_NWB) - ((vx1_NWB*vx1_NWB-vx3_NWB*vx3_NWB)));
-
- //////////////////////////////////////////////////////////////////////////
- //NWT//
- //////////////////////////////////////////////////////////////////////////
- //index
- kzero= kb;
- kw = kbw;
- ks = kbs;
- kb = neighborFZ[kb];
- ksw = kbsw;
- kbw = neighborFZ[kbw];
- kbs = neighborFZ[kbs];
- kbsw = neighborFZ[kbsw];
- ////////////////////////////////////////////////////////////////////////////////
- f_E = feF[kzero];
- f_W = fwF[kw];
- f_N = fnF[kzero];
- f_S = fsF[ks];
- f_T = ftF[kzero];
- f_B = fbF[kb];
- f_NE = fneF[kzero];
- f_SW = fswF[ksw];
- f_SE = fseF[ks];
- f_NW = fnwF[kw];
- f_TE = fteF[kzero];
- f_BW = fbwF[kbw];
- f_BE = fbeF[kb];
- f_TW = ftwF[kw];
- f_TN = ftnF[kzero];
- f_BS = fbsF[kbs];
- f_BN = fbnF[kb];
- f_TS = ftsF[ks];
- f_ZERO = fzeroF[kzero];
- f_TNE = ftneF[kzero];
- f_TSW = ftswF[ksw];
- f_TSE = ftseF[ks];
- f_TNW = ftnwF[kw];
- f_BNE = fbneF[kb];
- f_BSW = fbswF[kbsw];
- f_BSE = fbseF[kbs];
- f_BNW = fbnwF[kbw];
-
- drho_NWT = f_E+f_W+f_N+f_S+f_T+f_B+f_NE+f_SW+f_SE+f_NW+f_TE+f_BW+f_BE+f_TW+f_TN+f_BS+f_BN+f_TS+f_ZERO+f_TNE+f_TSW+f_TSE+f_TNW+f_BNE+f_BSW+f_BSE+f_BNW;
- vx1_NWT = (((f_TNE-f_BSW)+(f_TSE-f_BNW)+(f_BNE-f_TSW)+(f_BSE-f_TNW)) + (((f_NE-f_SW)+(f_TE-f_BW))+((f_SE-f_NW)+(f_BE-f_TW))) + (f_E-f_W))/(c1o1 + drho_NWT);
- vx2_NWT = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_BNE-f_TSW)+(f_BNW-f_TSE)) + (((f_NE-f_SW)+(f_TN-f_BS))+((f_BN-f_TS)+(f_NW-f_SE))) + (f_N-f_S))/(c1o1 + drho_NWT);
- vx3_NWT = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_TSE-f_BNW)+(f_TSW-f_BNE)) + (((f_TE-f_BW)+(f_TN-f_BS))+((f_TW-f_BE)+(f_TS-f_BN))) + (f_T-f_B))/(c1o1 + drho_NWT);
-
- kxyFromfcNEQ_NWT = -c3o1*omegaS*((f_SW+f_BSW+f_TSW-f_NW-f_BNW-f_TNW-f_SE-f_BSE-f_TSE+f_NE+f_BNE+f_TNE ) / (c1o1 + drho_NWT) - ((vx1_NWT*vx2_NWT)));
- kyzFromfcNEQ_NWT = -c3o1*omegaS*((f_BS+f_BSE+f_BSW-f_TS-f_TSE-f_TSW-f_BN-f_BNE-f_BNW+f_TN+f_TNE+f_TNW ) / (c1o1 + drho_NWT) - ((vx2_NWT*vx3_NWT)));
- kxzFromfcNEQ_NWT = -c3o1*omegaS*((f_BW+f_BSW+f_BNW-f_TW-f_TSW-f_TNW-f_BE-f_BSE-f_BNE+f_TE+f_TSE+f_TNE ) / (c1o1 + drho_NWT) - ((vx1_NWT*vx3_NWT)));
- kxxMyyFromfcNEQ_NWT = -c3o2*omegaS *((f_BW+f_W+f_TW-f_BS-f_S-f_TS-f_BN-f_N-f_TN+f_BE+f_E+f_TE ) / (c1o1 + drho_NWT) - ((vx1_NWT*vx1_NWT-vx2_NWT*vx2_NWT)));
- kxxMzzFromfcNEQ_NWT = -c3o2*omegaS *((f_SW+f_W+f_NW-f_BS-f_TS-f_B-f_T-f_BN-f_TN+f_SE+f_E+f_NE ) / (c1o1 + drho_NWT) - ((vx1_NWT*vx1_NWT-vx3_NWT*vx3_NWT)));
-
- //////////////////////////////////////////////////////////////////////////
- //NET//
- //////////////////////////////////////////////////////////////////////////
- //index
- kzero= kw;
- kw = neighborFX[kw];
- ks = ksw;
- kb = kbw;
- ksw = neighborFX[ksw];
- kbw = neighborFX[kbw];
- kbs = kbsw;
- kbsw = neighborFX[kbsw];
- ////////////////////////////////////////////////////////////////////////////////
- f_E = feF[kzero];
- f_W = fwF[kw];
- f_N = fnF[kzero];
- f_S = fsF[ks];
- f_T = ftF[kzero];
- f_B = fbF[kb];
- f_NE = fneF[kzero];
- f_SW = fswF[ksw];
- f_SE = fseF[ks];
- f_NW = fnwF[kw];
- f_TE = fteF[kzero];
- f_BW = fbwF[kbw];
- f_BE = fbeF[kb];
- f_TW = ftwF[kw];
- f_TN = ftnF[kzero];
- f_BS = fbsF[kbs];
- f_BN = fbnF[kb];
- f_TS = ftsF[ks];
- f_ZERO = fzeroF[kzero];
- f_TNE = ftneF[kzero];
- f_TSW = ftswF[ksw];
- f_TSE = ftseF[ks];
- f_TNW = ftnwF[kw];
- f_BNE = fbneF[kb];
- f_BSW = fbswF[kbsw];
- f_BSE = fbseF[kbs];
- f_BNW = fbnwF[kbw];
-
- drho_NET = f_E+f_W+f_N+f_S+f_T+f_B+f_NE+f_SW+f_SE+f_NW+f_TE+f_BW+f_BE+f_TW+f_TN+f_BS+f_BN+f_TS+f_ZERO+f_TNE+f_TSW+f_TSE+f_TNW+f_BNE+f_BSW+f_BSE+f_BNW;
- vx1_NET = (((f_TNE-f_BSW)+(f_TSE-f_BNW)+(f_BNE-f_TSW)+(f_BSE-f_TNW)) + (((f_NE-f_SW)+(f_TE-f_BW))+((f_SE-f_NW)+(f_BE-f_TW))) + (f_E-f_W))/(c1o1 + drho_NET);
- vx2_NET = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_BNE-f_TSW)+(f_BNW-f_TSE)) + (((f_NE-f_SW)+(f_TN-f_BS))+((f_BN-f_TS)+(f_NW-f_SE))) + (f_N-f_S))/(c1o1 + drho_NET);
- vx3_NET = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_TSE-f_BNW)+(f_TSW-f_BNE)) + (((f_TE-f_BW)+(f_TN-f_BS))+((f_TW-f_BE)+(f_TS-f_BN))) + (f_T-f_B))/(c1o1 + drho_NET);
-
- kxyFromfcNEQ_NET = -c3o1*omegaS*((f_SW+f_BSW+f_TSW-f_NW-f_BNW-f_TNW-f_SE-f_BSE-f_TSE+f_NE+f_BNE+f_TNE ) / (c1o1 + drho_NET) - ((vx1_NET*vx2_NET)));
- kyzFromfcNEQ_NET = -c3o1*omegaS*((f_BS+f_BSE+f_BSW-f_TS-f_TSE-f_TSW-f_BN-f_BNE-f_BNW+f_TN+f_TNE+f_TNW ) / (c1o1 + drho_NET) - ((vx2_NET*vx3_NET)));
- kxzFromfcNEQ_NET = -c3o1*omegaS*((f_BW+f_BSW+f_BNW-f_TW-f_TSW-f_TNW-f_BE-f_BSE-f_BNE+f_TE+f_TSE+f_TNE ) / (c1o1 + drho_NET) - ((vx1_NET*vx3_NET)));
- kxxMyyFromfcNEQ_NET = -c3o2*omegaS *((f_BW+f_W+f_TW-f_BS-f_S-f_TS-f_BN-f_N-f_TN+f_BE+f_E+f_TE ) / (c1o1 + drho_NET) - ((vx1_NET*vx1_NET-vx2_NET*vx2_NET)));
- kxxMzzFromfcNEQ_NET = -c3o2*omegaS *((f_SW+f_W+f_NW-f_BS-f_TS-f_B-f_T-f_BN-f_TN+f_SE+f_E+f_NE ) / (c1o1 + drho_NET) - ((vx1_NET*vx1_NET-vx3_NET*vx3_NET)));
-
- //////////////////////////////////////////////////////////////////////////
- //NEB//
- //////////////////////////////////////////////////////////////////////////
- //index
- kb = kzero;
- kbw = kw;
- kbs = ks;
- kbsw = ksw;
- kzero= k0w;
- kw = neighborFX[k0w];
- ks = k0sw;
- ksw = neighborFX[k0sw];
- ////////////////////////////////////////////////////////////////////////////////
- f_E = feF[kzero];
- f_W = fwF[kw];
- f_N = fnF[kzero];
- f_S = fsF[ks];
- f_T = ftF[kzero];
- f_B = fbF[kb];
- f_NE = fneF[kzero];
- f_SW = fswF[ksw];
- f_SE = fseF[ks];
- f_NW = fnwF[kw];
- f_TE = fteF[kzero];
- f_BW = fbwF[kbw];
- f_BE = fbeF[kb];
- f_TW = ftwF[kw];
- f_TN = ftnF[kzero];
- f_BS = fbsF[kbs];
- f_BN = fbnF[kb];
- f_TS = ftsF[ks];
- f_ZERO = fzeroF[kzero];
- f_TNE = ftneF[kzero];
- f_TSW = ftswF[ksw];
- f_TSE = ftseF[ks];
- f_TNW = ftnwF[kw];
- f_BNE = fbneF[kb];
- f_BSW = fbswF[kbsw];
- f_BSE = fbseF[kbs];
- f_BNW = fbnwF[kbw];
-
- drho_NEB = f_E+f_W+f_N+f_S+f_T+f_B+f_NE+f_SW+f_SE+f_NW+f_TE+f_BW+f_BE+f_TW+f_TN+f_BS+f_BN+f_TS+f_ZERO+f_TNE+f_TSW+f_TSE+f_TNW+f_BNE+f_BSW+f_BSE+f_BNW;
- vx1_NEB = (((f_TNE-f_BSW)+(f_TSE-f_BNW)+(f_BNE-f_TSW)+(f_BSE-f_TNW)) + (((f_NE-f_SW)+(f_TE-f_BW))+((f_SE-f_NW)+(f_BE-f_TW))) + (f_E-f_W))/(c1o1 + drho_NEB);
- vx2_NEB = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_BNE-f_TSW)+(f_BNW-f_TSE)) + (((f_NE-f_SW)+(f_TN-f_BS))+((f_BN-f_TS)+(f_NW-f_SE))) + (f_N-f_S))/(c1o1 + drho_NEB);
- vx3_NEB = (((f_TNE-f_BSW)+(f_TNW-f_BSE)+(f_TSE-f_BNW)+(f_TSW-f_BNE)) + (((f_TE-f_BW)+(f_TN-f_BS))+((f_TW-f_BE)+(f_TS-f_BN))) + (f_T-f_B))/(c1o1 + drho_NEB);
-
- kxyFromfcNEQ_NEB = -c3o1*omegaS*((f_SW+f_BSW+f_TSW-f_NW-f_BNW-f_TNW-f_SE-f_BSE-f_TSE+f_NE+f_BNE+f_TNE ) / (c1o1 + drho_NEB) - ((vx1_NEB*vx2_NEB)));
- kyzFromfcNEQ_NEB = -c3o1*omegaS*((f_BS+f_BSE+f_BSW-f_TS-f_TSE-f_TSW-f_BN-f_BNE-f_BNW+f_TN+f_TNE+f_TNW ) / (c1o1 + drho_NEB) - ((vx2_NEB*vx3_NEB)));
- kxzFromfcNEQ_NEB = -c3o1*omegaS*((f_BW+f_BSW+f_BNW-f_TW-f_TSW-f_TNW-f_BE-f_BSE-f_BNE+f_TE+f_TSE+f_TNE ) / (c1o1 + drho_NEB) - ((vx1_NEB*vx3_NEB)));
- kxxMyyFromfcNEQ_NEB = -c3o2*omegaS *((f_BW+f_W+f_TW-f_BS-f_S-f_TS-f_BN-f_N-f_TN+f_BE+f_E+f_TE ) / (c1o1 + drho_NEB) - ((vx1_NEB*vx1_NEB-vx2_NEB*vx2_NEB)));
- kxxMzzFromfcNEQ_NEB = -c3o2*omegaS *((f_SW+f_W+f_NW-f_BS-f_TS-f_B-f_T-f_BN-f_TN+f_SE+f_E+f_NE ) / (c1o1 + drho_NEB) - ((vx1_NEB*vx1_NEB-vx3_NEB*vx3_NEB)));
-
- //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- //kxyFromfcNEQ_SWB = zero;
- //kyzFromfcNEQ_SWB = zero;
- //kxzFromfcNEQ_SWB = zero;
- //kxxMyyFromfcNEQ_SWB = zero;
- //kxxMzzFromfcNEQ_SWB = zero;
- //kxyFromfcNEQ_SWT = zero;
- //kyzFromfcNEQ_SWT = zero;
- //kxzFromfcNEQ_SWT = zero;
- //kxxMyyFromfcNEQ_SWT = zero;
- //kxxMzzFromfcNEQ_SWT = zero;
- //kxyFromfcNEQ_SET = zero;
- //kyzFromfcNEQ_SET = zero;
- //kxzFromfcNEQ_SET = zero;
- //kxxMyyFromfcNEQ_SET = zero;
- //kxxMzzFromfcNEQ_SET = zero;
- //kxyFromfcNEQ_SEB = zero;
- //kyzFromfcNEQ_SEB = zero;
- //kxzFromfcNEQ_SEB = zero;
- //kxxMyyFromfcNEQ_SEB = zero;
- //kxxMzzFromfcNEQ_SEB = zero;
- //kxyFromfcNEQ_NWB = zero;
- //kyzFromfcNEQ_NWB = zero;
- //kxzFromfcNEQ_NWB = zero;
- //kxxMyyFromfcNEQ_NWB = zero;
- //kxxMzzFromfcNEQ_NWB = zero;
- //kxyFromfcNEQ_NWT = zero;
- //kyzFromfcNEQ_NWT = zero;
- //kxzFromfcNEQ_NWT = zero;
- //kxxMyyFromfcNEQ_NWT = zero;
- //kxxMzzFromfcNEQ_NWT = zero;
- //kxyFromfcNEQ_NET = zero;
- //kyzFromfcNEQ_NET = zero;
- //kxzFromfcNEQ_NET = zero;
- //kxxMyyFromfcNEQ_NET = zero;
- //kxxMzzFromfcNEQ_NET = zero;
- //kxyFromfcNEQ_NEB = zero;
- //kyzFromfcNEQ_NEB = zero;
- //kxzFromfcNEQ_NEB = zero;
- //kxxMyyFromfcNEQ_NEB = zero;
- //kxxMzzFromfcNEQ_NEB = zero;
- //////////////////////////////////////////////////////////////////////////
- //3
- //////////////////////////////////////////////////////////////////////////
- a0 = (-kxxMyyFromfcNEQ_NEB - kxxMyyFromfcNEQ_NET + kxxMyyFromfcNEQ_NWB + kxxMyyFromfcNEQ_NWT -
- kxxMyyFromfcNEQ_SEB - kxxMyyFromfcNEQ_SET + kxxMyyFromfcNEQ_SWB + kxxMyyFromfcNEQ_SWT -
- kxxMzzFromfcNEQ_NEB - kxxMzzFromfcNEQ_NET + kxxMzzFromfcNEQ_NWB + kxxMzzFromfcNEQ_NWT -
- kxxMzzFromfcNEQ_SEB - kxxMzzFromfcNEQ_SET + kxxMzzFromfcNEQ_SWB + kxxMzzFromfcNEQ_SWT -
- c2o1*kxyFromfcNEQ_NEB - c2o1*kxyFromfcNEQ_NET - c2o1*kxyFromfcNEQ_NWB - c2o1*kxyFromfcNEQ_NWT +
- c2o1*kxyFromfcNEQ_SEB + c2o1*kxyFromfcNEQ_SET + c2o1*kxyFromfcNEQ_SWB + c2o1*kxyFromfcNEQ_SWT +
- c2o1*kxzFromfcNEQ_NEB - c2o1*kxzFromfcNEQ_NET + c2o1*kxzFromfcNEQ_NWB - c2o1*kxzFromfcNEQ_NWT +
- c2o1*kxzFromfcNEQ_SEB - c2o1*kxzFromfcNEQ_SET + c2o1*kxzFromfcNEQ_SWB - c2o1*kxzFromfcNEQ_SWT +
- c8o1*vx1_NEB + c8o1*vx1_NET + c8o1*vx1_NWB + c8o1*vx1_NWT + c8o1*vx1_SEB +
- c8o1*vx1_SET + c8o1*vx1_SWB + c8o1*vx1_SWT + c2o1*vx2_NEB + c2o1*vx2_NET -
- c2o1*vx2_NWB - c2o1*vx2_NWT - c2o1*vx2_SEB - c2o1*vx2_SET + c2o1*vx2_SWB +
- c2o1*vx2_SWT - c2o1*vx3_NEB + c2o1*vx3_NET + c2o1*vx3_NWB - c2o1*vx3_NWT -
- c2o1*vx3_SEB + c2o1*vx3_SET + c2o1*vx3_SWB - c2o1*vx3_SWT)/c64o1;
- b0 = (c2o1*kxxMyyFromfcNEQ_NEB + c2o1*kxxMyyFromfcNEQ_NET + c2o1*kxxMyyFromfcNEQ_NWB + c2o1*kxxMyyFromfcNEQ_NWT -
- c2o1*kxxMyyFromfcNEQ_SEB - c2o1*kxxMyyFromfcNEQ_SET - c2o1*kxxMyyFromfcNEQ_SWB - c2o1*kxxMyyFromfcNEQ_SWT -
- kxxMzzFromfcNEQ_NEB - kxxMzzFromfcNEQ_NET - kxxMzzFromfcNEQ_NWB - kxxMzzFromfcNEQ_NWT +
- kxxMzzFromfcNEQ_SEB + kxxMzzFromfcNEQ_SET + kxxMzzFromfcNEQ_SWB + kxxMzzFromfcNEQ_SWT -
- c2o1*kxyFromfcNEQ_NEB - c2o1*kxyFromfcNEQ_NET + c2o1*kxyFromfcNEQ_NWB + c2o1*kxyFromfcNEQ_NWT -
- c2o1*kxyFromfcNEQ_SEB - c2o1*kxyFromfcNEQ_SET + c2o1*kxyFromfcNEQ_SWB + c2o1*kxyFromfcNEQ_SWT +
- c2o1*kyzFromfcNEQ_NEB - c2o1*kyzFromfcNEQ_NET + c2o1*kyzFromfcNEQ_NWB - c2o1*kyzFromfcNEQ_NWT +
- c2o1*kyzFromfcNEQ_SEB - c2o1*kyzFromfcNEQ_SET + c2o1*kyzFromfcNEQ_SWB - c2o1*kyzFromfcNEQ_SWT +
- c2o1*vx1_NEB + c2o1*vx1_NET - c2o1*vx1_NWB - c2o1*vx1_NWT -
- c2o1*vx1_SEB - c2o1*vx1_SET + c2o1*vx1_SWB + c2o1*vx1_SWT +
- c8o1*vx2_NEB + c8o1*vx2_NET + c8o1*vx2_NWB + c8o1*vx2_NWT +
- c8o1*vx2_SEB + c8o1*vx2_SET + c8o1*vx2_SWB + c8o1*vx2_SWT -
- c2o1*vx3_NEB + c2o1*vx3_NET - c2o1*vx3_NWB + c2o1*vx3_NWT +
- c2o1*vx3_SEB - c2o1*vx3_SET + c2o1*vx3_SWB - c2o1*vx3_SWT)/c64o1;
- c0 = (kxxMyyFromfcNEQ_NEB - kxxMyyFromfcNEQ_NET + kxxMyyFromfcNEQ_NWB - kxxMyyFromfcNEQ_NWT +
- kxxMyyFromfcNEQ_SEB - kxxMyyFromfcNEQ_SET + kxxMyyFromfcNEQ_SWB - kxxMyyFromfcNEQ_SWT -
- c2o1*kxxMzzFromfcNEQ_NEB + c2o1*kxxMzzFromfcNEQ_NET - c2o1*kxxMzzFromfcNEQ_NWB + c2o1*kxxMzzFromfcNEQ_NWT -
- c2o1*kxxMzzFromfcNEQ_SEB + c2o1*kxxMzzFromfcNEQ_SET - c2o1*kxxMzzFromfcNEQ_SWB + c2o1*kxxMzzFromfcNEQ_SWT -
- c2o1*kxzFromfcNEQ_NEB - c2o1*kxzFromfcNEQ_NET + c2o1*kxzFromfcNEQ_NWB + c2o1*kxzFromfcNEQ_NWT -
- c2o1*kxzFromfcNEQ_SEB - c2o1*kxzFromfcNEQ_SET + c2o1*kxzFromfcNEQ_SWB + c2o1*kxzFromfcNEQ_SWT -
- c2o1*kyzFromfcNEQ_NEB - c2o1*kyzFromfcNEQ_NET - c2o1*kyzFromfcNEQ_NWB - c2o1*kyzFromfcNEQ_NWT +
- c2o1*kyzFromfcNEQ_SEB + c2o1*kyzFromfcNEQ_SET + c2o1*kyzFromfcNEQ_SWB + c2o1*kyzFromfcNEQ_SWT -
- c2o1*vx1_NEB + c2o1*vx1_NET + c2o1*vx1_NWB - c2o1*vx1_NWT -
- c2o1*vx1_SEB + c2o1*vx1_SET + c2o1*vx1_SWB - c2o1*vx1_SWT -
- c2o1*vx2_NEB + c2o1*vx2_NET - c2o1*vx2_NWB + c2o1*vx2_NWT +
- c2o1*vx2_SEB - c2o1*vx2_SET + c2o1*vx2_SWB - c2o1*vx2_SWT +
- c8o1*vx3_NEB + c8o1*vx3_NET + c8o1*vx3_NWB + c8o1*vx3_NWT +
- c8o1*vx3_SEB + c8o1*vx3_SET + c8o1*vx3_SWB + c8o1*vx3_SWT)/c64o1;
- ax = (vx1_NEB + vx1_NET - vx1_NWB - vx1_NWT + vx1_SEB + vx1_SET - vx1_SWB - vx1_SWT)/c4o1;
- bx = (vx2_NEB + vx2_NET - vx2_NWB - vx2_NWT + vx2_SEB + vx2_SET - vx2_SWB - vx2_SWT)/c4o1;
- cx = (vx3_NEB + vx3_NET - vx3_NWB - vx3_NWT + vx3_SEB + vx3_SET - vx3_SWB - vx3_SWT)/c4o1;
- axx= (kxxMyyFromfcNEQ_NEB + kxxMyyFromfcNEQ_NET - kxxMyyFromfcNEQ_NWB - kxxMyyFromfcNEQ_NWT +
- kxxMyyFromfcNEQ_SEB + kxxMyyFromfcNEQ_SET - kxxMyyFromfcNEQ_SWB - kxxMyyFromfcNEQ_SWT +
- kxxMzzFromfcNEQ_NEB + kxxMzzFromfcNEQ_NET - kxxMzzFromfcNEQ_NWB - kxxMzzFromfcNEQ_NWT +
- kxxMzzFromfcNEQ_SEB + kxxMzzFromfcNEQ_SET - kxxMzzFromfcNEQ_SWB - kxxMzzFromfcNEQ_SWT +
- c2o1*vx2_NEB + c2o1*vx2_NET - c2o1*vx2_NWB - c2o1*vx2_NWT -
- c2o1*vx2_SEB - c2o1*vx2_SET + c2o1*vx2_SWB + c2o1*vx2_SWT -
- c2o1*vx3_NEB + c2o1*vx3_NET + c2o1*vx3_NWB - c2o1*vx3_NWT -
- c2o1*vx3_SEB + c2o1*vx3_SET + c2o1*vx3_SWB - c2o1*vx3_SWT)/c16o1;
- bxx= (kxyFromfcNEQ_NEB + kxyFromfcNEQ_NET - kxyFromfcNEQ_NWB - kxyFromfcNEQ_NWT +
- kxyFromfcNEQ_SEB + kxyFromfcNEQ_SET - kxyFromfcNEQ_SWB - kxyFromfcNEQ_SWT -
- c2o1*vx1_NEB - c2o1*vx1_NET + c2o1*vx1_NWB + c2o1*vx1_NWT +
- c2o1*vx1_SEB + c2o1*vx1_SET - c2o1*vx1_SWB - c2o1*vx1_SWT)/c8o1;
- cxx= (kxzFromfcNEQ_NEB + kxzFromfcNEQ_NET - kxzFromfcNEQ_NWB - kxzFromfcNEQ_NWT +
- kxzFromfcNEQ_SEB + kxzFromfcNEQ_SET - kxzFromfcNEQ_SWB - kxzFromfcNEQ_SWT +
- c2o1*vx1_NEB - c2o1*vx1_NET - c2o1*vx1_NWB + c2o1*vx1_NWT +
- c2o1*vx1_SEB - c2o1*vx1_SET - c2o1*vx1_SWB + c2o1*vx1_SWT)/c8o1;
- ay = (vx1_NEB + vx1_NET + vx1_NWB + vx1_NWT - vx1_SEB - vx1_SET - vx1_SWB - vx1_SWT)/c4o1;
- by = (vx2_NEB + vx2_NET + vx2_NWB + vx2_NWT - vx2_SEB - vx2_SET - vx2_SWB - vx2_SWT)/c4o1;
- cy = (vx3_NEB + vx3_NET + vx3_NWB + vx3_NWT - vx3_SEB - vx3_SET - vx3_SWB - vx3_SWT)/c4o1;
- ayy= (kxyFromfcNEQ_NEB + kxyFromfcNEQ_NET + kxyFromfcNEQ_NWB + kxyFromfcNEQ_NWT -
- kxyFromfcNEQ_SEB - kxyFromfcNEQ_SET - kxyFromfcNEQ_SWB - kxyFromfcNEQ_SWT -
- c2o1*vx2_NEB - c2o1*vx2_NET + c2o1*vx2_NWB + c2o1*vx2_NWT +
- c2o1*vx2_SEB + c2o1*vx2_SET - c2o1*vx2_SWB - c2o1*vx2_SWT)/c8o1;
- byy= (-c2o1*kxxMyyFromfcNEQ_NEB - c2o1*kxxMyyFromfcNEQ_NET - c2o1*kxxMyyFromfcNEQ_NWB - c2o1*kxxMyyFromfcNEQ_NWT +
- c2o1*kxxMyyFromfcNEQ_SEB + c2o1*kxxMyyFromfcNEQ_SET + c2o1*kxxMyyFromfcNEQ_SWB + c2o1*kxxMyyFromfcNEQ_SWT +
- kxxMzzFromfcNEQ_NEB + kxxMzzFromfcNEQ_NET + kxxMzzFromfcNEQ_NWB + kxxMzzFromfcNEQ_NWT -
- kxxMzzFromfcNEQ_SEB - kxxMzzFromfcNEQ_SET - kxxMzzFromfcNEQ_SWB - kxxMzzFromfcNEQ_SWT +
- c2o1*vx1_NEB + c2o1*vx1_NET - c2o1*vx1_NWB - c2o1*vx1_NWT -
- c2o1*vx1_SEB - c2o1*vx1_SET + c2o1*vx1_SWB + c2o1*vx1_SWT -
- c2o1*vx3_NEB + c2o1*vx3_NET - c2o1*vx3_NWB + c2o1*vx3_NWT +
- c2o1*vx3_SEB - c2o1*vx3_SET + c2o1*vx3_SWB - c2o1*vx3_SWT)/c16o1;
- cyy= (kyzFromfcNEQ_NEB + kyzFromfcNEQ_NET + kyzFromfcNEQ_NWB + kyzFromfcNEQ_NWT -
- kyzFromfcNEQ_SEB - kyzFromfcNEQ_SET - kyzFromfcNEQ_SWB - kyzFromfcNEQ_SWT +
- c2o1*vx2_NEB - c2o1*vx2_NET + c2o1*vx2_NWB - c2o1*vx2_NWT -
- c2o1*vx2_SEB + c2o1*vx2_SET - c2o1*vx2_SWB + c2o1*vx2_SWT)/c8o1;
- az = (-vx1_NEB + vx1_NET - vx1_NWB + vx1_NWT - vx1_SEB + vx1_SET - vx1_SWB + vx1_SWT)/c4o1;
- bz = (-vx2_NEB + vx2_NET - vx2_NWB + vx2_NWT - vx2_SEB + vx2_SET - vx2_SWB + vx2_SWT)/c4o1;
- cz = (-vx3_NEB + vx3_NET - vx3_NWB + vx3_NWT - vx3_SEB + vx3_SET - vx3_SWB + vx3_SWT)/c4o1;
- azz= (-kxzFromfcNEQ_NEB + kxzFromfcNEQ_NET - kxzFromfcNEQ_NWB + kxzFromfcNEQ_NWT -
- kxzFromfcNEQ_SEB + kxzFromfcNEQ_SET - kxzFromfcNEQ_SWB + kxzFromfcNEQ_SWT +
- c2o1*vx3_NEB - c2o1*vx3_NET - c2o1*vx3_NWB + c2o1*vx3_NWT +
- c2o1*vx3_SEB - c2o1*vx3_SET - c2o1*vx3_SWB + c2o1*vx3_SWT)/c8o1;
- bzz= (-kyzFromfcNEQ_NEB + kyzFromfcNEQ_NET - kyzFromfcNEQ_NWB + kyzFromfcNEQ_NWT -
- kyzFromfcNEQ_SEB + kyzFromfcNEQ_SET - kyzFromfcNEQ_SWB + kyzFromfcNEQ_SWT +
- c2o1*vx3_NEB - c2o1*vx3_NET + c2o1*vx3_NWB - c2o1*vx3_NWT -
- c2o1*vx3_SEB + c2o1*vx3_SET - c2o1*vx3_SWB + c2o1*vx3_SWT)/c8o1;
- czz= (-kxxMyyFromfcNEQ_NEB + kxxMyyFromfcNEQ_NET - kxxMyyFromfcNEQ_NWB + kxxMyyFromfcNEQ_NWT -
- kxxMyyFromfcNEQ_SEB + kxxMyyFromfcNEQ_SET - kxxMyyFromfcNEQ_SWB + kxxMyyFromfcNEQ_SWT +
- c2o1*kxxMzzFromfcNEQ_NEB - c2o1*kxxMzzFromfcNEQ_NET + c2o1*kxxMzzFromfcNEQ_NWB - c2o1*kxxMzzFromfcNEQ_NWT +
- c2o1*kxxMzzFromfcNEQ_SEB - c2o1*kxxMzzFromfcNEQ_SET + c2o1*kxxMzzFromfcNEQ_SWB - c2o1*kxxMzzFromfcNEQ_SWT -
- c2o1*vx1_NEB + c2o1*vx1_NET + c2o1*vx1_NWB - c2o1*vx1_NWT -
- c2o1*vx1_SEB + c2o1*vx1_SET + c2o1*vx1_SWB - c2o1*vx1_SWT -
- c2o1*vx2_NEB + c2o1*vx2_NET - c2o1*vx2_NWB + c2o1*vx2_NWT +
- c2o1*vx2_SEB - c2o1*vx2_SET + c2o1*vx2_SWB - c2o1*vx2_SWT)/c16o1;
- axy= (vx1_NEB + vx1_NET - vx1_NWB - vx1_NWT - vx1_SEB - vx1_SET + vx1_SWB + vx1_SWT)/c2o1;
- bxy= (vx2_NEB + vx2_NET - vx2_NWB - vx2_NWT - vx2_SEB - vx2_SET + vx2_SWB + vx2_SWT)/c2o1;
- cxy= (vx3_NEB + vx3_NET - vx3_NWB - vx3_NWT - vx3_SEB - vx3_SET + vx3_SWB + vx3_SWT)/c2o1;
- axz= (-vx1_NEB + vx1_NET + vx1_NWB - vx1_NWT - vx1_SEB + vx1_SET + vx1_SWB - vx1_SWT)/c2o1;
- bxz= (-vx2_NEB + vx2_NET + vx2_NWB - vx2_NWT - vx2_SEB + vx2_SET + vx2_SWB - vx2_SWT)/c2o1;
- cxz= (-vx3_NEB + vx3_NET + vx3_NWB - vx3_NWT - vx3_SEB + vx3_SET + vx3_SWB - vx3_SWT)/c2o1;
- ayz= (-vx1_NEB + vx1_NET - vx1_NWB + vx1_NWT + vx1_SEB - vx1_SET + vx1_SWB - vx1_SWT)/c2o1;
- byz= (-vx2_NEB + vx2_NET - vx2_NWB + vx2_NWT + vx2_SEB - vx2_SET + vx2_SWB - vx2_SWT)/c2o1;
- cyz= (-vx3_NEB + vx3_NET - vx3_NWB + vx3_NWT + vx3_SEB - vx3_SET + vx3_SWB - vx3_SWT)/c2o1;
- //axyz=-vx1_NEB + vx1_NET + vx1_NWB - vx1_NWT + vx1_SEB - vx1_SET - vx1_SWB + vx1_SWT;
- //bxyz=-vx2_NEB + vx2_NET + vx2_NWB - vx2_NWT + vx2_SEB - vx2_SET - vx2_SWB + vx2_SWT;
- //cxyz=-vx3_NEB + vx3_NET + vx3_NWB - vx3_NWT + vx3_SEB - vx3_SET - vx3_SWB + vx3_SWT;
- //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- real kxyAverage = c0o1;
- real kyzAverage = c0o1;
- real kxzAverage = c0o1;
- real kxxMyyAverage = c0o1;
- real kxxMzzAverage = c0o1;
- //real kxyAverage =(kxyFromfcNEQ_SWB+
- // kxyFromfcNEQ_SWT+
- // kxyFromfcNEQ_SET+
- // kxyFromfcNEQ_SEB+
- // kxyFromfcNEQ_NWB+
- // kxyFromfcNEQ_NWT+
- // kxyFromfcNEQ_NET+
- // kxyFromfcNEQ_NEB)*c1o8-(ay+bx);
- //real kyzAverage =(kyzFromfcNEQ_SWB+
- // kyzFromfcNEQ_SWT+
- // kyzFromfcNEQ_SET+
- // kyzFromfcNEQ_SEB+
- // kyzFromfcNEQ_NWB+
- // kyzFromfcNEQ_NWT+
- // kyzFromfcNEQ_NET+
- // kyzFromfcNEQ_NEB)*c1o8-(bz+cy);
- //real kxzAverage =(kxzFromfcNEQ_SWB+
- // kxzFromfcNEQ_SWT+
- // kxzFromfcNEQ_SET+
- // kxzFromfcNEQ_SEB+
- // kxzFromfcNEQ_NWB+
- // kxzFromfcNEQ_NWT+
- // kxzFromfcNEQ_NET+
- // kxzFromfcNEQ_NEB)*c1o8-(az+cx);
- //real kxxMyyAverage =(kxxMyyFromfcNEQ_SWB+
- // kxxMyyFromfcNEQ_SWT+
- // kxxMyyFromfcNEQ_SET+
- // kxxMyyFromfcNEQ_SEB+
- // kxxMyyFromfcNEQ_NWB+
- // kxxMyyFromfcNEQ_NWT+
- // kxxMyyFromfcNEQ_NET+
- // kxxMyyFromfcNEQ_NEB)*c1o8-(ax-by);
- //real kxxMzzAverage =(kxxMzzFromfcNEQ_SWB+
- // kxxMzzFromfcNEQ_SWT+
- // kxxMzzFromfcNEQ_SET+
- // kxxMzzFromfcNEQ_SEB+
- // kxxMzzFromfcNEQ_NWB+
- // kxxMzzFromfcNEQ_NWT+
- // kxxMzzFromfcNEQ_NET+
- // kxxMzzFromfcNEQ_NEB)*c1o8-(ax-cz);
-
-
-
- //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- ////Press
- //d0 = ( press_NEB + press_NET + press_NWB + press_NWT + press_SEB + press_SET + press_SWB + press_SWT) * c1o8;
- //dx = ( press_NEB + press_NET - press_NWB - press_NWT + press_SEB + press_SET - press_SWB - press_SWT) * c1o4;
- //dy = ( press_NEB + press_NET + press_NWB + press_NWT - press_SEB - press_SET - press_SWB - press_SWT) * c1o4;
- //dz = (-press_NEB + press_NET - press_NWB + press_NWT - press_SEB + press_SET - press_SWB + press_SWT) * c1o4;
- //dxy = ( press_NEB + press_NET - press_NWB - press_NWT - press_SEB - press_SET + press_SWB + press_SWT) * c1o2;
- //dxz = (-press_NEB + press_NET + press_NWB - press_NWT - press_SEB + press_SET + press_SWB - press_SWT) * c1o2;
- //dyz = (-press_NEB + press_NET - press_NWB + press_NWT + press_SEB - press_SET + press_SWB - press_SWT) * c1o2;
- //dxyz = -press_NEB + press_NET + press_NWB - press_NWT + press_SEB - press_SET - press_SWB + press_SWT;
- //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- //drho
- real LapRho = ((xoff != c0o1) || (yoff != c0o1) || (zoff != c0o1)) ? c0o1 : -c3o1*(ax*ax + by*by + cz*cz) - c6o1 * (bx*ay + cx*az + cy*bz);
- d0 = ( drho_NEB + drho_NET + drho_NWB + drho_NWT + drho_SEB + drho_SET + drho_SWB + drho_SWT - c2o1*LapRho) * c1o8;
- dx = ( drho_NEB + drho_NET - drho_NWB - drho_NWT + drho_SEB + drho_SET - drho_SWB - drho_SWT) * c1o4;
- dy = ( drho_NEB + drho_NET + drho_NWB + drho_NWT - drho_SEB - drho_SET - drho_SWB - drho_SWT) * c1o4;
- dz = (-drho_NEB + drho_NET - drho_NWB + drho_NWT - drho_SEB + drho_SET - drho_SWB + drho_SWT) * c1o4;
- dxy = ( drho_NEB + drho_NET - drho_NWB - drho_NWT - drho_SEB - drho_SET + drho_SWB + drho_SWT) * c1o2;
- dxz = (-drho_NEB + drho_NET + drho_NWB - drho_NWT - drho_SEB + drho_SET + drho_SWB - drho_SWT) * c1o2;
- dyz = (-drho_NEB + drho_NET - drho_NWB + drho_NWT + drho_SEB - drho_SET + drho_SWB - drho_SWT) * c1o2;
- //dxyz = -drho_NEB + drho_NET + drho_NWB - drho_NWT + drho_SEB - drho_SET - drho_SWB + drho_SWT;
- //d0 = zero;
- //dx = zero;
- //dy = zero;
- //dz = zero;
- //dxy = zero;
- //dxz = zero;
- //dyz = zero;
- //dxyz = zero;
- ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- //
- // Bernd das Brot
- //
- //
- // x------x
- // | |
- // | ---+--->X
- // | | \
- // x------x \
- // off-vector
- ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- a0 = a0 + xoff * ax + yoff * ay + zoff * az + xoff_sq * axx + yoff_sq * ayy + zoff_sq * azz + xoff*yoff*axy + xoff*zoff*axz + yoff*zoff*ayz;
- ax = ax + c2o1 * xoff * axx + yoff * axy + zoff * axz;
- ay = ay + c2o1 * yoff * ayy + xoff * axy + zoff * ayz;
- az = az + c2o1 * zoff * azz + xoff * axz + yoff * ayz;
- b0 = b0 + xoff * bx + yoff * by + zoff * bz + xoff_sq * bxx + yoff_sq * byy + zoff_sq * bzz + xoff*yoff*bxy + xoff*zoff*bxz + yoff*zoff*byz;
- bx = bx + c2o1 * xoff * bxx + yoff * bxy + zoff * bxz;
- by = by + c2o1 * yoff * byy + xoff * bxy + zoff * byz;
- bz = bz + c2o1 * zoff * bzz + xoff * bxz + yoff * byz;
- c0 = c0 + xoff * cx + yoff * cy + zoff * cz + xoff_sq * cxx + yoff_sq * cyy + zoff_sq * czz + xoff*yoff*cxy + xoff*zoff*cxz + yoff*zoff*cyz;
- cx = cx + c2o1 * xoff * cxx + yoff * cxy + zoff * cxz;
- cy = cy + c2o1 * yoff * cyy + xoff * cxy + zoff * cyz;
- cz = cz + c2o1 * zoff * czz + xoff * cxz + yoff * cyz;
- d0 = d0 + xoff * dx + yoff * dy + zoff * dz + xoff*yoff*dxy + xoff*zoff*dxz + yoff*zoff*dyz;
- ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- // FIX ///////////////////////////////////////////////////////////////////////////////////////////////////////////
- ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- //AAAAAAAAAAAAHHHHHHHHHHHH!!!!! Mieser Test!!!
- //b0= bx= by= bz= bxx= byy= bzz= bxy= bxz= byz= c0= cx= cy= cz= cxx= cyy= czz= cxy= cxz= cyz= axyz= bxyz= cxyz=zero;
- //b0=zero;
- //bx=zero;
- //by=zero;
- //bz=zero;
- //bxx=zero;
- //byy=zero;
- //bzz=zero;
- //bxy=zero;
- //bxz=zero;
- //byz=zero;
- //c0=zero;
- //cx=zero;
- //cy=zero;
- //cz=zero;
- //cxx=zero;
- //cyy=zero;
- //czz=zero;
- //cxy=zero;
- //cxz=zero;
- //cyz=zero;
- ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- real mfcbb = c0o1;
- real mfabb = c0o1;
- real mfbcb = c0o1;
- real mfbab = c0o1;
- real mfbbc = c0o1;
- real mfbba = c0o1;
- real mfccb = c0o1;
- real mfaab = c0o1;
- real mfcab = c0o1;
- real mfacb = c0o1;
- real mfcbc = c0o1;
- real mfaba = c0o1;
- real mfcba = c0o1;
- real mfabc = c0o1;
- real mfbcc = c0o1;
- real mfbaa = c0o1;
- real mfbca = c0o1;
- real mfbac = c0o1;
- real mfbbb = c0o1;
- real mfccc = c0o1;
- real mfaac = c0o1;
- real mfcac = c0o1;
- real mfacc = c0o1;
- real mfcca = c0o1;
- real mfaaa = c0o1;
- real mfcaa = c0o1;
- real mfaca = c0o1;
- ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- real m0, m1, m2, vvx, vvy, vvz, vx2, vy2, vz2, oMdrho;
- real mxxPyyPzz, mxxMyy, mxxMzz, mxxyPyzz, mxxyMyzz, mxxzPyyz, mxxzMyyz, mxyyPxzz, mxyyMxzz;
- //real qudricLimit = c1o100;//ganz schlechte Idee -> muss global sein
- //real O3 = c2o1 - o;
- //real residu, residutmp;
- //residutmp = c0o1;///*-*/ c2o9 * (1./o - c1o2) * eps_new * eps_new;
- real NeqOn = c1o1;//zero;//one; //.... one = on ..... zero = off
- ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
-
- ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- //
- //Position C 0., 0., 0.
- //
- ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- //x = 0.;
- //y = 0.;
- //z = 0.;
- ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- //real mxoff = -xoff;
- //real myoff = -yoff;
- //real mzoff = -zoff;
- ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
- //press = press_NET * (c1o8 - c1o4 * mxoff - c1o4 * myoff - c1o4 * mzoff) +
- // press_NWT * (c1o8 + c1o4 * mxoff - c1o4 * myoff - c1o4 * mzoff) +
- // press_SET * (c1o8 - c1o4 * mxoff + c1o4 * myoff - c1o4 * mzoff) +
- // press_SWT * (c1o8 + c1o4 * mxoff + c1o4 * myoff - c1o4 * mzoff) +
- // press_NEB * (c1o8 - c1o4 * mxoff - c1o4 * myoff + c1o4 * mzoff) +
- // press_NWB * (c1o8 + c1o4 * mxoff - c1o4 * myoff + c1o4 * mzoff) +
- // press_SEB * (c1o8 - c1o4 * mxoff + c1o4 * myoff + c1o4 * mzoff) +
- // press_SWB * (c1o8 + c1o4 * mxoff + c1o4 * myoff + c1o4 * mzoff);
- //drho = drho_NET * (c1o8 - c1o4 * xoff - c1o4 * yoff - c1o4 * zoff) +
- // drho_NWT * (c1o8 + c1o4 * xoff - c1o4 * yoff - c1o4 * zoff) +
- // drho_SET * (c1o8 - c1o4 * xoff + c1o4 * yoff - c1o4 * zoff) +
- // drho_SWT * (c1o8 + c1o4 * xoff + c1o4 * yoff - c1o4 * zoff) +
- // drho_NEB * (c1o8 - c1o4 * xoff - c1o4 * yoff + c1o4 * zoff) +
- // drho_NWB * (c1o8 + c1o4 * xoff - c1o4 * yoff + c1o4 * zoff) +
- // drho_SEB * (c1o8 - c1o4 * xoff + c1o4 * yoff + c1o4 * zoff) +
- // drho_SWB * (c1o8 + c1o4 * xoff + c1o4 * yoff + c1o4 * zoff);
- press = d0;
- vvx = a0;
- vvy = b0;
- vvz = c0;
-
- //mfaaa = drho;
- //mfaaa = press + (ax+by+cz)/three; // 1/3 = 2/3*(1/op-1/2)
- mfaaa = press; // if drho is interpolated directly
-
- vx2 = vvx*vvx;
- vy2 = vvy*vvy;
- vz2 = vvz*vvz;
- oMdrho = c1o1;
- //oMdrho = one - mfaaa;
-
- //two
- // linear combinations
- mxxPyyPzz = mfaaa;
- //mxxMyy = -c2o3*(ax - by)*eps_new/o;
- //mxxMzz = -c2o3*(ax - cz)*eps_new/o;
-
- //mfabb = -c1o3 * (bz + cy)*eps_new/o;
- //mfbab = -c1o3 * (az + cx)*eps_new/o;
- //mfbba = -c1o3 * (ay + bx)*eps_new/o;
- mxxMyy = -c2o3*((ax - by)+kxxMyyAverage)*eps_new/o * (c1o1 + press);
- mxxMzz = -c2o3*((ax - cz)+kxxMzzAverage)*eps_new/o * (c1o1 + press);
-
- mfabb = -c1o3 * ((bz + cy)+kyzAverage)*eps_new/o * (c1o1 + press);
- mfbab = -c1o3 * ((az + cx)+kxzAverage)*eps_new/o * (c1o1 + press);
- mfbba = -c1o3 * ((ay + bx)+kxyAverage)*eps_new/o * (c1o1 + press);
-
-
- // linear combinations back
- mfcaa = c1o3 * ( mxxMyy + mxxMzz + mxxPyyPzz) * NeqOn;
- mfaca = c1o3 * (-c2o1 * mxxMyy + mxxMzz + mxxPyyPzz) * NeqOn;
- mfaac = c1o3 * ( mxxMyy - c2o1 * mxxMzz + mxxPyyPzz) * NeqOn;
-
- //3.
- // linear combinations
- //residu = residutmp * (ayz + bxz + cxy );
- //mfbbb = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
- mfbbb = c0o1;
-
- //residu = residutmp * (axy + two*bxx + two*bzz + cyz );
- //residu = -(c1o9*(axy - 2*bxx - 2*bzz + cyz ));
- //mxxyPyzz = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
- mxxyPyzz = c0o1;
-
- //residu = residutmp * (axy + two*bxx - two*bzz - cyz );
- //residu = c1o9*(axy - 2*bxx + 2*bzz - cyz );
- //mxxyMyzz = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
- mxxyMyzz = c0o1;
-
- //residu = residutmp * (axz + byz + two*cxx + two*cyy );
- //residu = -(c1o9*(axz + byz - 2*cxx - 2*cyy ));
- //mxxzPyyz = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
- mxxzPyyz = c0o1;
-
- //residu = residutmp * (axz - byz + two*cxx - two*cyy );
- //residu = c1o9*(axz - byz - 2*cxx + 2*cyy );
- //mxxzMyyz = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
- mxxzMyyz = c0o1;
-
- //residu = residutmp * (two*ayy + two*azz + bxy + cxz );
- //residu = c1o9*(2*ayy + 2*azz - bxy - cxz );
- //mxyyPxzz = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
- mxyyPxzz = c0o1;
-
- //residu = residutmp * (two*ayy - two*azz + bxy - cxz );
- //residu = c1o9*(-2*ayy + 2*azz + bxy - cxz );
- //mxyyMxzz = (abs(residu)+qudricLimit) * residu / (qudricLimit * O3 + abs(residu));
- mxyyMxzz = c0o1;
-
- // linear combinations back
- mfcba = ( mxxyMyzz + mxxyPyzz) * c1o2;
- mfabc = (-mxxyMyzz + mxxyPyzz) * c1o2;
- mfcab = ( mxxzMyyz + mxxzPyyz) * c1o2;
- mfacb = (-mxxzMyyz + mxxzPyyz) * c1o2;
- mfbca = ( mxyyMxzz + mxyyPxzz) * c1o2;
- mfbac = (-mxyyMxzz + mxyyPxzz) * c1o2;
-
- //4.
- mfacc = mfaaa*c1o9;
- mfcac = mfacc;
- mfcca = mfacc;
- //5.
-
- //6.
- mfccc = mfaaa*c1o27;
- ////////////////////////////////////////////////////////////////////////////////////
- //back
- ////////////////////////////////////////////////////////////////////////////////////
- //mit 1, 0, 1/3, 0, 0, 0, 1/3, 0, 1/9 Konditionieren
- ////////////////////////////////////////////////////////////////////////////////////
- // Z - Dir
- m0 = mfaac * c1o2 + mfaab * (vvz - c1o2) + (mfaaa + c1o1 * oMdrho) * ( vz2 - vvz) * c1o2;
- m1 = -mfaac - c2o1 * mfaab * vvz + mfaaa * (c1o1 - vz2) - c1o1 * oMdrho * vz2;
- m2 = mfaac * c1o2 + mfaab * (vvz + c1o2) + (mfaaa + c1o1 * oMdrho) * ( vz2 + vvz) * c1o2;
- mfaaa = m0;
- mfaab = m1;
- mfaac = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- m0 = mfabc * c1o2 + mfabb * (vvz - c1o2) + mfaba * ( vz2 - vvz) * c1o2;
- m1 = -mfabc - c2o1 * mfabb * vvz + mfaba * (c1o1 - vz2);
- m2 = mfabc * c1o2 + mfabb * (vvz + c1o2) + mfaba * ( vz2 + vvz) * c1o2;
- mfaba = m0;
- mfabb = m1;
- mfabc = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- m0 = mfacc * c1o2 + mfacb * (vvz - c1o2) + (mfaca + c1o3 * oMdrho) * ( vz2 - vvz) * c1o2;
- m1 = -mfacc - c2o1 * mfacb * vvz + mfaca * (c1o1 - vz2) - c1o3 * oMdrho * vz2;
- m2 = mfacc * c1o2 + mfacb * (vvz + c1o2) + (mfaca + c1o3 * oMdrho) * ( vz2 + vvz) * c1o2;
- mfaca = m0;
- mfacb = m1;
- mfacc = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- ////////////////////////////////////////////////////////////////////////////////////
- m0 = mfbac * c1o2 + mfbab * (vvz - c1o2) + mfbaa * ( vz2 - vvz) * c1o2;
- m1 = -mfbac - c2o1 * mfbab * vvz + mfbaa * (c1o1 - vz2);
- m2 = mfbac * c1o2 + mfbab * (vvz + c1o2) + mfbaa * ( vz2 + vvz) * c1o2;
- mfbaa = m0;
- mfbab = m1;
- mfbac = m2;
- /////////b//////////////////////////////////////////////////////////////////////////
- m0 = mfbbc * c1o2 + mfbbb * (vvz - c1o2) + mfbba * ( vz2 - vvz) * c1o2;
- m1 = -mfbbc - c2o1 * mfbbb * vvz + mfbba * (c1o1 - vz2);
- m2 = mfbbc * c1o2 + mfbbb * (vvz + c1o2) + mfbba * ( vz2 + vvz) * c1o2;
- mfbba = m0;
- mfbbb = m1;
- mfbbc = m2;
- /////////b//////////////////////////////////////////////////////////////////////////
- m0 = mfbcc * c1o2 + mfbcb * (vvz - c1o2) + mfbca * ( vz2 - vvz) * c1o2;
- m1 = -mfbcc - c2o1 * mfbcb * vvz + mfbca * (c1o1 - vz2);
- m2 = mfbcc * c1o2 + mfbcb * (vvz + c1o2) + mfbca * ( vz2 + vvz) * c1o2;
- mfbca = m0;
- mfbcb = m1;
- mfbcc = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- ////////////////////////////////////////////////////////////////////////////////////
- m0 = mfcac * c1o2 + mfcab * (vvz - c1o2) + (mfcaa + c1o3 * oMdrho) * ( vz2 - vvz) * c1o2;
- m1 = -mfcac - c2o1 * mfcab * vvz + mfcaa * (c1o1 - vz2) - c1o3 * oMdrho * vz2;
- m2 = mfcac * c1o2 + mfcab * (vvz + c1o2) + (mfcaa + c1o3 * oMdrho) * ( vz2 + vvz) * c1o2;
- mfcaa = m0;
- mfcab = m1;
- mfcac = m2;
- /////////c//////////////////////////////////////////////////////////////////////////
- m0 = mfcbc * c1o2 + mfcbb * (vvz - c1o2) + mfcba * ( vz2 - vvz) * c1o2;
- m1 = -mfcbc - c2o1 * mfcbb * vvz + mfcba * (c1o1 - vz2);
- m2 = mfcbc * c1o2 + mfcbb * (vvz + c1o2) + mfcba * ( vz2 + vvz) * c1o2;
- mfcba = m0;
- mfcbb = m1;
- mfcbc = m2;
- /////////c//////////////////////////////////////////////////////////////////////////
- m0 = mfccc * c1o2 + mfccb * (vvz - c1o2) + (mfcca + c1o9 * oMdrho) * ( vz2 - vvz) * c1o2;
- m1 = -mfccc - c2o1 * mfccb * vvz + mfcca * (c1o1 - vz2) - c1o9 * oMdrho * vz2;
- m2 = mfccc * c1o2 + mfccb * (vvz + c1o2) + (mfcca + c1o9 * oMdrho) * ( vz2 + vvz) * c1o2;
- mfcca = m0;
- mfccb = m1;
- mfccc = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- ////////////////////////////////////////////////////////////////////////////////////
- //mit 1/6, 2/3, 1/6, 0, 0, 0, 1/18, 2/9, 1/18 Konditionieren
- ////////////////////////////////////////////////////////////////////////////////////
- // Y - Dir
- m0 = mfaca * c1o2 + mfaba * (vvy - c1o2) + (mfaaa + c1o6 * oMdrho) * ( vy2 - vvy) * c1o2;
- m1 = -mfaca - c2o1 * mfaba * vvy + mfaaa * (c1o1 - vy2) - c1o6 * oMdrho * vy2;
- m2 = mfaca * c1o2 + mfaba * (vvy + c1o2) + (mfaaa + c1o6 * oMdrho) * ( vy2 + vvy) * c1o2;
- mfaaa = m0;
- mfaba = m1;
- mfaca = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- m0 = mfacb * c1o2 + mfabb * (vvy - c1o2) + (mfaab + c2o3 * oMdrho) * ( vy2 - vvy) * c1o2;
- m1 = -mfacb - c2o1 * mfabb * vvy + mfaab * (c1o1 - vy2) - c2o3 * oMdrho * vy2;
- m2 = mfacb * c1o2 + mfabb * (vvy + c1o2) + (mfaab + c2o3 * oMdrho) * ( vy2 + vvy) * c1o2;
- mfaab = m0;
- mfabb = m1;
- mfacb = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- m0 = mfacc * c1o2 + mfabc * (vvy - c1o2) + (mfaac + c1o6 * oMdrho) * ( vy2 - vvy) * c1o2;
- m1 = -mfacc - c2o1 * mfabc * vvy + mfaac * (c1o1 - vy2) - c1o6 * oMdrho * vy2;
- m2 = mfacc * c1o2 + mfabc * (vvy + c1o2) + (mfaac + c1o6 * oMdrho) * ( vy2 + vvy) * c1o2;
- mfaac = m0;
- mfabc = m1;
- mfacc = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- ////////////////////////////////////////////////////////////////////////////////////
- m0 = mfbca * c1o2 + mfbba * (vvy - c1o2) + mfbaa * ( vy2 - vvy) * c1o2;
- m1 = -mfbca - c2o1 * mfbba * vvy + mfbaa * (c1o1 - vy2);
- m2 = mfbca * c1o2 + mfbba * (vvy + c1o2) + mfbaa * ( vy2 + vvy) * c1o2;
- mfbaa = m0;
- mfbba = m1;
- mfbca = m2;
- /////////b//////////////////////////////////////////////////////////////////////////
- m0 = mfbcb * c1o2 + mfbbb * (vvy - c1o2) + mfbab * ( vy2 - vvy) * c1o2;
- m1 = -mfbcb - c2o1 * mfbbb * vvy + mfbab * (c1o1 - vy2);
- m2 = mfbcb * c1o2 + mfbbb * (vvy + c1o2) + mfbab * ( vy2 + vvy) * c1o2;
- mfbab = m0;
- mfbbb = m1;
- mfbcb = m2;
- /////////b//////////////////////////////////////////////////////////////////////////
- m0 = mfbcc * c1o2 + mfbbc * (vvy - c1o2) + mfbac * ( vy2 - vvy) * c1o2;
- m1 = -mfbcc - c2o1 * mfbbc * vvy + mfbac * (c1o1 - vy2);
- m2 = mfbcc * c1o2 + mfbbc * (vvy + c1o2) + mfbac * ( vy2 + vvy) * c1o2;
- mfbac = m0;
- mfbbc = m1;
- mfbcc = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- ////////////////////////////////////////////////////////////////////////////////////
- m0 = mfcca * c1o2 + mfcba * (vvy - c1o2) + (mfcaa + c1o18 * oMdrho) * ( vy2 - vvy) * c1o2;
- m1 = -mfcca - c2o1 * mfcba * vvy + mfcaa * (c1o1 - vy2) - c1o18 * oMdrho * vy2;
- m2 = mfcca * c1o2 + mfcba * (vvy + c1o2) + (mfcaa + c1o18 * oMdrho) * ( vy2 + vvy) * c1o2;
- mfcaa = m0;
- mfcba = m1;
- mfcca = m2;
- /////////c//////////////////////////////////////////////////////////////////////////
- m0 = mfccb * c1o2 + mfcbb * (vvy - c1o2) + (mfcab + c2o9 * oMdrho) * ( vy2 - vvy) * c1o2;
- m1 = -mfccb - c2o1 * mfcbb * vvy + mfcab * (c1o1 - vy2) - c2o9 * oMdrho * vy2;
- m2 = mfccb * c1o2 + mfcbb * (vvy + c1o2) + (mfcab + c2o9 * oMdrho) * ( vy2 + vvy) * c1o2;
- mfcab = m0;
- mfcbb = m1;
- mfccb = m2;
- /////////c//////////////////////////////////////////////////////////////////////////
- m0 = mfccc * c1o2 + mfcbc * (vvy - c1o2) + (mfcac + c1o18 * oMdrho) * ( vy2 - vvy) * c1o2;
- m1 = -mfccc - c2o1 * mfcbc * vvy + mfcac * (c1o1 - vy2) - c1o18 * oMdrho * vy2;
- m2 = mfccc * c1o2 + mfcbc * (vvy + c1o2) + (mfcac + c1o18 * oMdrho) * ( vy2 + vvy) * c1o2;
- mfcac = m0;
- mfcbc = m1;
- mfccc = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- ////////////////////////////////////////////////////////////////////////////////////
- //mit 1/36, 1/9, 1/36, 1/9, 4/9, 1/9, 1/36, 1/9, 1/36 Konditionieren
- ////////////////////////////////////////////////////////////////////////////////////
- // X - Dir
- m0 = mfcaa * c1o2 + mfbaa * (vvx - c1o2) + (mfaaa + c1o36 * oMdrho) * ( vx2 - vvx) * c1o2;
- m1 = -mfcaa - c2o1 * mfbaa * vvx + mfaaa * (c1o1 - vx2) - c1o36 * oMdrho * vx2;
- m2 = mfcaa * c1o2 + mfbaa * (vvx + c1o2) + (mfaaa + c1o36 * oMdrho) * ( vx2 + vvx) * c1o2;
- mfaaa = m0;
- mfbaa = m1;
- mfcaa = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- m0 = mfcba * c1o2 + mfbba * (vvx - c1o2) + (mfaba + c1o9 * oMdrho) * ( vx2 - vvx) * c1o2;
- m1 = -mfcba - c2o1 * mfbba * vvx + mfaba * (c1o1 - vx2) - c1o9 * oMdrho * vx2;
- m2 = mfcba * c1o2 + mfbba * (vvx + c1o2) + (mfaba + c1o9 * oMdrho) * ( vx2 + vvx) * c1o2;
- mfaba = m0;
- mfbba = m1;
- mfcba = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- m0 = mfcca * c1o2 + mfbca * (vvx - c1o2) + (mfaca + c1o36 * oMdrho) * ( vx2 - vvx) * c1o2;
- m1 = -mfcca - c2o1 * mfbca * vvx + mfaca * (c1o1 - vx2) - c1o36 * oMdrho * vx2;
- m2 = mfcca * c1o2 + mfbca * (vvx + c1o2) + (mfaca + c1o36 * oMdrho) * ( vx2 + vvx) * c1o2;
- mfaca = m0;
- mfbca = m1;
- mfcca = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- ////////////////////////////////////////////////////////////////////////////////////
- m0 = mfcab * c1o2 + mfbab * (vvx - c1o2) + (mfaab + c1o9 * oMdrho) * ( vx2 - vvx) * c1o2;
- m1 = -mfcab - c2o1 * mfbab * vvx + mfaab * (c1o1 - vx2) - c1o9 * oMdrho * vx2;
- m2 = mfcab * c1o2 + mfbab * (vvx + c1o2) + (mfaab + c1o9 * oMdrho) * ( vx2 + vvx) * c1o2;
- mfaab = m0;
- mfbab = m1;
- mfcab = m2;
- ///////////b////////////////////////////////////////////////////////////////////////
- m0 = mfcbb * c1o2 + mfbbb * (vvx - c1o2) + (mfabb + c4o9 * oMdrho) * ( vx2 - vvx) * c1o2;
- m1 = -mfcbb - c2o1 * mfbbb * vvx + mfabb * (c1o1 - vx2) - c4o9 * oMdrho * vx2;
- m2 = mfcbb * c1o2 + mfbbb * (vvx + c1o2) + (mfabb + c4o9 * oMdrho) * ( vx2 + vvx) * c1o2;
- mfabb = m0;
- mfbbb = m1;
- mfcbb = m2;
- ///////////b////////////////////////////////////////////////////////////////////////
- m0 = mfccb * c1o2 + mfbcb * (vvx - c1o2) + (mfacb + c1o9 * oMdrho) * ( vx2 - vvx) * c1o2;
- m1 = -mfccb - c2o1 * mfbcb * vvx + mfacb * (c1o1 - vx2) - c1o9 * oMdrho * vx2;
- m2 = mfccb * c1o2 + mfbcb * (vvx + c1o2) + (mfacb + c1o9 * oMdrho) * ( vx2 + vvx) * c1o2;
- mfacb = m0;
- mfbcb = m1;
- mfccb = m2;
- ////////////////////////////////////////////////////////////////////////////////////
- ////////////////////////////////////////////////////////////////////////////////////
- m0 = mfcac * c1o2 + mfbac * (vvx - c1o2) + (mfaac + c1o36 * oMdrho) * ( vx2 - vvx) * c1o2;
- m1 = -mfcac - c2o1 * mfbac * vvx + mfaac * (c1o1 - vx2) - c1o36 * oMdrho * vx2;
- m2 = mfcac * c1o2 + mfbac * (vvx + c1o2) + (mfaac + c1o36 * oMdrho) * ( vx2 + vvx) * c1o2;
- mfaac = m0;
- mfbac = m1;
- mfcac = m2;
- ///////////c////////////////////////////////////////////////////////////////////////
- m0 = mfcbc * c1o2 + mfbbc * (vvx - c1o2) + (mfabc + c1o9 * oMdrho) * ( vx2 - vvx) * c1o2;
- m1 = -mfcbc - c2o1 * mfbbc * vvx + mfabc * (c1o1 - vx2) - c1o9 * oMdrho * vx2;
- m2 = mfcbc * c1o2 + mfbbc * (vvx + c1o2) + (mfabc + c1o9 * oMdrho) * ( vx2 + vvx) * c1o2;
- mfabc = m0;
- mfbbc = m1;
- mfcbc = m2;
- ///////////c////////////////////////////////////////////////////////////////////////
- m0 = mfccc * c1o2 + mfbcc * (vvx - c1o2) + (mfacc + c1o36 * oMdrho) * ( vx2 - vvx) * c1o2;
- m1 = -mfccc - c2o1 * mfbcc * vvx + mfacc * (c1o1 - vx2) - c1o36 * oMdrho * vx2;
- m2 = mfccc * c1o2 + mfbcc * (vvx + c1o2) + (mfacc + c1o36 * oMdrho) * ( vx2 + vvx) * c1o2;
- mfacc = m0;
- mfbcc = m1;
- mfccc = m2;
- ////////////////////////////////////////////////////////////////////////////////////
-
- ////////////////////////////////////////////////////////////////////////////////////
- //index 0
- kzero= posC[k];
- kw = neighborCX[kzero];
- ks = neighborCY[kzero];
- kb = neighborCZ[kzero];
- ksw = neighborCY[kw];
- kbw = neighborCZ[kw];
- kbs = neighborCZ[ks];
- kbsw = neighborCZ[ksw];
- ////////////////////////////////////////////////////////////////////////////////////
-
-
- ////////////////////////////////////////////////////////////////////////////////////
- feC[kzero] = mfcbb;
- fwC[kw] = mfabb;
- fnC[kzero] = mfbcb;
- fsC[ks] = mfbab;
- ftC[kzero] = mfbbc;
- fbC[kb] = mfbba;
- fneC[kzero] = mfccb;
- fswC[ksw] = mfaab;
- fseC[ks] = mfcab;
- fnwC[kw] = mfacb;
- fteC[kzero] = mfcbc;
- fbwC[kbw] = mfaba;
- fbeC[kb] = mfcba;
- ftwC[kw] = mfabc;
- ftnC[kzero] = mfbcc;
- fbsC[kbs] = mfbaa;
- fbnC[kb] = mfbca;
- ftsC[ks] = mfbac;
- fzeroC[kzero] = mfbbb;
- ftneC[kzero] = mfccc;
- ftseC[ks] = mfcac;
- fbneC[kb] = mfcca;
- fbseC[kbs] = mfcaa;
- ftnwC[kw] = mfacc;
- ftswC[ksw] = mfaac;
- fbnwC[kbw] = mfaca;
- fbswC[kbsw] = mfaaa;
- ////////////////////////////////////////////////////////////////////////////////////
- }
+extern "C" __global__ void scaleFC_RhoSq_comp_27_Stream(real *DC, real *DF, unsigned int *neighborCX, unsigned int *neighborCY,
+ unsigned int *neighborCZ, unsigned int *neighborFX,
+ unsigned int *neighborFY, unsigned int *neighborFZ,
+ unsigned int size_MatC, unsigned int size_MatF, bool evenOrOdd,
+ unsigned int *posC, unsigned int *posFSWB, unsigned int kFC,
+ real omCoarse, real omFine, real nu, unsigned int nxC, unsigned int nyC, unsigned int nxF, unsigned int nyF,
+ OffFC offFC, const unsigned int *fluidNodeIndices, unsigned int numberOfFluidNodes)
+{
+ ////////////////////////////////////////////////////////////////////////////////
+ const unsigned ix = threadIdx.x; // Globaler x-Index
+ const unsigned iy = blockIdx.x; // Globaler y-Index
+ const unsigned iz = blockIdx.y; // Globaler z-Index
+
+ const unsigned nx = blockDim.x;
+ const unsigned ny = gridDim.x;
+
+ const unsigned k_thread = nx * (ny * iz + iy) + ix;
+
+ if (k_thread < numberOfFluidNodes) {
+ const unsigned k = fluidNodeIndices[k_thread];
+ scaleFC_RhoSq_comp_27_Calculation(DC, DF, neighborCX, neighborCY, neighborCZ, neighborFX, neighborFY,
+ neighborFZ, size_MatC, size_MatF, evenOrOdd, posC, posFSWB, kFC, omCoarse,
+ omFine, nu, nxC, nyC, nxF, nyF, offFC, k);
+ }
}
-//////////////////////////////////////////////////////////////////////////
-