Alsvinn  0.5.3
The fast FVM simulator with UQ support
structure_common.hpp
Go to the documentation of this file.
1 /* Copyright (c) 2018 ETH Zurich, Kjetil Olsen Lye
2  * This program is free software: you can redistribute it and/or modify
3  * it under the terms of the GNU General Public License as published by
4  * the Free Software Foundation, either version 3 of the License, or
5  * (at your option) any later version.
6  *
7  * This program is distributed in the hope that it will be useful,
8  * but WITHOUT ANY WARRANTY; without even the implied warranty of
9  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
10  * GNU General Public License for more details.
11  *
12  * You should have received a copy of the GNU General Public License
13  * along with this program. If not, see <http://www.gnu.org/licenses/>.
14  */
15 
16 #pragma once
17 #include "alsfvm/types.hpp"
18 #include "alsfvm/memory/Memory.hpp"
19 #include "alsfvm/volume/Volume.hpp"
23 
24 namespace alsfvm {
25 namespace functional {
26 
27 inline __device__ __host__ int makePositive(int position, int N) {
28  if (position < 0) {
29  position += N;
30  }
31 
32  return position;
33 
34 }
35 
36 
37 template<alsfvm::boundary::Type BoundaryType, class Function>
39  Function f,
41  int i, int j, int k, int h, int nx, int ny, int nz,
42  int ngx, int ngy, int ngz, int dimensions) {
43  const auto u = input.at(i + ngx, j + ngy, k + ngz);
44 
45  const auto numberOfCellsWithoutGhostCells = ivec3{nx, ny, nz};
46  const auto numberOfGhostCells = ivec3{ngx, ngy, ngz};
47 
48  for (int d = 0; d < dimensions; d++) {
49  // side = 0 represents bottom, side = 1 represents top
50  for (int side = 0; side < 2; side++) {
51  const bool zDir = (d == 2);
52  const bool yDir = (d == 1);
53  const bool xDir = (d == 0);
54  // Either we start on the left (i == 0), or on the right(i==1)
55  const int zStart = zDir ?
56  (side == 0 ? k - h : k + h) : (dimensions > 2 ? k - h : 0);
57 
58  const int zEnd = zDir ?
59  (zStart + 1) : (dimensions > 2 ? k + h + 1 : 1);
60 
61  int yStart = yDir ?
62  (side == 0 ? j - h : j + h) : (dimensions == 2 ? j - h + 1 :
63  (dimensions == 3 ? j - h + 1 : 0));
64 
65  int yEnd = yDir ?
66  (yStart + 1) : (dimensions == 2 ? j + h : (dimensions == 3 ? j + h : 1));
67 
68  int xStart = xDir ?
69  (side == 0 ? i - h : i + h) : i - h;
70 
71  int xEnd = xDir ?
72  (xStart + 1) : i + h + 1;
73 
74  if (zDir) {
75  xStart += 1;
76  xEnd -= 1;
77 
78  //yStart += 1;
79  //yEnd -= 1;
80  }
81 
82  for (int z = zStart; z < zEnd; z++) {
83  for (int y = yStart; y < yEnd; y++) {
84  for (int x = xStart; x < xEnd; x++) {
85  const auto discretePosition = ivec3{i, j, k};
86  const auto discretePositionPlusH = ivec3{x, y, z};
87 
88 
89  const auto u_ijk_h =
91  input,
92  discretePositionPlusH,
93  numberOfCellsWithoutGhostCells,
95 
96 
97  f(u, u_ijk_h);
98  }
99  }
100  }
101  }
102  }
103 }
104 
105 template<alsfvm::boundary::Type BoundaryType, class PowerClass>
108  output,
110  int i, int j, int k, int h, int nx, int ny, int nz,
111  int ngx, int ngy, int ngz, int dimensions, real p) {
112 
113  forEachPointInComputeStructureCube<BoundaryType>([&](double u, double u_h) {
114  output.at(h) += PowerClass::power(fabs(u - u_h), p) / (nx * ny * nz);
115  }, input, i, j, k, h, nx, ny, nz, ngx, ngy, ngz, dimensions);
116 }
117 
118 
119 
120 
121 template<alsfvm::boundary::Type BoundaryType, class PowerClass>
123  const alsfvm::volume::Volume& input, int numberOfH, double p) {
124  for (size_t var = 0; var < input.getNumberOfVariables(); ++var) {
125  auto inputView = input[var]->getView();
126  auto outputView = output[var]->getView();
127 
128  int ngx = int(input.getNumberOfXGhostCells());
129  int ngy = int(input.getNumberOfYGhostCells());
130  int ngz = int(input.getNumberOfZGhostCells());
131 
132  int nx = int(input.getNumberOfXCells());
133  int ny = int(input.getNumberOfYCells());
134  int nz = int(input.getNumberOfZCells());
135 
136  for (int k = 0; k < nz; ++k) {
137  for (int j = 0; j < ny; ++j) {
138  for (int i = 0; i < nx; ++i) {
139  for (int h = 1; h < numberOfH; ++h) {
140 
141  computeStructureCube<BoundaryType, PowerClass>(outputView, inputView, i, j, k,
142  h, nx, ny, nz,
143  ngx, ngy, ngz, input.getDimensions(), p);
144 
145  }
146  }
147  }
148  }
149 
150 
151  }
152 
153 
154 }
155 
156 template<alsfvm::boundary::Type BoundaryType>
158  const alsfvm::volume::Volume& input, int numberOfH, double p) {
159  if (p == 1.0) {
160  computeStructureCubeCPU <BoundaryType, alsutils::math::FastPower<1>>
161  (output, input, numberOfH, p);
162  } else if (p == 2.0) {
163  computeStructureCubeCPU <BoundaryType, alsutils::math::FastPower<2>>
164  (output, input, numberOfH, p);
165  }
166 
167  else if (p == 3.0) {
168  computeStructureCubeCPU <BoundaryType, alsutils::math::FastPower<3>>
169  (output, input, numberOfH, p);
170  } else if (p == 4.0) {
171  computeStructureCubeCPU <BoundaryType, alsutils::math::FastPower<4>>
172  (output, input, numberOfH, p);
173  } else if (p == 5.0) {
174  computeStructureCubeCPU <BoundaryType, alsutils::math::FastPower<5>>
175  (output, input, numberOfH, p);
176  } else {
177  computeStructureCubeCPU <BoundaryType, alsutils::math::PowPower>
178  (output, input, numberOfH, p);
179  }
180 }
181 
182 }
183 }
void computeStructureCubeCPU(alsfvm::volume::Volume &output, const alsfvm::volume::Volume &input, int numberOfH, double p)
Definition: structure_common.hpp:122
void dispatchComputeStructureCubeCPU(alsfvm::volume::Volume &output, const alsfvm::volume::Volume &input, int numberOfH, double p)
Definition: structure_common.hpp:157
__device__ __host__ void forEachPointInComputeStructureCube(Function f, const alsfvm::memory::View< const real > &input, int i, int j, int k, int h, int nx, int ny, int nz, int ngx, int ngy, int ngz, int dimensions)
Definition: structure_common.hpp:38
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
#define __host__
Definition: types.hpp:46
__device__ __host__ T & at(size_t x, size_t y, size_t z)
at returns a reference to the element at the given location
Definition: View.hpp:65
size_t nx
Definition: VolumeFactory.cpp:87
double real
Definition: types.hpp:65
virtual size_t getNumberOfZCells() const
Definition: Volume.cpp:215
virtual size_t getNumberOfXGhostCells() const
Definition: Volume.cpp:238
virtual size_t getNumberOfZGhostCells() const
Definition: Volume.cpp:254
__device__ __host__ void computeStructureCube(alsfvm::memory::View< real > &output, const alsfvm::memory::View< const real > &input, int i, int j, int k, int h, int nx, int ny, int nz, int ngx, int ngy, int ngz, int dimensions, real p)
Definition: structure_common.hpp:106
size_t ny
Definition: VolumeFactory.cpp:88
virtual size_t getNumberOfYGhostCells() const
Definition: Volume.cpp:246
__device__ __host__ int makePositive(int position, int N)
Definition: structure_common.hpp:27
virtual size_t getNumberOfXCells() const
Definition: Volume.cpp:201
Simple templated class to get the value at the boundary.
Definition: ValueAtBoundary.hpp:14
size_t numberOfGhostCells
Definition: VolumeFactory.cpp:90
virtual size_t getDimensions() const
Gets the number of space dimensions.
Definition: Volume.cpp:279
#define __device__
Definition: types.hpp:45
Definition: View.hpp:32
size_t nz
Definition: VolumeFactory.cpp:89
Various utility functions to implement the tecno flux.
Definition: types.hpp:30
float p
Definition: sodshocktube.py:5
functional
Definition: FunctionalStatistics.cpp:76