2 #include <thrust/device_vector.h> 16 template<alsfvm::boundary::Type BoundaryType,
class PowerClass>
20 int nx,
int ny,
int nz,
int ngx,
int ngy,
int ngz,
21 real p,
int dimensions) {
22 const int index = threadIdx.x + blockIdx.x * blockDim.x;
24 if (index >= nx * ny * nz) {
28 const int i = index %
nx;
29 const int j = (index /
nx) % ny;
30 const int k = (index / nx /
ny);
34 forEachPointInComputeStructureCube<BoundaryType>([&] (
double u,
double u_h) {
35 output[index] += PowerClass::power(fabs(u - u_h), p);
36 }, input, i, j, k, h,
nx,
ny,
nz, ngx, ngy, ngz, dimensions);
41 template<alsfvm::boundary::Type BoundaryType,
class PowerClass,
class BufferClass>
46 size_t numberOfH,
double p) {
48 auto inputView = input[var]->getView();
50 auto outputMemory = output[var]->getHostMemory();
51 auto outputView = outputMemory->getView();
61 buffer.resize(nx * ny * nz);
64 for (
int h = 1; h < numberOfH; ++h) {
66 const int threads = (p < 6 && p == int(p)) ? 1024 : 512;
67 const int size = nx * ny *
nz;
68 const int blockNumber = (size + threads - 1) / threads;
71 computeStructureCubeKernel<BoundaryType, PowerClass>
72 <<< blockNumber, threads>>>
73 (thrust::raw_pointer_cast(
76 h, nx, ny, nz, ngx, ngy, ngz, p, dimensions);
79 real structureResult = thrust::reduce(buffer.begin(),
81 0.0, thrust::plus<real>());
83 outputView.at(h) += structureResult / (nx * ny *
nz);
86 if (!output[var]->isOnHost()) {
87 output[var]->copyFrom(*outputMemory);
94 template<alsfvm::boundary::Type BoundaryType>
97 int numberOfH,
double p) {
99 computeStructureCubeCUDA < BoundaryType, alsutils::math::FastPower<1>>
100 (output, input, buffer, numberOfH,
p);
101 }
else if (p == 2.0) {
102 computeStructureCubeCUDA <BoundaryType, alsutils::math::FastPower<2>>
103 (output, input, buffer, numberOfH,
p);
107 computeStructureCubeCUDA < BoundaryType, alsutils::math::FastPower<3>>
108 (output, input, buffer, numberOfH,
p);
109 }
else if (p == 4.0) {
110 computeStructureCubeCUDA < BoundaryType, alsutils::math::FastPower<4>>
111 (output, input, buffer, numberOfH,
p);
112 }
else if (p == 5.0) {
113 computeStructureCubeCUDA < BoundaryType, alsutils::math::FastPower<5>>
114 (output, input, buffer, numberOfH,
p);
116 computeStructureCubeCUDA < BoundaryType, alsutils::math::PowPower>
117 (output, input, buffer, numberOfH,
p);
virtual size_t getNumberOfYCells() const
Definition: Volume.cpp:208
The Volume class represents a volume (a collection of cells with values for each cell (eg...
Definition: Volume.hpp:30
virtual size_t getNumberOfVariables() const
getNumberOfVariables gets the number of variables used
Definition: Volume.cpp:83
size_t nx
Definition: VolumeFactory.cpp:87
double real
Definition: types.hpp:65
void dispatchComputeStructureCubeCUDA(alsfvm::volume::Volume &output, const alsfvm::volume::Volume &input, thrust::device_vector< real > &buffer, int numberOfH, double p)
Definition: structure_common_cuda.hpp:95
virtual size_t getNumberOfZCells() const
Definition: Volume.cpp:215
virtual size_t getNumberOfXGhostCells() const
Definition: Volume.cpp:238
__global__ void computeStructureCubeKernel(real *output, alsfvm::memory::View< const real > input, int h, int nx, int ny, int nz, int ngx, int ngy, int ngz, real p, int dimensions)
Definition: structure_common_cuda.hpp:17
virtual size_t getNumberOfZGhostCells() const
Definition: Volume.cpp:254
size_t ny
Definition: VolumeFactory.cpp:88
virtual size_t getNumberOfYGhostCells() const
Definition: Volume.cpp:246
virtual size_t getNumberOfXCells() const
Definition: Volume.cpp:201
virtual size_t getDimensions() const
Gets the number of space dimensions.
Definition: Volume.cpp:279
size_t nz
Definition: VolumeFactory.cpp:89
void computeStructureCubeCUDA(alsfvm::volume::Volume &output, const alsfvm::volume::Volume &input, BufferClass &buffer, size_t numberOfH, double p)
Definition: structure_common_cuda.hpp:42
Various utility functions to implement the tecno flux.
Definition: types.hpp:30
float p
Definition: sodshocktube.py:5
functional
Definition: FunctionalStatistics.cpp:76