diff --git a/src/gpu/VirtualFluids_GPU/CMakeLists.txt b/src/gpu/VirtualFluids_GPU/CMakeLists.txt
index cd8abaf95b1904e835d47601bb80006ed0651fe9..4cd0d1df5db24385c2b528294bb665d405e73937 100644
--- a/src/gpu/VirtualFluids_GPU/CMakeLists.txt
+++ b/src/gpu/VirtualFluids_GPU/CMakeLists.txt
@@ -27,4 +27,5 @@ if(BUILD_VF_UNIT_TESTS)
set_source_files_properties(Kernel/Kernels/BasicKernels/FluidFlow/Compressible/CumulantK17/CumulantK17Test.cpp PROPERTIES LANGUAGE CUDA)
set_source_files_properties(Kernel/Utilities/DistributionHelperTests.cpp PROPERTIES LANGUAGE CUDA)
set_source_files_properties(Parameter/ParameterTest.cpp PROPERTIES LANGUAGE CUDA)
+ set_source_files_properties(PreCollisionInteractor/ActuatorFarmInlinesTest.cpp PROPERTIES LANGUAGE CUDA)
endif()
diff --git a/src/gpu/VirtualFluids_GPU/PreCollisionInteractor/ActuatorFarm.cu b/src/gpu/VirtualFluids_GPU/PreCollisionInteractor/ActuatorFarm.cu
index 3f3cfe384f513141d4b8bda36337ac0d96dc8aea..0dc3bb78e57cbdf3641c479a9ca27be80625f78a 100644
--- a/src/gpu/VirtualFluids_GPU/PreCollisionInteractor/ActuatorFarm.cu
+++ b/src/gpu/VirtualFluids_GPU/PreCollisionInteractor/ActuatorFarm.cu
@@ -31,6 +31,7 @@
//! \author Henrik Asmuth, Henry Korb
//======================================================================================
#include "ActuatorFarm.h"
+#include "ActuatorFarmInlines.h"
#include <cuda.h>
#include <cuda_runtime.h>
@@ -49,22 +50,6 @@
using namespace vf::basics::constant;
-
-__host__ __device__ __inline__ uint calcNode(uint bladeNode, uint numberOfBladeNodes, uint blade, uint numberOfBlades, uint turbine, uint numberOfTurbines)
-{
-
- return bladeNode+numberOfBladeNodes*(blade+numberOfBlades*turbine);
-}
-
-__host__ __device__ __inline__ void calcTurbineBladeAndBladeNode(uint node, uint& bladeNode, uint numberOfBladeNodes, uint& blade, uint numberOfBlades, uint& turbine, uint numberOfTurbines)
-{
- turbine = node/(numberOfBladeNodes*numberOfBlades);
- uint x_off = turbine*numberOfBladeNodes*numberOfBlades;
- blade = (node - x_off)/numberOfBlades;
- uint y_off = numberOfBladeNodes*blade+x_off;
- bladeNode = (node - y_off);
-}
-
__host__ __device__ __forceinline__ real distSqrd(real distX, real distY, real distZ)
{
return distX*distX+distY*distY+distZ*distZ;
@@ -100,6 +85,24 @@ __host__ __device__ __inline__ void rotateFromGlobalToBlade(
invRotateAboutX3D(azimuth, tmpX, tmpY, tmpZ, bladeCoordX_BF, bladeCoordY_BF, bladeCoordZ_BF);
}
+__host__ __device__ __inline__ void calcCoordinateOfBladeNodeInGridFrameOfReference(uint bladeNodeIndex, uint turbine, real localAzimuth, real yaw,
+ const real* bladeCoordsX, const real* bladeCoordsY, const real* bladeCoordsZ,
+ const real* turbPosX, const real* turbPosY, const real* turbPosZ,
+ real& bladeCoordX_GF, real& bladeCoordY_GF, real& bladeCoordZ_GF)
+{
+
+ const real bladeCoordX_BF = bladeCoordsX[bladeNodeIndex];
+ const real bladeCoordY_BF = bladeCoordsY[bladeNodeIndex];
+ const real bladeCoordZ_BF = bladeCoordsZ[bladeNodeIndex];
+
+ rotateFromBladeToGlobal(bladeCoordX_BF, bladeCoordY_BF, bladeCoordZ_BF, bladeCoordX_GF, bladeCoordY_GF,
+ bladeCoordZ_GF, localAzimuth, yaw);
+
+ bladeCoordX_GF += turbPosX[turbine];
+ bladeCoordY_GF += turbPosY[turbine];
+ bladeCoordZ_GF += turbPosZ[turbine];
+}
+
__global__ void interpolateVelocities(real* gridCoordsX, real* gridCoordsY, real* gridCoordsZ,
uint* neighborsX, uint* neighborsY, uint* neighborsZ, uint* neighborsWSB,
real* vx, real* vy, real* vz,
@@ -116,29 +119,19 @@ __global__ void interpolateVelocities(real* gridCoordsX, real* gridCoordsY, real
//!
const unsigned nodeIndex = vf::gpu::getNodeIndex();
- if(nodeIndex>=numberOfBladeNodes*numberOfBlades*numberOfTurbines) return;
-
- uint turbine, bladeNode, blade;
+ if (nodeIndex >= numberOfBladeNodes * numberOfBlades * numberOfTurbines) return;
- calcTurbineBladeAndBladeNode(nodeIndex, bladeNode, numberOfBladeNodes, blade, numberOfBlades, turbine, numberOfTurbines);
+ uint turbine, bladeNodeIndexOnBlade, blade;
+ calcTurbineBladeAndBladeNode(nodeIndex, bladeNodeIndexOnBlade, numberOfBladeNodes, blade, numberOfBlades, turbine, numberOfTurbines);
- real bladeCoordX_BF = bladeCoordsX[nodeIndex];
- real bladeCoordY_BF = bladeCoordsY[nodeIndex];
- real bladeCoordZ_BF = bladeCoordsZ[nodeIndex];
+ const real localAzimuth = azimuths[turbine] + blade * c2Pi / numberOfBlades;
+ const real yaw = yaws[turbine];
real bladeCoordX_GF, bladeCoordY_GF, bladeCoordZ_GF;
-
- real localAzimuth = azimuths[turbine]+blade*c2Pi/numberOfBlades;
- real yaw = yaws[turbine];
-
-
- rotateFromBladeToGlobal(bladeCoordX_BF, bladeCoordY_BF, bladeCoordZ_BF,
- bladeCoordX_GF, bladeCoordY_GF, bladeCoordZ_GF,
- localAzimuth, yaw);
-
- bladeCoordX_GF += turbPosX[turbine];
- bladeCoordY_GF += turbPosY[turbine];
- bladeCoordZ_GF += turbPosZ[turbine];
+ calcCoordinateOfBladeNodeInGridFrameOfReference(nodeIndex, turbine, localAzimuth, yaw,
+ bladeCoordsX, bladeCoordsY, bladeCoordsZ,
+ turbPosX, turbPosY, turbPosZ,
+ bladeCoordX_GF, bladeCoordY_GF, bladeCoordZ_GF);
uint k, ke, kn, kt;
uint kne, kte, ktn, ktne;
@@ -147,7 +140,7 @@ __global__ void interpolateVelocities(real* gridCoordsX, real* gridCoordsY, real
gridCoordsX, gridCoordsY, gridCoordsZ,
bladeCoordX_GF, bladeCoordY_GF, bladeCoordZ_GF,
neighborsX, neighborsY, neighborsZ, neighborsWSB);
-
+
bladeIndices[nodeIndex] = k;
getNeighborIndicesOfBSW(k, ke, kn, kt, kne, kte, ktn, ktne, neighborsX, neighborsY, neighborsZ);
@@ -171,7 +164,7 @@ __global__ void interpolateVelocities(real* gridCoordsX, real* gridCoordsY, real
localAzimuth, yaw);
bladeVelocitiesX[nodeIndex] = bladeVelX_BF;
- bladeVelocitiesY[nodeIndex] = bladeVelY_BF+omegas[turbine]*bladeCoordZ_BF;
+ bladeVelocitiesY[nodeIndex] = bladeVelY_BF + omegas[turbine] * bladeCoordsZ[nodeIndex];
bladeVelocitiesZ[nodeIndex] = bladeVelZ_BF;
}
@@ -191,7 +184,6 @@ __global__ void applyBodyForces(real* gridCoordsX, real* gridCoordsY, real* grid
if(index>=nIndices) return;
-
uint gridIndex = gridIndices[index];
real gridCoordX_GF = gridCoordsX[gridIndex];
@@ -229,7 +221,7 @@ __global__ void applyBodyForces(real* gridCoordsX, real* gridCoordsY, real* grid
localAzimuth, yaw);
uint node;
- uint nextNode = calcNode(0, numberOfBladeNodes, blade, numberOfBlades, turbine, numberOfTurbines);
+ uint nextNode = calcNodeIndexInBladeArrays(0, numberOfBladeNodes, blade, numberOfBlades, turbine, numberOfTurbines);
real last_z = c0o1;
real current_z = c0o1;
@@ -240,7 +232,7 @@ __global__ void applyBodyForces(real* gridCoordsX, real* gridCoordsY, real* grid
for( uint bladeNode=0; bladeNode<numberOfBladeNodes-1; bladeNode++)
{
node = nextNode;
- nextNode = calcNode(bladeNode+1, numberOfBladeNodes, blade, numberOfBlades, turbine, numberOfTurbines);
+ nextNode = calcNodeIndexInBladeArrays(bladeNode+1, numberOfBladeNodes, blade, numberOfBlades, turbine, numberOfTurbines);
x = bladeCoordsX[node];
y = bladeCoordsY[node];
@@ -254,7 +246,7 @@ __global__ void applyBodyForces(real* gridCoordsX, real* gridCoordsY, real* grid
rotateFromBladeToGlobal(bladeForcesX[node], bladeForcesY[node], bladeForcesZ[node],
forceX_RF, forceY_RF, forceZ_RF,
localAzimuth, yaw);
-
+
gridForceX_RF += forceX_RF*eta;
gridForceY_RF += forceY_RF*eta;
gridForceZ_RF += forceZ_RF*eta;
@@ -276,7 +268,7 @@ __global__ void applyBodyForces(real* gridCoordsX, real* gridCoordsY, real* grid
rotateFromBladeToGlobal(bladeForcesX[node], bladeForcesY[node], bladeForcesZ[node],
forceX_RF, forceY_RF, forceZ_RF,
localAzimuth, yaw);
-
+
gridForceX_RF += forceX_RF*eta;
gridForceY_RF += forceY_RF*eta;
gridForceZ_RF += forceZ_RF*eta;
@@ -319,7 +311,7 @@ void ActuatorFarm::interact(Parameter* para, CudaMemoryManager* cudaMemoryManage
cudaStream_t stream = para->getStreamManager()->getStream(CudaStreamIndex::ActuatorFarm, this->streamIndex);
- if(useHostArrays) cudaMemoryManager->cudaCopyBladeCoordsHtoD(this);
+ if (useHostArrays) cudaMemoryManager->cudaCopyBladeCoordsHtoD(this);
vf::cuda::CudaGrid bladeGrid = vf::cuda::CudaGrid(para->getParH(level)->numberofthreads, this->numberOfGridNodes);
@@ -386,7 +378,7 @@ void ActuatorFarm::calcForcesEllipticWing()
{
for( uint bladeNode=0; bladeNode<this->numberOfBladeNodes; bladeNode++)
{
- uint node = calcNode(bladeNode, this->numberOfBladeNodes, blade, this->numberOfBlades, turbine, this->numberOfTurbines);
+ uint node = calcNodeIndexInBladeArrays(bladeNode, this->numberOfBladeNodes, blade, this->numberOfBlades, turbine, this->numberOfTurbines);
u_rel = this->bladeVelocitiesXH[node];
v_rel = this->bladeVelocitiesYH[node];
@@ -433,7 +425,7 @@ void ActuatorFarm::initTurbineGeometries(CudaMemoryManager* cudaMemoryManager)
{
for(uint node=0; node<this->numberOfBladeNodes; node++)
{
- this->bladeRadiiH[calcNode(node, numberOfBladeNodes, 0, 1, turbine, numberOfTurbines)] = this->preInitBladeRadii[turbine][node];
+ this->bladeRadiiH[calcNodeIndexInBladeArrays(node, numberOfBladeNodes, 0, 1, turbine, numberOfTurbines)] = this->preInitBladeRadii[turbine][node];
}
}
@@ -461,11 +453,11 @@ void ActuatorFarm::initBladeCoords(CudaMemoryManager* cudaMemoryManager)
{
for(uint bladeNode=0; bladeNode<this->numberOfBladeNodes; bladeNode++)
{
- uint node = calcNode(bladeNode, this->numberOfBladeNodes, blade, this->numberOfBlades, turbine, this->numberOfTurbines);
+ uint node = calcNodeIndexInBladeArrays(bladeNode, this->numberOfBladeNodes, blade, this->numberOfBlades, turbine, this->numberOfTurbines);
this->bladeCoordsXH[node] = c0o1;
this->bladeCoordsYH[node] = c0o1;
- this->bladeCoordsZH[node] = this->bladeRadiiH[calcNode(bladeNode, numberOfBladeNodes, 0, 1, turbine, numberOfTurbines)];
+ this->bladeCoordsZH[node] = this->bladeRadiiH[calcNodeIndexInBladeArrays(bladeNode, numberOfBladeNodes, 0, 1, turbine, numberOfTurbines)];
}
}
}
@@ -593,7 +585,8 @@ void ActuatorFarm::setAllBladeForces(real* _bladeForcesX, real* _bladeForcesY, r
std::copy_n(_bladeForcesY, this->numberOfGridNodes, this->bladeForcesYH);
std::copy_n(_bladeForcesZ, this->numberOfGridNodes, this->bladeForcesZH);
-}void ActuatorFarm::setTurbineBladeCoords(uint turbine, real* _bladeCoordsX, real* _bladeCoordsY, real* _bladeCoordsZ)
+}
+void ActuatorFarm::setTurbineBladeCoords(uint turbine, real* _bladeCoordsX, real* _bladeCoordsY, real* _bladeCoordsZ)
{
std::copy_n(_bladeCoordsX, numberOfBladeNodes*numberOfBlades, &this->bladeCoordsXH[turbine*numberOfBladeNodes*numberOfBlades]);
std::copy_n(_bladeCoordsY, numberOfBladeNodes*numberOfBlades, &this->bladeCoordsYH[turbine*numberOfBladeNodes*numberOfBlades]);
diff --git a/src/gpu/VirtualFluids_GPU/PreCollisionInteractor/ActuatorFarmInlines.h b/src/gpu/VirtualFluids_GPU/PreCollisionInteractor/ActuatorFarmInlines.h
new file mode 100644
index 0000000000000000000000000000000000000000..bbe38b7dffb06f6a5b84d601db6549acabbc0cdf
--- /dev/null
+++ b/src/gpu/VirtualFluids_GPU/PreCollisionInteractor/ActuatorFarmInlines.h
@@ -0,0 +1,55 @@
+//=======================================================================================
+// ____ ____ __ ______ __________ __ __ __ __
+// \ \ | | | | | _ \ |___ ___| | | | | / \ | |
+// \ \ | | | | | |_) | | | | | | | / \ | |
+// \ \ | | | | | _ / | | | | | | / /\ \ | |
+// \ \ | | | | | | \ \ | | | \__/ | / ____ \ | |____
+// \ \ | | |__| |__| \__\ |__| \________/ /__/ \__\ |_______|
+// \ \ | | ________________________________________________________________
+// \ \ | | | ______________________________________________________________|
+// \ \| | | | __ __ __ __ ______ _______
+// \ | | |_____ | | | | | | | | | _ \ / _____)
+// \ | | _____| | | | | | | | | | | \ \ \_______
+// \ | | | | |_____ | \_/ | | | | |_/ / _____ |
+// \ _____| |__| |________| \_______/ |__| |______/ (_______/
+//
+// This file is part of VirtualFluids. VirtualFluids is free software: you can
+// redistribute it and/or modify it under the terms of the GNU General Public
+// License as published by the Free Software Foundation, either version 3 of
+// the License, or (at your option) any later version.
+//
+// VirtualFluids is distributed in the hope that it will be useful, but WITHOUT
+// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+// for more details.
+//
+// You should have received a copy of the GNU General Public License along
+// with VirtualFluids (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>.
+//
+//! \file ActuatorFarm.cu
+//! \ingroup PreCollisionInteractor
+//! \author Henrik Asmuth, Henry Korb, Anna Wellmann
+//======================================================================================
+
+#ifndef ACTUATOR_FARM_INLINES
+#define ACTUATOR_FARM_INLINES
+
+#include "basics/DataTypes.h"
+
+__host__ __device__ __inline__ uint calcNodeIndexInBladeArrays(uint bladeNode, uint numberOfNodesPerBlade, uint blade, uint numberOfBlades, uint turbine, uint numberOfTurbines)
+{
+ // see https://git.rz.tu-bs.de/irmb/VirtualFluids_dev/-/merge_requests/248 for visualization
+ return bladeNode + numberOfNodesPerBlade * (blade + numberOfBlades * turbine);
+}
+
+__host__ __device__ __inline__ void calcTurbineBladeAndBladeNode(uint node, uint &bladeNode, uint numberOfNodesPerBlade, uint &blade, uint numberOfBlades, uint &turbine, uint numberOfTurbines)
+{
+ // see https://git.rz.tu-bs.de/irmb/VirtualFluids_dev/-/merge_requests/248 for visualization
+ turbine = node / (numberOfNodesPerBlade * numberOfBlades);
+ uint x_off = turbine * numberOfNodesPerBlade * numberOfBlades;
+ blade = (node - x_off) / numberOfNodesPerBlade;
+ uint y_off = numberOfNodesPerBlade * blade + x_off;
+ bladeNode = node - y_off;
+}
+
+#endif
diff --git a/src/gpu/VirtualFluids_GPU/PreCollisionInteractor/ActuatorFarmInlinesTest.cpp b/src/gpu/VirtualFluids_GPU/PreCollisionInteractor/ActuatorFarmInlinesTest.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..720ff3c3c3217456ace0284e73d4953d08d7baa8
--- /dev/null
+++ b/src/gpu/VirtualFluids_GPU/PreCollisionInteractor/ActuatorFarmInlinesTest.cpp
@@ -0,0 +1,99 @@
+#include <gmock/gmock.h>
+
+#include "ActuatorFarmInlines.h"
+
+TEST(ActuatorFarmInlinesTest, calcNodeIndexInBladeArrays)
+{
+ const uint numberOfNodesPerBlade = 4;
+ const uint numberOfBlades = 3;
+ const uint numberOfTurbines = 2;
+
+ // first node on first blade
+ uint bladeNode = 0;
+ uint blade = 0;
+ uint turbine = 0;
+ auto nodeIndexInBladeArrays = calcNodeIndexInBladeArrays(bladeNode, numberOfNodesPerBlade, blade, numberOfBlades, turbine, numberOfTurbines);
+ EXPECT_THAT(nodeIndexInBladeArrays, testing::Eq(0));
+
+ // last node on first blade
+ bladeNode = 3;
+ blade = 0;
+ turbine = 0;
+ nodeIndexInBladeArrays = calcNodeIndexInBladeArrays(bladeNode, numberOfNodesPerBlade, blade, numberOfBlades, turbine, numberOfTurbines);
+ EXPECT_THAT(nodeIndexInBladeArrays, testing::Eq(3));
+
+ // first node on third blade
+ bladeNode = 0;
+ blade = 2;
+ turbine = 0;
+ nodeIndexInBladeArrays = calcNodeIndexInBladeArrays(bladeNode, numberOfNodesPerBlade, blade, numberOfBlades, turbine, numberOfTurbines);
+ EXPECT_THAT(nodeIndexInBladeArrays, testing::Eq(8));
+
+ // last node on third blade, also last node on first turbine
+ bladeNode = 3;
+ blade = 2;
+ turbine = 0;
+ nodeIndexInBladeArrays = calcNodeIndexInBladeArrays(bladeNode, numberOfNodesPerBlade, blade, numberOfBlades, turbine, numberOfTurbines);
+ EXPECT_THAT(nodeIndexInBladeArrays, testing::Eq(11));
+
+ // first node on second turbine
+ bladeNode = 0;
+ blade = 0;
+ turbine = 1;
+ nodeIndexInBladeArrays = calcNodeIndexInBladeArrays(bladeNode, numberOfNodesPerBlade, blade, numberOfBlades, turbine, numberOfTurbines);
+ EXPECT_THAT(nodeIndexInBladeArrays, testing::Eq(12));
+
+ // last node on second turbine
+ bladeNode = 3;
+ blade = 2;
+ turbine = 1;
+ nodeIndexInBladeArrays = calcNodeIndexInBladeArrays(bladeNode, numberOfNodesPerBlade, blade, numberOfBlades, turbine, numberOfTurbines);
+ EXPECT_THAT(nodeIndexInBladeArrays, testing::Eq(23));
+}
+
+TEST(ActuatorFarmInlinesTest, calcTurbineBladeAndBladeNode)
+{
+ const uint numberOfNodesPerBlade = 4;
+ const uint numberOfBlades = 3;
+ const uint numberOfTurbines = 2;
+
+ uint bladeNode;
+ uint blade;
+ uint turbine;
+
+ uint node = 0; // first node on first blade
+ calcTurbineBladeAndBladeNode(node, bladeNode, numberOfNodesPerBlade, blade, numberOfBlades, turbine, numberOfTurbines);
+ EXPECT_THAT(bladeNode, testing::Eq(0));
+ EXPECT_THAT(blade, testing::Eq(0));
+ EXPECT_THAT(turbine, testing::Eq(0));
+
+ node = 3; // last node on first blade
+ calcTurbineBladeAndBladeNode(node, bladeNode, numberOfNodesPerBlade, blade, numberOfBlades, turbine, numberOfTurbines);
+ EXPECT_THAT(bladeNode, testing::Eq(3));
+ EXPECT_THAT(blade, testing::Eq(0));
+ EXPECT_THAT(turbine, testing::Eq(0));
+
+ node = 8; // first node on third blade
+ calcTurbineBladeAndBladeNode(node, bladeNode, numberOfNodesPerBlade, blade, numberOfBlades, turbine, numberOfTurbines);
+ EXPECT_THAT(bladeNode, testing::Eq(0));
+ EXPECT_THAT(blade, testing::Eq(2));
+ EXPECT_THAT(turbine, testing::Eq(0));
+
+ node = 11; // last node on third blade, also last node on first turbine
+ calcTurbineBladeAndBladeNode(node, bladeNode, numberOfNodesPerBlade, blade, numberOfBlades, turbine, numberOfTurbines);
+ EXPECT_THAT(bladeNode, testing::Eq(3));
+ EXPECT_THAT(blade, testing::Eq(2));
+ EXPECT_THAT(turbine, testing::Eq(0));
+
+ node = 12; // first node on second turbine
+ calcTurbineBladeAndBladeNode(node, bladeNode, numberOfNodesPerBlade, blade, numberOfBlades, turbine, numberOfTurbines);
+ EXPECT_THAT(bladeNode, testing::Eq(0));
+ EXPECT_THAT(blade, testing::Eq(0));
+ EXPECT_THAT(turbine, testing::Eq(1));
+
+ node = 23; // last node on second turbine
+ calcTurbineBladeAndBladeNode(node, bladeNode, numberOfNodesPerBlade, blade, numberOfBlades, turbine, numberOfTurbines);
+ EXPECT_THAT(bladeNode, testing::Eq(3));
+ EXPECT_THAT(blade, testing::Eq(2));
+ EXPECT_THAT(turbine, testing::Eq(1));
+}