diff --git a/apps/gpu/LBM/MusselOyster/MusselOyster.cpp b/apps/gpu/LBM/MusselOyster/MusselOyster.cpp
index 065adf2af96c2e0a31b785b6d1bc62f3116c04f2..f593891f328cb9254857b8e05dbbb4212ec63970 100644
--- a/apps/gpu/LBM/MusselOyster/MusselOyster.cpp
+++ b/apps/gpu/LBM/MusselOyster/MusselOyster.cpp
@@ -120,6 +120,7 @@ void multipleLevel(const std::string& configPath)
bool useStreams = false;
bool useLevels = true;
para->useReducedCommunicationAfterFtoC = false;
+ para->setCalcTurbulenceIntensity(true);
if (para->getNumprocs() == 1) {
useStreams = false;
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/apps/gpu/LBM/MusselOyster2x/MusselOyster2x.cpp b/apps/gpu/LBM/MusselOyster2x/MusselOyster2x.cpp
index ca8cba6c3e9d77743ac105f31b4823af3530f077..c6e309f5ecdbeba85e96e7b7abe2f93e11faeb6f 100644
--- a/apps/gpu/LBM/MusselOyster2x/MusselOyster2x.cpp
+++ b/apps/gpu/LBM/MusselOyster2x/MusselOyster2x.cpp
@@ -110,10 +110,11 @@ void multipleLevel(const std::string& configPath)
////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
bool useGridGenerator = true;
- bool useMultiGPU = true;
- bool useStreams = true;
+ bool useMultiGPU = false;
+ bool useStreams = false;
bool useLevels = true;
para->useReducedCommunicationAfterFtoC = true;
+ para->setCalcTurbulenceIntensity(true);
////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
diff --git a/apps/gpu/LBM/MusselOyster3z/MusselOyster3z.cpp b/apps/gpu/LBM/MusselOyster3z/MusselOyster3z.cpp
index 0bc6d2e872f41e80a910fc910eb4d1b1eb07bd66..224896fb6d82336df9a1522efe2679b3760aeaec 100644
--- a/apps/gpu/LBM/MusselOyster3z/MusselOyster3z.cpp
+++ b/apps/gpu/LBM/MusselOyster3z/MusselOyster3z.cpp
@@ -114,6 +114,8 @@ void multipleLevel(const std::string& configPath)
bool useStreams = true;
bool useLevels = true;
para->useReducedCommunicationAfterFtoC = true;
+ para->setCalcTurbulenceIntensity(true);
+
////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
diff --git a/src/gpu/VirtualFluids_GPU/Calculation/CalcTurbulenceIntensity.cpp b/src/gpu/VirtualFluids_GPU/Calculation/CalcTurbulenceIntensity.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..c33249c920ab9a9907b25fcf99ba0c7e0c4f1f27
--- /dev/null
+++ b/src/gpu/VirtualFluids_GPU/Calculation/CalcTurbulenceIntensity.cpp
@@ -0,0 +1,182 @@
+// _ ___ __ __________ _ __ ______________ __
+// | | / (_)____/ /___ ______ _/ / ____/ /_ __(_)___/ /____ / ___/ __ / / / /
+// | | / / / ___/ __/ / / / __ `/ / /_ / / / / / / __ / ___/ / /___/ /_/ / / / /
+// | |/ / / / / /_/ /_/ / /_/ / / __/ / / /_/ / / /_/ (__ ) / /_) / ____/ /__/ /
+// |___/_/_/ \__/\__,_/\__,_/_/_/ /_/\__,_/_/\__,_/____/ \____/_/ \_____/
+//
+//////////////////////////////////////////////////////////////////////////
+#include "Calculation/CalcTurbulenceIntensity.h"
+#include <cuda_runtime.h>
+#include <helper_cuda.h>
+#include <basics/Core/StringUtilities/StringUtil.h>
+
+void allocTurbulenceIntensity(Parameter *para, CudaMemoryManager *cudaManager)
+{
+ for (int lev=para->getCoarse(); lev <= para->getFine(); lev++) {
+ cudaManager->cudaAllocTurbulenceIntensity(lev, para->getParH(lev)->size_Mat_SP);
+ para->getParH(lev)->turbulenceIntensity.resize(para->getParH(lev)->size_Mat_SP);
+ }
+ resetVelocityFluctuationsAndMeans(para, cudaManager);
+}
+
+
+void calcVelocityAndFluctuations(Parameter *para, CudaMemoryManager *cudaManager, uint tdiff)
+{
+ for (int lev = para->getCoarse(); lev <= para->getFine(); lev++) {
+ cudaManager->cudaCopyTurbulenceIntensityDH(lev, para->getParH(lev)->size_Mat_SP);
+
+ for (uint i = 0; i < para->getParH(lev)->size_Mat_SP; 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)->vxy[i] = para->getParH(lev)->vxy[i] / (real)tdiff;
+ para->getParH(lev)->vxz[i] = para->getParH(lev)->vxz[i] / (real)tdiff;
+ para->getParH(lev)->vyz[i] = para->getParH(lev)->vyz[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];
+ para->getParH(lev)->vxy[i] =
+ para->getParH(lev)->vxy[i] - para->getParH(lev)->vx_mean[i] * para->getParH(lev)->vy_mean[i];
+ para->getParH(lev)->vxz[i] =
+ para->getParH(lev)->vxz[i] - para->getParH(lev)->vx_mean[i] * para->getParH(lev)->vz_mean[i];
+ para->getParH(lev)->vyz[i] =
+ para->getParH(lev)->vyz[i] - para->getParH(lev)->vy_mean[i] * para->getParH(lev)->vz_mean[i];
+ }
+ }
+}
+
+
+void calcTurbulenceIntensity(Parameter *para, CudaMemoryManager *cudaManager, uint tdiff) {
+
+
+ real fluc_squared;
+ real v_mean_squared;
+
+ for (int lev = para->getCoarse(); lev <= para->getFine(); lev++) {
+ calcVelocityAndFluctuations(para, cudaManager, tdiff);
+
+ for (uint i = 0; i < para->getParH(lev)->size_Mat_SP; 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 resetVelocityFluctuationsAndMeans(Parameter *para, CudaMemoryManager *cudaManager)
+{
+ for (int lev = para->getCoarse(); lev <= para->getFine(); lev++) {
+ for (unsigned int i = 0; i < para->getParH(lev)->size_Mat_SP; 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)->vxy[i] = (real)0.0;
+ para->getParH(lev)->vxz[i] = (real)0.0;
+ para->getParH(lev)->vyz[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, para->getParH(lev)->size_Mat_SP);
+ }
+}
+
+void cudaFreeTurbulenceIntensityArrays(Parameter *para, CudaMemoryManager *cudaManager)
+{
+ for (int lev = para->getCoarse(); lev <= para->getFine(); lev++) {
+ cudaManager->cudaFreeTurbulenceIntensity(lev);
+ }
+}
+
+void writeTurbulenceIntensityToFile(Parameter *para, uint timestep)
+{
+ for (int lev = para->getCoarse(); lev <= para->getFine(); lev++) {
+ std::vector<real *> data = { para->getParH(lev)->turbulenceIntensity.data() };
+ std::vector<std::string> datanames = { "ti" };
+ writeTiStuffToFile(para, timestep, para->getParH(lev)->size_Mat_SP, data, datanames);
+ }
+}
+
+void writeVeloFluctuationToFile(Parameter *para, uint timestep)
+{
+ for (int lev = para->getCoarse(); lev <= para->getFine(); lev++) {
+ std::vector<real *> data = { para->getParH(lev)->vxx, para->getParH(lev)->vyy, para->getParH(lev)->vzz };
+ std::vector<std::string> datanames = { "vxx", "vyy", "vzz" };
+ writeTiStuffToFile(para, timestep, para->getParH(lev)->size_Mat_SP, data, datanames);
+ }
+}
+
+void writeVeloMeansToFile(Parameter *para, uint timestep) {
+ for (int lev = para->getCoarse(); lev <= para->getFine(); lev++) {
+ std::vector<real *> data = { para->getParH(lev)->vx_mean,
+ para->getParH(lev)->vy_mean,
+ para->getParH(lev)->vz_mean };
+ std::vector<std::string> datanames = { "vx_mean", "vy_mean", "vz_mean" };
+ writeTiStuffToFile(para, timestep, para->getParH(lev)->size_Mat_SP, data, datanames);
+ }
+}
+
+void writeAllTiDatafToFile(Parameter *para, uint timestep)
+{
+ for (int lev = para->getCoarse(); lev <= para->getFine(); lev++) {
+ std::vector<real *> data = { para->getParH(lev)->vxx,
+ para->getParH(lev)->vyy,
+ para->getParH(lev)->vzz,
+ para->getParH(lev)->vx_mean,
+ para->getParH(lev)->vy_mean,
+ para->getParH(lev)->vz_mean,
+ para->getParH(lev)->turbulenceIntensity.data() };
+ std::vector<std::string> datanames = { "vxx", "vyy", "vzz", "vx_mean", "vy_mean", "vz_mean", "ti" };
+ writeTiStuffToFile(para, timestep, para->getParH(lev)->size_Mat_SP, data, datanames);
+ }
+}
+
+void writeTiStuffToFile(Parameter *para, uint timestep, int sizeOfTiArray, std::vector<real *> &data,
+ std::vector<std::string> &datanames)
+{
+ //////////////////////////////////////////////////////////////////////////
+ // set filename
+ std::string names;
+ std::for_each(datanames.begin(), datanames.end(), [&names](const std::string &s) { return names += "_" + s; });
+ std::string ffname = para->getFName() + StringUtil::toString<int>(para->getMyID()) + "_" +
+ StringUtil::toString<int>(timestep) + names + "_ti.txt";
+ const char *fname = ffname.c_str();
+ ////////////////////////////////////////////////////////////////////////
+ // set ofstream
+ std::ofstream ostr;
+ //////////////////////////////////////////////////////////////////////////
+ // open file
+ ostr.open(fname);
+ //////////////////////////////////////////////////////////////////////////
+ // add header
+ ostr << "index_sp";
+ for (auto name : datanames) ostr << "\t" << name;
+ ostr << std::endl;
+ ////////////////////////////////////////////////////////////////////////
+ // fill file with data
+ for (int i = 0; i < sizeOfTiArray; i++) {
+ ostr << i;
+ for (auto dataset : data)
+ ostr << "\t" << dataset[i];
+ ostr << std::endl;
+ }
+ //////////////////////////////////////////////////////////////////////////
+ // close file
+ ostr.close();
+ //////////////////////////////////////////////////////////////////////////
+}
\ No newline at end of file
diff --git a/src/gpu/VirtualFluids_GPU/Calculation/CalcTurbulenceIntensity.h b/src/gpu/VirtualFluids_GPU/Calculation/CalcTurbulenceIntensity.h
new file mode 100644
index 0000000000000000000000000000000000000000..4a2d539f3ae31f3975d03cbc0ea73dad90c20f73
--- /dev/null
+++ b/src/gpu/VirtualFluids_GPU/Calculation/CalcTurbulenceIntensity.h
@@ -0,0 +1,24 @@
+#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);
+extern "C" void calcVelocityAndFluctuations(Parameter *para, CudaMemoryManager *cudaManager, uint tdiff);
+extern "C" void calcTurbulenceIntensity(Parameter *para, CudaMemoryManager *cudaManager, uint tdiff);
+extern "C" void resetVelocityFluctuationsAndMeans(Parameter *para, CudaMemoryManager *cudaManager);
+extern "C" void cudaFreeTurbulenceIntensityArrays(Parameter *para, CudaMemoryManager *cudaManager);
+
+
+void writeTurbulenceIntensityToFile(Parameter *para, uint timestep);
+void writeVeloFluctuationToFile(Parameter *para, uint timeste);
+void writeVeloMeansToFile(Parameter *para, uint timestep);
+void writeAllTiDatafToFile(Parameter *para, uint timestep);
+
+void writeTiStuffToFile(Parameter *para, uint timestep, int sizeOfTiArray, std::vector<real *> &data,
+ std::vector<std::string> &datanames);
+
+#endif
diff --git a/src/gpu/VirtualFluids_GPU/GPU/CalcMac27.cu b/src/gpu/VirtualFluids_GPU/GPU/CalcMac27.cu
index ce6e034d1c1eee57e062b736cfcea97e07306f3c..42d0a6f4aff4a85d4912d0dcdf73e30b6bee9eb9 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/CalcMac27.cu
+++ b/src/gpu/VirtualFluids_GPU/GPU/CalcMac27.cu
@@ -273,15 +273,15 @@ extern "C" __global__ void LBCalcMacCompSP27(real *vxD, real *vyD, real *vzD, re
if(k >= size_Mat)
return;
- if (!vf::gpu::isValidFluidNode(geoD[k]))
- return;
-
pressD[k] = c0o1;
rhoD[k] = c0o1;
vxD[k] = c0o1;
vyD[k] = c0o1;
vzD[k] = c0o1;
+ if (!vf::gpu::isValidFluidNode(geoD[k]))
+ return;
+
vf::gpu::DistributionWrapper distr_wrapper(distributions, size_Mat, isEvenTimestep, k, neighborX, neighborY,
neighborZ);
const auto &distribution = distr_wrapper.distribution;
diff --git a/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.cpp b/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.cpp
index 8f00dce4b9decd84fecb74d1e7f96849b163ca0f..140927be10537d592ec7f46550d2942caf206503 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.cpp
+++ b/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.cpp
@@ -975,6 +975,74 @@ 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)->vxy ), mem_size));
+ checkCudaErrors( cudaMallocHost((void**) &(parameter->getParH(lev)->vxz ), mem_size));
+ checkCudaErrors( cudaMallocHost((void**) &(parameter->getParH(lev)->vyz ), 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)->vxy ), mem_size));
+ checkCudaErrors( cudaMalloc((void**) &(parameter->getParD(lev)->vxz ), mem_size));
+ checkCudaErrors( cudaMalloc((void**) &(parameter->getParD(lev)->vyz ), 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 = 9. * (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)->vxy , parameter->getParH(lev)->vxy , mem_size , cudaMemcpyHostToDevice));
+ checkCudaErrors( cudaMemcpy(parameter->getParD(lev)->vxz , parameter->getParH(lev)->vxz , mem_size , cudaMemcpyHostToDevice));
+ checkCudaErrors( cudaMemcpy(parameter->getParD(lev)->vyz , parameter->getParH(lev)->vyz , 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)->vxy , parameter->getParD(lev)->vxy , mem_size , cudaMemcpyDeviceToHost));
+ checkCudaErrors( cudaMemcpy(parameter->getParH(lev)->vxz , parameter->getParD(lev)->vxz , mem_size , cudaMemcpyDeviceToHost));
+ checkCudaErrors( cudaMemcpy(parameter->getParH(lev)->vyz , parameter->getParD(lev)->vyz , 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)->vxy ));
+ checkCudaErrors( cudaFreeHost(parameter->getParH(lev)->vxz ));
+ checkCudaErrors( cudaFreeHost(parameter->getParH(lev)->vyz ));
+ 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 7e81ecbe9bd2c2bc39cc9a5dc3b01cfffbc16d1e..3d7d39d85138b2661960cbf822756421380bfd0c 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.h
+++ b/src/gpu/VirtualFluids_GPU/GPU/CudaMemoryManager.h
@@ -139,6 +139,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 c2f01e4f57c25f219d517fdf2a111ecd25f0f794..93e78c82bf149a5b5398f145dbed86cde2302ea5 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/GPU_Interface.h
+++ b/src/gpu/VirtualFluids_GPU/GPU/GPU_Interface.h
@@ -2526,4 +2526,23 @@ extern "C" void generateRandomValuesDevice(curandState* state,
real* randArray,
unsigned int numberOfThreads);
+extern "C" void CalcTurbulenceIntensityDevice(
+ real* vxx,
+ real* vyy,
+ real* vzz,
+ real* vxy,
+ real* vxz,
+ real* vyz,
+ real* vx_mean,
+ real* vy_mean,
+ real* vz_mean,
+ real* DD,
+ uint *typeOfGridNode,
+ 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..c1e4ce45b8e4ccf70f084eb7742fcf9d1bfdb47f 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/GPU_Kernels.cuh
+++ b/src/gpu/VirtualFluids_GPU/GPU/GPU_Kernels.cuh
@@ -2407,6 +2407,23 @@ 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* vxy,
+ real* vxz,
+ real* vyz,
+ real* vx_mean,
+ real* vy_mean,
+ real* vz_mean,
+ real* DD,
+ uint *typeOfGridNode,
+ 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 57194fe3bf36c6357454d1e6cef6e80fd0f80b60..d3bfb5d125028ac94a8bbdc87192351d90617989 100644
--- a/src/gpu/VirtualFluids_GPU/GPU/LBMKernel.cu
+++ b/src/gpu/VirtualFluids_GPU/GPU/LBMKernel.cu
@@ -7345,6 +7345,61 @@ 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* vxy,
+ real* vxz,
+ real* vyz,
+ real* vx_mean,
+ real* vy_mean,
+ real* vz_mean,
+ real* DD,
+ uint* typeOfGridNode,
+ 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,
+ vxy,
+ vxz,
+ vyz,
+ vx_mean,
+ vy_mean,
+ vz_mean,
+ DD,
+ typeOfGridNode,
+ 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..2282411cbc4ecb7bfa2cd306e40f8bb214dead39
--- /dev/null
+++ b/src/gpu/VirtualFluids_GPU/GPU/TurbulenceIntensity.cu
@@ -0,0 +1,79 @@
+// _ ___ __ __________ _ __ ______________ __
+// | | / (_)____/ /___ ______ _/ / ____/ /_ __(_)___/ /____ / ___/ __ / / / /
+// | | / / / ___/ __/ / / / __ `/ / /_ / / / / / / __ / ___/ / /___/ /_/ / / / /
+// | |/ / / / / /_/ /_/ / /_/ / / __/ / / /_/ / / /_/ (__ ) / /_) / ____/ /__/ /
+// |___/_/_/ \__/\__,_/\__,_/_/_/ /_/\__,_/_/\__,_/____/ \____/_/ \_____/
+//
+//////////////////////////////////////////////////////////////////////////
+
+/* Device code */
+#include "LBM/LB.h"
+#include "LBM/D3Q27.h"
+#include <lbm/constants/NumericConstants.h>
+
+using namespace vf::lbm::constant;
+#include "lbm/MacroscopicQuantities.h"
+#include "../Kernel/Utilities/DistributionHelper.cuh"
+
+
+using namespace vf::lbm::constant;
+
+//////////////////////////////////////////////////////////////////////////////
+extern "C" __global__ void CalcTurbulenceIntensity(
+ real* vxx,
+ real* vyy,
+ real* vzz,
+ real* vxy,
+ real* vxz,
+ real* vyz,
+ real* vx_mean,
+ real* vy_mean,
+ real* vz_mean,
+ real *distributions,
+ uint* typeOfGridNode,
+ unsigned int* neighborX,
+ unsigned int* neighborY,
+ unsigned int* neighborZ,
+ unsigned int size_Mat,
+ bool isEvenTimestep)
+{
+
+ //////////////////////////////////////////////////////////////////////////
+
+
+ ////////////////////////////////////////////////////////////////////////////////
+ const unsigned k = vf::gpu::getNodeIndex();
+
+ if (k >= size_Mat)
+ return;
+
+ if (!vf::gpu::isValidFluidNode(typeOfGridNode[k]))
+ return;
+
+ vf::gpu::DistributionWrapper distr_wrapper(distributions, size_Mat, isEvenTimestep, k, neighborX, neighborY,
+ neighborZ);
+ const auto &distribution = distr_wrapper.distribution;
+
+ // analogue to LBCalcMacCompSP27
+ real rho = vf::lbm::getDensity(distribution.f);
+ real vx = vf::lbm::getCompressibleVelocityX1(distribution.f, rho);
+ real vy = vf::lbm::getCompressibleVelocityX2(distribution.f, rho);
+ real vz = vf::lbm::getCompressibleVelocityX3(distribution.f, rho);
+
+
+ // compute subtotals:
+ // fluctuations
+ vxx[k] = vxx[k] + vx * vx;
+ vyy[k] = vyy[k] + vy * vy;
+ vzz[k] = vzz[k] + vz * vz;
+ vxy[k] = vxy[k] + vx * vy;
+ vxz[k] = vxz[k] + vx * vz;
+ vyz[k] = vyz[k] + vy * 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 ddcfa298ef4bcc121470ddd8c94c4d7fd1c7833b..b503a1601e3176ca9c2b0760bcf98c4985ebcb47 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());
+ }
+
//////////////////////////////////////////////////////////////////////////
//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);
@@ -506,6 +516,30 @@ void Simulation::run()
}
}
+
+ if (para->getCalcTurbulenceIntensity()) {
+ for (int lev = para->getCoarse(); lev <= para->getFine(); lev++) {
+ CalcTurbulenceIntensityDevice(
+ para->getParD(lev)->vxx,
+ para->getParD(lev)->vyy,
+ para->getParD(lev)->vzz,
+ para->getParD(lev)->vxy,
+ para->getParD(lev)->vxz,
+ para->getParD(lev)->vyz,
+ para->getParD(lev)->vx_mean,
+ para->getParD(lev)->vy_mean,
+ para->getParD(lev)->vz_mean,
+ para->getParD(lev)->d0SP.f[0],
+ para->getParD(lev)->geoSP,
+ para->getParD(lev)->neighborX_SP,
+ para->getParD(lev)->neighborY_SP,
+ para->getParD(lev)->neighborZ_SP,
+ para->getParD(lev)->size_Mat_SP,
+ para->getParD(lev)->evenOrOdd,
+ para->getParD(lev)->numberofthreads
+ );
+ }
+ }
////////////////////////////////////////////////////////////////////////////////
@@ -952,9 +986,20 @@ void Simulation::run()
resetMedian(para.get());
/////////////////////////////////
}
+ if (para->getCalcTurbulenceIntensity())
+ {
+ uint t_diff = t - t_turbulenceIntensity;
+ calcTurbulenceIntensity(para.get(), cudaManager.get(), t_diff);
+ //writeAllTiDatafToFile(para.get(), t);
+ }
////////////////////////////////////////////////////////////////////////
dataWriter->writeTimestep(para, t);
////////////////////////////////////////////////////////////////////////
+ if (para->getCalcTurbulenceIntensity()) {
+ t_turbulenceIntensity = t;
+ resetVelocityFluctuationsAndMeans(para.get(), cudaManager.get());
+ }
+ ////////////////////////////////////////////////////////////////////////
if (para->getCalcDragLift()) printDragLift(para.get(), cudaManager.get(), t);
////////////////////////////////////////////////////////////////////////
if (para->getCalcParticle()) copyAndPrintParticles(para.get(), cudaManager.get(), t, false);
@@ -1334,6 +1379,10 @@ void Simulation::free()
}
}
//////////////////////////////////////////////////////////////////////////
+ // Turbulence Intensity
+ if (para->getCalcTurbulenceIntensity()) {
+ cudaFreeTurbulenceIntensityArrays(para.get(), cudaManager.get());
+ }
delete comm;
diff --git a/src/gpu/VirtualFluids_GPU/Output/FileWriter.cpp b/src/gpu/VirtualFluids_GPU/Output/FileWriter.cpp
index 0a9e34166652613659b5b1d609bc40d07f6cee76..06694ca352249ed2ec0907504f21c66a3164fc32 100644
--- a/src/gpu/VirtualFluids_GPU/Output/FileWriter.cpp
+++ b/src/gpu/VirtualFluids_GPU/Output/FileWriter.cpp
@@ -183,6 +183,16 @@ void FileWriter::writeUnstrucuredGridLT(std::shared_ptr<Parameter> para, int lev
nodedatanames.push_back("geo");
//nodedatanames.push_back("sendNodes");
//nodedatanames.push_back("sparseIndex");
+
+ uint firstTurbNode = nodedatanames.size();
+ if (para->getCalcTurbulenceIntensity()) {
+ nodedatanames.push_back("vxx");
+ nodedatanames.push_back("vyy");
+ nodedatanames.push_back("vzz");
+ nodedatanames.push_back("vxy");
+ nodedatanames.push_back("vxz");
+ nodedatanames.push_back("vyz");
+ }
unsigned int number1, number2, number3, number4, number5, number6, number7, number8;
uint dn1, dn2, dn3, dn4, dn5, dn6, dn7, dn8;
bool neighborsAreFluid;
@@ -229,10 +239,20 @@ void FileWriter::writeUnstrucuredGridLT(std::shared_ptr<Parameter> para, int lev
nodedata[4][dn1] = (double)para->getParH(level)->vz_SP[pos] * (double)para->getVelocityRatio();
nodedata[5][dn1] = (double)para->getParH(level)->geoSP[pos];
+ //nodedata[6][dn1] = (double) pos;
+
//int sendNode = 0; // 0 - not a sendNode; 1 - sendNode; 2 - sendNode in communication after fine to coarse
// testForSendNodeZ(para, level, pos, sendNode); // slow and should not be done multiple times --> use for debugging only!
//nodedata[6][dn1] = (double) sendNode;
- // nodedata[7][dn1] = (double) pos;
+
+ if (para->getCalcTurbulenceIntensity()) {
+ nodedata[firstTurbNode ][dn1] = (double)para->getParH(level)->vxx[pos];
+ nodedata[firstTurbNode + 1][dn1] = (double)para->getParH(level)->vyy[pos];
+ nodedata[firstTurbNode + 2][dn1] = (double)para->getParH(level)->vzz[pos];
+ nodedata[firstTurbNode + 3][dn1] = (double)para->getParH(level)->vxy[pos];
+ nodedata[firstTurbNode + 4][dn1] = (double)para->getParH(level)->vxz[pos];
+ nodedata[firstTurbNode + 5][dn1] = (double)para->getParH(level)->vyz[pos];
+ }
//////////////////////////////////////////////////////////////////////////
number2 = para->getParH(level)->neighborX_SP[number1];
@@ -641,5 +661,3 @@ void FileWriter::writeUnstrucuredGridMedianLTConc(std::shared_ptr<Parameter> par
-
-
diff --git a/src/gpu/VirtualFluids_GPU/Parameter/Parameter.cpp b/src/gpu/VirtualFluids_GPU/Parameter/Parameter.cpp
index ab810a5e38740e1d9b6953f2bbec85864276e482..cc7d3b1941e3d2ecbf885b6d239709bb158168a8 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 2fc896e3e18e6c4e9c72163a748a040e4a826a0e..0ed19933314560f7f32f814d03e764e52d0d8432 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, *vxy, *vxz, *vyz; // fluctuations
+ std::vector<real> turbulenceIntensity;
+
// macroscopic values//////
real *vx, *vy, *vz, *rho;
real *vx_SP, *vy_SP, *vz_SP, *rho_SP, *press_SP;
@@ -382,6 +387,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);
@@ -585,6 +591,7 @@ public:
bool getCompOn();
bool getPrintFiles();
bool getReadGeo();
+ bool getCalcTurbulenceIntensity();
bool getCalcMedian();
bool getCalcDragLift();
bool getCalcCp();
@@ -808,6 +815,7 @@ private:
bool calcCp { false };
bool writeVeloASCII { false };
bool calcPlaneConc { false };
+ bool calcVelocityAndFluctuations{ false };
bool isBodyForce;
int diffMod {27};
int maxlevel {0};