GpuBuffer.hpp
Go to the documentation of this file.
1/*
2 Copyright 2024 SINTEF AS
3
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 3 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18*/
19#ifndef OPM_GPUBUFFER_HEADER_HPP
20#define OPM_GPUBUFFER_HEADER_HPP
21#include <dune/common/fvector.hh>
22#include <dune/istl/bvector.hh>
23#include <exception>
24#include <fmt/core.h>
25#include <opm/common/ErrorMacros.hpp>
28#include <vector>
29#include <string>
30#include <cuda_runtime.h>
31
32
33namespace Opm::gpuistl
34{
35
55template <typename T>
56class GpuBuffer
57{
58public:
59 using field_type = T;
60 using size_type = size_t;
61 using value_type = T;
62
71 GpuBuffer(const GpuBuffer<T>& other);
72
82 explicit GpuBuffer(const std::vector<T>& data);
83
87 GpuBuffer() = default;
88
94 explicit GpuBuffer(const size_t numberOfElements);
95
96
106 GpuBuffer(const T* dataOnHost, const size_t numberOfElements);
107
111 virtual ~GpuBuffer();
112
116 T* data();
117
121 const T* data() const;
122
130 template <int BlockDimension>
131 void copyFromHost(const Dune::BlockVector<Dune::FieldVector<T, BlockDimension>>& bvector)
132 {
133 // TODO: [perf] vector.size() can be replaced by bvector.N() * BlockDimension
134 if (m_numberOfElements != bvector.size()) {
135 OPM_THROW(std::runtime_error,
136 fmt::format("Given incompatible vector size. GpuBuffer has size {}, \n"
137 "however, BlockVector has N() = {}, and size = {}.",
138 m_numberOfElements,
139 bvector.N(),
140 bvector.size()));
141 }
142 const auto dataPointer = static_cast<const T*>(&(bvector[0][0]));
143 copyFromHost(dataPointer, m_numberOfElements);
144 }
145
153 template <int BlockDimension>
154 void copyToHost(Dune::BlockVector<Dune::FieldVector<T, BlockDimension>>& bvector) const
155 {
156 // TODO: [perf] vector.size() can be replaced by bvector.N() * BlockDimension
157 if (m_numberOfElements != bvector.size()) {
158 OPM_THROW(std::runtime_error,
159 fmt::format("Given incompatible vector size. GpuBuffer has size {},\n however, the BlockVector "
160 "has has N() = {}, and size() = {}.",
161 m_numberOfElements,
162 bvector.N(),
163 bvector.size()));
164 }
165 const auto dataPointer = static_cast<T*>(&(bvector[0][0]));
166 copyToHost(dataPointer, m_numberOfElements);
167 }
168
176 void copyFromHost(const T* dataPointer, size_t numberOfElements);
177
185 void copyToHost(T* dataPointer, size_t numberOfElements) const;
186
194 void copyFromHost(const std::vector<T>& data);
195
203 void copyToHost(std::vector<T>& data) const;
204
209 size_type size() const;
210
215 void resize(size_t);
216
221 std::vector<T> asStdVector() const;
222
223private:
224 T* m_dataOnDevice = nullptr;
225 size_t m_numberOfElements = 0;
226
227 void assertSameSize(const GpuBuffer<T>& other) const;
228 void assertSameSize(size_t size) const;
229
230 void assertHasElements() const;
231};
232
233template <class T>
234GpuView<T> make_view(GpuBuffer<T>&);
235
236template <class T>
237GpuView<const T> make_view(const GpuBuffer<T>&);
238
239} // namespace Opm::gpuistl
240#endif
Definition: autotuner.hpp:30
PointerView< T > make_view(const std::shared_ptr< T > &ptr)
Definition: gpu_smart_pointer.hpp:271