diff --git a/apps/gpu/LBM/MusselOyster2x/CMakeLists.txt b/apps/gpu/LBM/MusselOyster2x/CMakeLists.txt
index 47f36604929442374a26f0cb7d51014e75c60a32..881b64c4434ae32c1fe4b12cff6de48009e38726 100644
--- a/apps/gpu/LBM/MusselOyster2x/CMakeLists.txt
+++ b/apps/gpu/LBM/MusselOyster2x/CMakeLists.txt
@@ -5,6 +5,6 @@ vf_add_library(BUILDTYPE binary PRIVATE_LINK basics VirtualFluids_GPU GridGenera
set_source_files_properties(MusselOyster2x.cpp PROPERTIES LANGUAGE CUDA)
set_target_properties(MusselOyster2x PROPERTIES
- CUDA_SEPARABLE_COMPILATION ON
- VS_DEBUGGER_COMMAND "C:/Program Files/Microsoft MPI/Bin/mpiexec.exe"
- VS_DEBUGGER_COMMAND_ARGUMENTS "-n 2 \"$<TARGET_FILE:MusselOyster2x>\"")
\ No newline at end of file
+ CUDA_SEPARABLE_COMPILATION ON)
+# VS_DEBUGGER_COMMAND "C:/Program Files/Microsoft MPI/Bin/mpiexec.exe"
+# VS_DEBUGGER_COMMAND_ARGUMENTS "-n 2 \"$<TARGET_FILE:MusselOyster2x>\"")
\ No newline at end of file
diff --git a/src/gpu/VirtualFluids_GPU/Calculation/CalcTurbulenceIntensity.cpp b/src/gpu/VirtualFluids_GPU/Calculation/CalcTurbulenceIntensity.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..e9a986daf3b489827f11d8d2526dc99e08057b97
--- /dev/null
+++ b/src/gpu/VirtualFluids_GPU/Calculation/CalcTurbulenceIntensity.cpp
@@ -0,0 +1,78 @@
+// _ ___ __ __________ _ __ ______________ __
+// | | / (_)____/ /___ ______ _/ / ____/ /_ __(_)___/ /____ / ___/ __ / / / /
+// | | / / / ___/ __/ / / / __ `/ / /_ / / / / / / __ / ___/ / /___/ /_/ / / / /
+// | |/ / / / / /_/ /_/ / /_/ / / __/ / / /_/ / / /_/ (__ ) / /_) / ____/ /__/ /
+// |___/_/_/ \__/\__,_/\__,_/_/_/ /_/\__,_/_/\__,_/____/ \____/_/ \_____/
+//
+//////////////////////////////////////////////////////////////////////////
+#include "Calculation/CalcTurbulenceIntensity.h"
+#include <cuda_runtime.h>
+#include <helper_cuda.h>
+
+void allocTurbulenceIntensity(Parameter *para, CudaMemoryManager *cudaManager, uint size)
+{
+ int lev = para->getFine();
+ cudaManager->cudaAllocTurbulenceIntensity(lev, size);
+ resetTurbulenceIntensity(para, cudaManager, size);
+}
+
+
+void calcVelocityAndFluctuations(Parameter *para, CudaMemoryManager *cudaManager, uint tdiff, uint size)
+{
+ int lev = para->getFine();
+ cudaManager->cudaCopyTurbulenceIntensityDH(lev, size);
+
+ for (uint i = 0; i < size; i++)
+ {
+ // mean velocity
+ para->getParH(lev)->vx_mean[i] = para->getParH(lev)->vx_mean[i] / (real)tdiff;
+ para->getParH(lev)->vy_mean[i] = para->getParH(lev)->vy_mean[i] / (real)tdiff;
+ para->getParH(lev)->vz_mean[i] = para->getParH(lev)->vz_mean[i] / (real)tdiff;
+
+ // fluctuations
+ para->getParH(lev)->vxx[i] = para->getParH(lev)->vxx[i] / (real)tdiff;
+ para->getParH(lev)->vyy[i] = para->getParH(lev)->vyy[i] / (real)tdiff;
+ para->getParH(lev)->vzz[i] = para->getParH(lev)->vzz[i] / (real)tdiff;
+
+ para->getParH(lev)->vxx[i] = para->getParH(lev)->vxx[i] - para->getParH(lev)->vx_mean[i] * para->getParH(lev)->vx_mean[i];
+ para->getParH(lev)->vyy[i] = para->getParH(lev)->vyy[i] - para->getParH(lev)->vy_mean[i] * para->getParH(lev)->vy_mean[i];
+ para->getParH(lev)->vzz[i] = para->getParH(lev)->vzz[i] - para->getParH(lev)->vz_mean[i] * para->getParH(lev)->vz_mean[i];
+ }
+}
+
+
+void calcTurbulenceIntensity(Parameter *para, CudaMemoryManager *cudaManager, uint tdiff, uint size) {
+ calcVelocityAndFluctuations(para, cudaManager, tdiff, size);
+ int lev = para->getFine();
+
+ para->getParH(lev)->turbulenceIntensity.resize(size);
+
+ real fluc_squared;
+ real v_mean_squared;
+ for (uint i = 0; i < size; i++) {
+ fluc_squared =
+ (real)(1.0 / 3.0 * (para->getParH(lev)->vxx[i] + para->getParH(lev)->vyy[i] + para->getParH(lev)->vzz[i]));
+ v_mean_squared = para->getParH(lev)->vx_mean[i] * para->getParH(lev)->vx_mean[i] +
+ para->getParH(lev)->vy_mean[i] * para->getParH(lev)->vy_mean[i] +
+ para->getParH(lev)->vz_mean[i] * para->getParH(lev)->vz_mean[i];
+ para->getParH(lev)->turbulenceIntensity[i] = (real)sqrt(fluc_squared / v_mean_squared);
+ }
+}
+
+
+void resetTurbulenceIntensity(Parameter *para, CudaMemoryManager *cudaManager, uint size)
+{
+ int lev = para->getFine();
+
+ for (unsigned int i = 0; i < size; i++) {
+ para->getParH(lev)->vxx[i] = (real)0.0;
+ para->getParH(lev)->vyy[i] = (real)0.0;
+ para->getParH(lev)->vzz[i] = (real)0.0;
+ para->getParH(lev)->vx_mean[i] = (real)0.0;
+ para->getParH(lev)->vy_mean[i] = (real)0.0;
+ para->getParH(lev)->vz_mean[i] = (real)0.0;
+ }
+
+ cudaManager->cudaCopyTurbulenceIntensityHD(lev, size);
+}
+
diff --git a/src/gpu/VirtualFluids_GPU/Calculation/CalcTurbulenceIntensity.h b/src/gpu/VirtualFluids_GPU/Calculation/CalcTurbulenceIntensity.h
new file mode 100644
index 0000000000000000000000000000000000000000..d9e2575945f48f8890c31952b3d89cd111ada481
--- /dev/null
+++ b/src/gpu/VirtualFluids_GPU/Calculation/CalcTurbulenceIntensity.h
@@ -0,0 +1,14 @@
+#ifndef CalcTurbulenceIntensity_H
+#define CalcTurbulenceIntensity_H
+
+#include "LBM/LB.h"
+#include "GPU/GPU_Interface.h"
+#include "Parameter/Parameter.h"
+#include "GPU/CudaMemoryManager.h"
+
+extern "C" void allocTurbulenceIntensity(Parameter *para, CudaMemoryManager *cudaManager, uint size);
+extern "C" void calcVelocityAndFluctuations(Parameter *para, CudaMemoryManager *cudaManager, uint tdiff, uint size);
+extern "C" void calcTurbulenceIntensity(Parameter *para, CudaMemoryManager *cudaManager, uint tdiff, uint size);
+extern "C" void resetTurbulenceIntensity(Parameter *para, CudaMemoryManager *cudaManager, uint size);
+
+#endif
diff --git a/src/gpu/VirtualFluids_GPU/Calculation/UpdateGrid27.cpp b/src/gpu/VirtualFluids_GPU/Calculation/UpdateGrid27.cpp
index 01f8b4bdaf908ef5caf06f5f8305143b8178a7d0..b3df5ba4925b708ad97eb046753ac5178425114f 100644
--- a/src/gpu/VirtualFluids_GPU/Calculation/UpdateGrid27.cpp
+++ b/src/gpu/VirtualFluids_GPU/Calculation/UpdateGrid27.cpp
@@ -45,6 +45,27 @@ void UpdateGrid27::updateGrid(Parameter *para, vf::gpu::Communicator *comm, Cuda
if( level != para->getFine() )
{
this->refinementAndExchange(para, level, comm, cudaManager);
+ } else {
+ if (para->getCalcTurbulenceIntensity()) {
+ CalcTurbulenceIntensityDevice(
+ para->getParD(level)->vxx,
+ para->getParD(level)->vyy,
+ para->getParD(level)->vzz,
+ para->getParD(level)->vx_mean,
+ para->getParD(level)->vy_mean,
+ para->getParD(level)->vz_mean,
+ para->getParD(level)->d0SP.f[0],
+ para->getParD(level)->QGeom.k,
+ para->getParD(level)->QGeom.kQ,
+ para->getParD(level)->omega,
+ para->getParD(level)->neighborX_SP,
+ para->getParD(level)->neighborY_SP,
+ para->getParD(level)->neighborZ_SP,
+ para->getParD(level)->size_Mat_SP,
+ para->getParD(level)->evenOrOdd,
+ para->getParD(level)->numberofthreads
+ );
+ }
}
}
diff --git a/src/gpu/VirtualFluids_GPU/DataStructureInitializer/GridReaderGenerator/GridGenerator.cpp b/src/gpu/VirtualFluids_GPU/DataStructureInitializer/GridReaderGenerator/GridGenerator.cpp
index ab55abcec171b3bc135f367c23ae2ae877b8e239..e5c5cadba0b26d15ea7460d52efce603a895e93d 100644
--- a/src/gpu/VirtualFluids_GPU/DataStructureInitializer/GridReaderGenerator/GridGenerator.cpp
+++ b/src/gpu/VirtualFluids_GPU/DataStructureInitializer/GridReaderGenerator/GridGenerator.cpp
@@ -1380,7 +1380,7 @@ void GridGenerator::allocArrays_OffsetScale()
// split fine-to-coarse indices into border and bulk
splitFineToCoarseIntoBorderAndBulk(level);
// split coarse-to-fine indices into border and bulk
- splitCoarseToFineIntoBorderAndBulk(level);
+ //splitCoarseToFineIntoBorderAndBulk(level);
}
////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
//copy
diff --git a/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.cpp b/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.cpp
index d1f9eaaa002a761275c1c7b4e7ee3ddd7ee6d2a5..9918386bf18c19adc017082dca39673a5e8d893e 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.cpp
+++ b/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.cpp
@@ -971,6 +971,59 @@ void CudaMemoryManager::cudaFreeTurbulentViscosity(int lev)
checkCudaErrors(cudaFreeHost(parameter->getParH(lev)->gDyvz));
checkCudaErrors(cudaFreeHost(parameter->getParH(lev)->gDzvz));
}
+//turbulence intensity
+void CudaMemoryManager::cudaAllocTurbulenceIntensity(int lev, uint size)
+{
+ uint mem_size = sizeof(real) * size;
+ // Host
+ checkCudaErrors( cudaMallocHost((void**) &(parameter->getParH(lev)->vxx ), mem_size));
+ checkCudaErrors( cudaMallocHost((void**) &(parameter->getParH(lev)->vyy ), mem_size));
+ checkCudaErrors( cudaMallocHost((void**) &(parameter->getParH(lev)->vzz ), mem_size));
+ checkCudaErrors( cudaMallocHost((void**) &(parameter->getParH(lev)->vx_mean ), mem_size));
+ checkCudaErrors( cudaMallocHost((void**) &(parameter->getParH(lev)->vy_mean ), mem_size));
+ checkCudaErrors( cudaMallocHost((void**) &(parameter->getParH(lev)->vz_mean ), mem_size));
+ //Device
+ checkCudaErrors( cudaMalloc((void**) &(parameter->getParD(lev)->vxx ), mem_size));
+ checkCudaErrors( cudaMalloc((void**) &(parameter->getParD(lev)->vyy ), mem_size));
+ checkCudaErrors( cudaMalloc((void**) &(parameter->getParD(lev)->vzz ), mem_size));
+ checkCudaErrors( cudaMalloc((void**) &(parameter->getParD(lev)->vx_mean ), mem_size));
+ checkCudaErrors( cudaMalloc((void**) &(parameter->getParD(lev)->vy_mean ), mem_size));
+ checkCudaErrors( cudaMalloc((void**) &(parameter->getParD(lev)->vz_mean ), mem_size));
+ //////////////////////////////////////////////////////////////////////////
+ double tmp = 6. * (double)mem_size;
+ setMemsizeGPU(tmp, false);
+}
+void CudaMemoryManager::cudaCopyTurbulenceIntensityHD(int lev, uint size)
+{
+ uint mem_size = sizeof(real) * size;
+ //copy host to device
+ checkCudaErrors( cudaMemcpy(parameter->getParD(lev)->vxx , parameter->getParH(lev)->vxx , mem_size , cudaMemcpyHostToDevice));
+ checkCudaErrors( cudaMemcpy(parameter->getParD(lev)->vyy , parameter->getParH(lev)->vyy , mem_size , cudaMemcpyHostToDevice));
+ checkCudaErrors( cudaMemcpy(parameter->getParD(lev)->vzz , parameter->getParH(lev)->vzz , mem_size , cudaMemcpyHostToDevice));
+ checkCudaErrors( cudaMemcpy(parameter->getParD(lev)->vx_mean, parameter->getParH(lev)->vx_mean, mem_size , cudaMemcpyHostToDevice));
+ checkCudaErrors( cudaMemcpy(parameter->getParD(lev)->vy_mean, parameter->getParH(lev)->vy_mean, mem_size , cudaMemcpyHostToDevice));
+ checkCudaErrors( cudaMemcpy(parameter->getParD(lev)->vz_mean, parameter->getParH(lev)->vz_mean, mem_size , cudaMemcpyHostToDevice));
+}
+void CudaMemoryManager::cudaCopyTurbulenceIntensityDH(int lev, uint size)
+{
+ uint mem_size = sizeof(real) * size;
+ //copy device to host
+ checkCudaErrors( cudaMemcpy(parameter->getParH(lev)->vxx , parameter->getParD(lev)->vxx , mem_size , cudaMemcpyDeviceToHost));
+ checkCudaErrors( cudaMemcpy(parameter->getParH(lev)->vyy , parameter->getParD(lev)->vyy , mem_size , cudaMemcpyDeviceToHost));
+ checkCudaErrors( cudaMemcpy(parameter->getParH(lev)->vzz , parameter->getParD(lev)->vzz , mem_size , cudaMemcpyDeviceToHost));
+ checkCudaErrors( cudaMemcpy(parameter->getParH(lev)->vx_mean, parameter->getParD(lev)->vx_mean, mem_size , cudaMemcpyDeviceToHost));
+ checkCudaErrors( cudaMemcpy(parameter->getParH(lev)->vy_mean, parameter->getParD(lev)->vy_mean, mem_size , cudaMemcpyDeviceToHost));
+ checkCudaErrors( cudaMemcpy(parameter->getParH(lev)->vz_mean, parameter->getParD(lev)->vz_mean, mem_size , cudaMemcpyDeviceToHost));
+}
+void CudaMemoryManager::cudaFreeTurbulenceIntensity(int lev)
+{
+ checkCudaErrors( cudaFreeHost(parameter->getParH(lev)->vxx ));
+ checkCudaErrors( cudaFreeHost(parameter->getParH(lev)->vyy ));
+ checkCudaErrors( cudaFreeHost(parameter->getParH(lev)->vzz ));
+ checkCudaErrors( cudaFreeHost(parameter->getParH(lev)->vx_mean ));
+ checkCudaErrors( cudaFreeHost(parameter->getParH(lev)->vy_mean ));
+ checkCudaErrors( cudaFreeHost(parameter->getParH(lev)->vz_mean ));
+}
//median
void CudaMemoryManager::cudaAllocMedianSP(int lev)
{
diff --git a/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.h b/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.h
index 0f8008ae491d9f7711f9056648bf1f2e8326e26a..0dd31d91c3d7b8dbe351d51da37d730560c9b5f1 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.h
+++ b/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.h
@@ -135,6 +135,11 @@ public:
void cudaCopyTurbulentViscosityHD(int lev);
void cudaCopyTurbulentViscosityDH(int lev);
void cudaFreeTurbulentViscosity(int lev);
+
+ void cudaAllocTurbulenceIntensity(int lev, uint size);
+ void cudaCopyTurbulenceIntensityHD(int lev, uint size);
+ void cudaCopyTurbulenceIntensityDH(int lev, uint size);
+ void cudaFreeTurbulenceIntensity(int lev);
void cudaAllocMedianSP(int lev);
void cudaCopyMedianSP(int lev);
diff --git a/src/gpu/VirtualFluids_GPU/GPU/GPU_Interface.h b/src/gpu/VirtualFluids_GPU/GPU/GPU_Interface.h
index 5b7a32481bd076c70eb62f7a7fdfe8afa32936a7..0583dcc35c189b1210de99dd136a2ad6c2298546 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/GPU_Interface.h
+++ b/src/gpu/VirtualFluids_GPU/GPU/GPU_Interface.h
@@ -2525,4 +2525,22 @@ extern "C" void generateRandomValuesDevice(curandState* state,
real* randArray,
unsigned int numberOfThreads);
+extern "C" void CalcTurbulenceIntensityDevice(
+ real* vxx,
+ real* vyy,
+ real* vzz,
+ real* vx_mean,
+ real* vy_mean,
+ real* vz_mean,
+ real* DD,
+ int* k_Q,
+ uint sizeQ,
+ real om1,
+ unsigned int* neighborX,
+ unsigned int* neighborY,
+ unsigned int* neighborZ,
+ unsigned int size_Mat,
+ bool evenOrOdd,
+ uint numberOfThreads);
+
#endif
\ No newline at end of file
diff --git a/src/gpu/VirtualFluids_GPU/GPU/GPU_Kernels.cuh b/src/gpu/VirtualFluids_GPU/GPU/GPU_Kernels.cuh
index 288db43e7bcd36dc4d187982b86178d345601094..ecb8a12a39f9df440b6c5323c83b24684d3b69a9 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/GPU_Kernels.cuh
+++ b/src/gpu/VirtualFluids_GPU/GPU/GPU_Kernels.cuh
@@ -2407,6 +2407,22 @@ extern "C" __global__ void initRandom(curandState* state);
extern "C" __global__ void generateRandomValues(curandState* state,
real* randArray);
+extern "C" __global__ void CalcTurbulenceIntensity(
+ real* vxx,
+ real* vyy,
+ real* vzz,
+ real* vx_mean,
+ real* vy_mean,
+ real* vz_mean,
+ real* DD,
+ int* k_Q,
+ uint sizeQ,
+ real om1,
+ unsigned int* neighborX,
+ unsigned int* neighborY,
+ unsigned int* neighborZ,
+ unsigned int size_Mat,
+ bool evenOrOdd);
#endif
\ No newline at end of file
diff --git a/src/gpu/VirtualFluids_GPU/GPU/LBMKernel.cu b/src/gpu/VirtualFluids_GPU/GPU/LBMKernel.cu
index 6e072c6dc490c7671aa81bcfe14df7d27683cd16..956bbec8beddb93d09f978c8c92277c2bd7f4d02 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/LBMKernel.cu
+++ b/src/gpu/VirtualFluids_GPU/GPU/LBMKernel.cu
@@ -7344,6 +7344,59 @@ extern "C" void generateRandomValuesDevice( curandState* state,
generateRandomValues<<< gridQ, threads >>> (state,randArray);
getLastCudaError("generateRandomValues execution failed");
}
+//////////////////////////////////////////////////////////////////////////
+extern "C" void CalcTurbulenceIntensityDevice(
+ real* vxx,
+ real* vyy,
+ real* vzz,
+ real* vx_mean,
+ real* vy_mean,
+ real* vz_mean,
+ real* DD,
+ int* k_Q,
+ uint sizeQ,
+ real om1,
+ unsigned int* neighborX,
+ unsigned int* neighborY,
+ unsigned int* neighborZ,
+ unsigned int size_Mat,
+ bool evenOrOdd,
+ uint numberOfThreads)
+{
+ int Grid = (size_Mat / numberOfThreads)+1;
+ int Grid1, Grid2;
+ if (Grid>512)
+ {
+ Grid1 = 512;
+ Grid2 = (Grid/Grid1)+1;
+ }
+ else
+ {
+ Grid1 = 1;
+ Grid2 = Grid;
+ }
+ dim3 gridQ(Grid1, Grid2);
+ dim3 threads(numberOfThreads, 1, 1 );
+
+ CalcTurbulenceIntensity<<<gridQ, threads>>>(
+ vxx,
+ vyy,
+ vzz,
+ vx_mean,
+ vy_mean,
+ vz_mean,
+ DD,
+ k_Q,
+ sizeQ,
+ om1,
+ neighborX,
+ neighborY,
+ neighborZ,
+ size_Mat,
+ evenOrOdd);
+
+ getLastCudaError("CalcTurbulenceIntensity execution failed");
+}
diff --git a/src/gpu/VirtualFluids_GPU/GPU/TurbulenceIntensity.cu b/src/gpu/VirtualFluids_GPU/GPU/TurbulenceIntensity.cu
new file mode 100644
index 0000000000000000000000000000000000000000..889be31590e0469e69c1f4c65dce7a51e6bdc843
--- /dev/null
+++ b/src/gpu/VirtualFluids_GPU/GPU/TurbulenceIntensity.cu
@@ -0,0 +1,196 @@
+// _ ___ __ __________ _ __ ______________ __
+// | | / (_)____/ /___ ______ _/ / ____/ /_ __(_)___/ /____ / ___/ __ / / / /
+// | | / / / ___/ __/ / / / __ `/ / /_ / / / / / / __ / ___/ / /___/ /_/ / / / /
+// | |/ / / / / /_/ /_/ / /_/ / / __/ / / /_/ / / /_/ (__ ) / /_) / ____/ /__/ /
+// |___/_/_/ \__/\__,_/\__,_/_/_/ /_/\__,_/_/\__,_/____/ \____/_/ \_____/
+//
+//////////////////////////////////////////////////////////////////////////
+
+/* Device code */
+#include "LBM/LB.h"
+#include "LBM/D3Q27.h"
+#include <lbm/constants/NumericConstants.h>
+
+using namespace vf::lbm::constant;
+
+//////////////////////////////////////////////////////////////////////////////
+extern "C" __global__ void CalcTurbulenceIntensity(
+ real* vxx,
+ real* vyy,
+ real* vzz,
+ real* vx_mean,
+ real* vy_mean,
+ real* vz_mean,
+ real* DD,
+ int* k_Q,
+ uint sizeQ,
+ real om1,
+ unsigned int* neighborX,
+ unsigned int* neighborY,
+ unsigned int* neighborZ,
+ unsigned int size_Mat,
+ bool evenOrOdd)
+{
+ Distributions27 D;
+ if (evenOrOdd==true)
+ {
+ D.f[dirE ] = &DD[dirE *size_Mat];
+ D.f[dirW ] = &DD[dirW *size_Mat];
+ D.f[dirN ] = &DD[dirN *size_Mat];
+ D.f[dirS ] = &DD[dirS *size_Mat];
+ D.f[dirT ] = &DD[dirT *size_Mat];
+ D.f[dirB ] = &DD[dirB *size_Mat];
+ D.f[dirNE ] = &DD[dirNE *size_Mat];
+ D.f[dirSW ] = &DD[dirSW *size_Mat];
+ D.f[dirSE ] = &DD[dirSE *size_Mat];
+ D.f[dirNW ] = &DD[dirNW *size_Mat];
+ D.f[dirTE ] = &DD[dirTE *size_Mat];
+ D.f[dirBW ] = &DD[dirBW *size_Mat];
+ D.f[dirBE ] = &DD[dirBE *size_Mat];
+ D.f[dirTW ] = &DD[dirTW *size_Mat];
+ D.f[dirTN ] = &DD[dirTN *size_Mat];
+ D.f[dirBS ] = &DD[dirBS *size_Mat];
+ D.f[dirBN ] = &DD[dirBN *size_Mat];
+ D.f[dirTS ] = &DD[dirTS *size_Mat];
+ D.f[dirZERO] = &DD[dirZERO*size_Mat];
+ D.f[dirTNE ] = &DD[dirTNE *size_Mat];
+ D.f[dirTSW ] = &DD[dirTSW *size_Mat];
+ D.f[dirTSE ] = &DD[dirTSE *size_Mat];
+ D.f[dirTNW ] = &DD[dirTNW *size_Mat];
+ D.f[dirBNE ] = &DD[dirBNE *size_Mat];
+ D.f[dirBSW ] = &DD[dirBSW *size_Mat];
+ D.f[dirBSE ] = &DD[dirBSE *size_Mat];
+ D.f[dirBNW ] = &DD[dirBNW *size_Mat];
+ }
+ else
+ {
+ D.f[dirW ] = &DD[dirE *size_Mat];
+ D.f[dirE ] = &DD[dirW *size_Mat];
+ D.f[dirS ] = &DD[dirN *size_Mat];
+ D.f[dirN ] = &DD[dirS *size_Mat];
+ D.f[dirB ] = &DD[dirT *size_Mat];
+ D.f[dirT ] = &DD[dirB *size_Mat];
+ D.f[dirSW ] = &DD[dirNE *size_Mat];
+ D.f[dirNE ] = &DD[dirSW *size_Mat];
+ D.f[dirNW ] = &DD[dirSE *size_Mat];
+ D.f[dirSE ] = &DD[dirNW *size_Mat];
+ D.f[dirBW ] = &DD[dirTE *size_Mat];
+ D.f[dirTE ] = &DD[dirBW *size_Mat];
+ D.f[dirTW ] = &DD[dirBE *size_Mat];
+ D.f[dirBE ] = &DD[dirTW *size_Mat];
+ D.f[dirBS ] = &DD[dirTN *size_Mat];
+ D.f[dirTN ] = &DD[dirBS *size_Mat];
+ D.f[dirTS ] = &DD[dirBN *size_Mat];
+ D.f[dirBN ] = &DD[dirTS *size_Mat];
+ D.f[dirZERO] = &DD[dirZERO*size_Mat];
+ D.f[dirTNE ] = &DD[dirBSW *size_Mat];
+ D.f[dirTSW ] = &DD[dirBNE *size_Mat];
+ D.f[dirTSE ] = &DD[dirBNW *size_Mat];
+ D.f[dirTNW ] = &DD[dirBSE *size_Mat];
+ D.f[dirBNE ] = &DD[dirTSW *size_Mat];
+ D.f[dirBSW ] = &DD[dirTNE *size_Mat];
+ D.f[dirBSE ] = &DD[dirTNW *size_Mat];
+ D.f[dirBNW ] = &DD[dirTSE *size_Mat];
+ }
+ ////////////////////////////////////////////////////////////////////////////////
+ const unsigned x = threadIdx.x; // Globaler x-Index
+ const unsigned y = blockIdx.x; // Globaler y-Index
+ const unsigned z = blockIdx.y; // Globaler z-Index
+
+ const unsigned nx = blockDim.x;
+ const unsigned ny = gridDim.x;
+
+ const unsigned k = nx*(ny*z + y) + x;
+ //////////////////////////////////////////////////////////////////////////
+
+ if (k < sizeQ) {
+ ////////////////////////////////////////////////////////////////////////////////
+ // index
+ unsigned int KQK = k_Q[k];
+ unsigned int kzero = KQK;
+ unsigned int ke = KQK;
+ unsigned int kw = neighborX[KQK];
+ unsigned int kn = KQK;
+ unsigned int ks = neighborY[KQK];
+ unsigned int kt = KQK;
+ unsigned int kb = neighborZ[KQK];
+ unsigned int ksw = neighborY[kw];
+ unsigned int kne = KQK;
+ unsigned int kse = ks;
+ unsigned int knw = kw;
+ unsigned int kbw = neighborZ[kw];
+ unsigned int kte = KQK;
+ unsigned int kbe = kb;
+ unsigned int ktw = kw;
+ unsigned int kbs = neighborZ[ks];
+ unsigned int ktn = KQK;
+ unsigned int kbn = kb;
+ unsigned int kts = ks;
+ unsigned int ktse = ks;
+ unsigned int kbnw = kbw;
+ unsigned int ktnw = kw;
+ unsigned int kbse = kbs;
+ unsigned int ktsw = ksw;
+ unsigned int kbne = kb;
+ unsigned int ktne = KQK;
+ unsigned int kbsw = neighborZ[ksw];
+ ////////////////////////////////////////////////////////////////////////////////
+ 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_TNE,
+ f_TSW, f_TSE, f_TNW, f_BNE, f_BSW, f_BSE, f_BNW;
+
+ f_W = (D.f[dirE])[ke];
+ f_E = (D.f[dirW])[kw];
+ f_S = (D.f[dirN])[kn];
+ f_N = (D.f[dirS])[ks];
+ f_B = (D.f[dirT])[kt];
+ f_T = (D.f[dirB])[kb];
+ f_SW = (D.f[dirNE])[kne];
+ f_NE = (D.f[dirSW])[ksw];
+ f_NW = (D.f[dirSE])[kse];
+ f_SE = (D.f[dirNW])[knw];
+ f_BW = (D.f[dirTE])[kte];
+ f_TE = (D.f[dirBW])[kbw];
+ f_TW = (D.f[dirBE])[kbe];
+ f_BE = (D.f[dirTW])[ktw];
+ f_BS = (D.f[dirTN])[ktn];
+ f_TN = (D.f[dirBS])[kbs];
+ f_TS = (D.f[dirBN])[kbn];
+ f_BN = (D.f[dirTS])[kts];
+ f_BSW = (D.f[dirTNE])[ktne];
+ f_BNE = (D.f[dirTSW])[ktsw];
+ f_BNW = (D.f[dirTSE])[ktse];
+ f_BSE = (D.f[dirTNW])[ktnw];
+ f_TSW = (D.f[dirBNE])[kbne];
+ f_TNE = (D.f[dirBSW])[kbsw];
+ f_TNW = (D.f[dirBSE])[kbse];
+ f_TSE = (D.f[dirBNW])[kbnw];
+ ////////////////////////////////////////////////////////////////////////////////
+ real vx, vy, vz, drho;
+ drho = f_TSE + f_TNW + f_TNE + f_TSW + f_BSE + f_BNW + f_BNE + f_BSW + f_BN + f_TS + f_TN + f_BS + f_BE + f_TW +
+ f_TE + f_BW + f_SE + f_NW + f_NE + f_SW + f_T + f_B + f_N + f_S + f_E + f_W + ((D.f[dirREST])[kzero]);
+
+ vx = (((f_TSE - f_BNW) - (f_TNW - f_BSE)) + ((f_TNE - f_BSW) - (f_TSW - f_BNE)) +
+ ((f_BE - f_TW) + (f_TE - f_BW)) + ((f_SE - f_NW) + (f_NE - f_SW)) + (f_E - f_W)) /
+ (c1o1 + drho);
+
+ vy = ((-(f_TSE - f_BNW) + (f_TNW - f_BSE)) + ((f_TNE - f_BSW) - (f_TSW - f_BNE)) +
+ ((f_BN - f_TS) + (f_TN - f_BS)) + (-(f_SE - f_NW) + (f_NE - f_SW)) + (f_N - f_S)) /
+ (c1o1 + drho);
+
+ vz = (((f_TSE - f_BNW) + (f_TNW - f_BSE)) + ((f_TNE - f_BSW) + (f_TSW - f_BNE)) +
+ (-(f_BN - f_TS) + (f_TN - f_BS)) + ((f_TE - f_BW) - (f_BE - f_TW)) + (f_T - f_B)) /
+ (c1o1 + drho);
+
+ // compute subtotals:
+ // fluctuations
+ vxx[k] = vxx[k] + vx * vx;
+ vyy[k] = vyy[k] + vy * vy;
+ vzz[k] = vzz[k] + vz * vz;
+
+ // velocity (for mean velocity)
+ vx_mean[k] = vx_mean[k] + vx;
+ vy_mean[k] = vy_mean[k] + vy;
+ vz_mean[k] = vz_mean[k] + vz;
+ }
+}
+////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
diff --git a/src/gpu/VirtualFluids_GPU/LBM/Simulation.cpp b/src/gpu/VirtualFluids_GPU/LBM/Simulation.cpp
index 3c3f159d4886c77dea0ad6b3d35538999e431a9e..c01c265feb7e785698d7d0ad62a0268c9941b800 100644
--- a/src/gpu/VirtualFluids_GPU/LBM/Simulation.cpp
+++ b/src/gpu/VirtualFluids_GPU/LBM/Simulation.cpp
@@ -36,6 +36,7 @@
#include "Calculation/Cp.h"
#include "Calculation/Calc2ndMoments.h"
#include "Calculation/CalcMedian.h"
+#include "Calculation/CalcTurbulenceIntensity.h"
#include "Calculation/ForceCalculations.h"
#include "Calculation/PorousMedia.h"
//////////////////////////////////////////////////////////////////////////
@@ -202,12 +203,20 @@ void Simulation::init(SPtr<Parameter> para, SPtr<GridProvider> gridProvider, std
//////////////////////////////////////////////////////////////////////////
if (para->getCalcMedian())
{
- output << "alloc Calculation for Mean Valus " << "\n";
+ output << "alloc Calculation for Mean Values " << "\n";
if (para->getDiffOn()) allocMedianAD(para.get(), cudaManager.get());
else allocMedian(para.get(), cudaManager.get());
}
+ //////////////////////////////////////////////////////////////////////////
+ // Turbulence Intensity
+ //////////////////////////////////////////////////////////////////////////
+ if (para->getCalcTurbulenceIntensity()) {
+ output << "alloc arrays for calculating Turbulence Intensity " << "\n";
+ allocTurbulenceIntensity(para.get(), cudaManager.get(), para->getParH(para->getFine())->QGeom.kQ);
+ }
+
//////////////////////////////////////////////////////////////////////////
//allocate memory and initialize 2nd, 3rd and higher order moments
//////////////////////////////////////////////////////////////////////////
@@ -402,6 +411,7 @@ void Simulation::run()
ftimeE = 0.0f;
ftimeS = 0.0f;
unsigned int t, t_prev;
+ uint t_turbulenceIntensity = 0;
unsigned int t_MP = 0;
//////////////////////////////////////////////////////////////////////////
para->setStepEnsight(0);
@@ -952,6 +962,14 @@ void Simulation::run()
resetMedian(para.get());
/////////////////////////////////
}
+ if (para->getCalcTurbulenceIntensity())
+ {
+ uint t_diff = t - t_turbulenceIntensity;
+ calcTurbulenceIntensity(para.get(), cudaManager.get(), t_diff, para->getParH(para->getFine())->QGeom.kQ);
+
+ t_turbulenceIntensity = t;
+ resetTurbulenceIntensity(para.get(), cudaManager.get(), para->getParH(para->getFine())->QGeom.kQ);
+ }
////////////////////////////////////////////////////////////////////////
dataWriter->writeTimestep(para, t);
////////////////////////////////////////////////////////////////////////
diff --git a/src/gpu/VirtualFluids_GPU/Parameter/Parameter.cpp b/src/gpu/VirtualFluids_GPU/Parameter/Parameter.cpp
index 78b065bffb8eb626e70c4d04da3a26e2103db388..226c9ad9041494f6af67d18fe106ae9377870683 100644
--- a/src/gpu/VirtualFluids_GPU/Parameter/Parameter.cpp
+++ b/src/gpu/VirtualFluids_GPU/Parameter/Parameter.cpp
@@ -663,7 +663,11 @@ void Parameter::setTStartOut(unsigned int tStartOut)
}
void Parameter::setTimestepOfCoarseLevel(unsigned int timestep)
{
- this->timestep = timestep;
+ this->timestep = timestep;
+}
+void Parameter::setCalcTurbulenceIntensity(bool calcVelocityAndFluctuations)
+{
+ this->calcVelocityAndFluctuations = calcVelocityAndFluctuations;
}
void Parameter::setCalcMedian(bool calcMedian)
{
@@ -1767,7 +1771,11 @@ bool Parameter::getPrintFiles()
}
bool Parameter::getReadGeo()
{
- return ic.readGeo;
+ return ic.readGeo;
+}
+bool Parameter::getCalcTurbulenceIntensity()
+{
+ return this->calcVelocityAndFluctuations;
}
real Parameter::getDiffusivity()
{
diff --git a/src/gpu/VirtualFluids_GPU/Parameter/Parameter.h b/src/gpu/VirtualFluids_GPU/Parameter/Parameter.h
index e468fe15bd09252f374065145545b6c5dc9f7ca9..f854b2bdb6f34802146aacb139abdf07d18da0da 100644
--- a/src/gpu/VirtualFluids_GPU/Parameter/Parameter.h
+++ b/src/gpu/VirtualFluids_GPU/Parameter/Parameter.h
@@ -118,6 +118,11 @@ struct LBMSimulationParameter
real *turbViscosity;
real *gSij, *gSDij, *gDxvx, *gDyvx, *gDzvx, *gDxvy, *gDyvy, *gDzvy, *gDxvz, *gDyvz, *gDzvz; // DebugInformation
+ // turbulence intensity //
+ real *vx_mean, *vy_mean, *vz_mean; // means
+ real *vxx, *vyy, *vzz; // fluctuations
+ std::vector<real> turbulenceIntensity;
+
// macroscopic values//////
real *vx, *vy, *vz, *rho;
real *vx_SP, *vy_SP, *vz_SP, *rho_SP, *press_SP;
@@ -373,6 +378,7 @@ public:
void setTOut(unsigned int tout);
void setTStartOut(unsigned int tStartOut);
void setTimestepOfCoarseLevel(unsigned int timestep);
+ void setCalcTurbulenceIntensity(bool calcVelocityAndFluctuations);
void setCalcMedian(bool calcMedian);
void setCalcDragLift(bool calcDragLift);
void setCalcCp(bool calcCp);
@@ -576,6 +582,7 @@ public:
bool getCompOn();
bool getPrintFiles();
bool getReadGeo();
+ bool getCalcTurbulenceIntensity();
bool getCalcMedian();
bool getCalcDragLift();
bool getCalcCp();
@@ -799,6 +806,7 @@ private:
bool calcCp { false };
bool writeVeloASCII { false };
bool calcPlaneConc { false };
+ bool calcVelocityAndFluctuations{ false };
bool isBodyForce;
int diffMod {27};
int maxlevel {0};