Alsvinn  0.5.3
The fast FVM simulator with UQ support
NumericalFluxCUDA.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
21 #include "alsfvm/grid/Grid.hpp"
22 
23 namespace alsfvm {
24 namespace numflux {
25 
26 template<class Flux, class Equation, size_t dimension>
28 public:
30  alsfvm::shared_ptr<reconstruction::Reconstruction>& reconstruction,
32  alsfvm::shared_ptr<DeviceConfiguration>& deviceConfiguration
33  );
34 
51  virtual void computeFlux(const volume::Volume& conservedVariables,
52  rvec3& waveSpeed, bool computeWaveSpeed,
53  volume::Volume& output, const ivec3& start = {0, 0, 0},
54  const ivec3& end = {0, 0, 0}
55  );
56 
60  virtual size_t getNumberOfGhostCells();
61 private:
62  alsfvm::shared_ptr<reconstruction::Reconstruction> reconstruction;
63  alsfvm::shared_ptr<volume::Volume> left;
64  alsfvm::shared_ptr<volume::Volume> right;
65  alsfvm::shared_ptr<volume::Volume> fluxOutput;
66  typename Equation::Parameters equationParameters;
67  Equation equation;
68 
69 
70  void callComputeFlux(const Equation& equation,
71  const volume::Volume& conservedVariables, volume::Volume& left,
72  volume::Volume& right, volume::Volume& output, volume::Volume& temporaryOutput,
73  size_t numberOfGhostCells, rvec3& waveSpeeds,
74  reconstruction::Reconstruction& reconstruction, const ivec3& start,
75  const ivec3& end);
76 
77  template<bool xDir, bool yDir, bool zDir, size_t direction>
78  void computeFlux(const Equation& equation,
79  const volume::Volume& left,
80  const volume::Volume& right,
81  volume::Volume& output,
82  int numberOfGhostCells,
83  real& waveSpeed, const ivec3& start,
84  const ivec3& end);
85 
87  std::unique_ptr<memory::Memory<real>> waveSpeedBuffer;
88 
90  std::unique_ptr<memory::Memory<real>> waveSpeedBufferOut;
91 
93  std::unique_ptr<memory::Memory<real>>temporaryReductionMemory;
94 
96  size_t temporaryReductionMemoryStorageSizeBytes = 0;
97 
98  void initializeWaveSpeedAndReductionMemory(int size);
99 };
100 
101 } // namespace alsfvm
102 } // namespace numflux
alsfvm::shared_ptr< DeviceConfiguration > & deviceConfiguration
Definition: NumericalFluxFactory.cpp:103
NumericalFluxCUDA(const grid::Grid &grid, alsfvm::shared_ptr< reconstruction::Reconstruction > &reconstruction, const simulator::SimulatorParameters &parameters, alsfvm::shared_ptr< DeviceConfiguration > &deviceConfiguration)
Definition: Grid.hpp:27
simulator::SimulatorParameters & parameters
Definition: CellComputerFactory.cpp:60
The Volume class represents a volume (a collection of cells with values for each cell (eg...
Definition: Volume.hpp:30
Definition: NumericalFlux.hpp:25
double real
Definition: types.hpp:65
Definition: Reconstruction.hpp:21
virtual size_t getNumberOfGhostCells()
const grid::Grid & grid
Definition: NumericalFluxFactory.cpp:104
size_t numberOfGhostCells
Definition: VolumeFactory.cpp:90
virtual void computeFlux(const volume::Volume &conservedVariables, rvec3 &waveSpeed, bool computeWaveSpeed, volume::Volume &output, const ivec3 &start={0, 0, 0}, const ivec3 &end={0, 0, 0})
Definition: NumericalFluxCUDA.hpp:27
Various utility functions to implement the tecno flux.
Definition: types.hpp:30
Definition: SimulatorParameters.hpp:22
alsutils::parameters::Parameters Parameters
Definition: Parameters.hpp:24