Opm::cuistl Namespace Reference
Function Documentation◆ setDevice()
Sets the correct CUDA device in the setting of MPI.
◆ setZeroAtIndexSet()
The CuVector class is a simple (arithmetic) vector class for the GPU.
Example usage: #include <opm/simulators/linalg/cuistl/CuVector.hpp>
void someFunction() {
auto someDataOnCPU = std::vector<double>({1.0, 2.0, 42.0, 59.9451743, 10.7132692});
auto dataOnGPU = CuVector<double>(someDataOnCPU);
// Multiply by 4.0:
dataOnGPU *= 4.0;
// Get data back on CPU in another vector:
auto stdVectorOnCPU = dataOnGPU.asStdVector();
}
@tparam T the type to store. Can be either float, double or int.
/
template <typename T>
class CuVector
{
public:
using field_type = T;
using size_type = size_t;
CuVector(const CuVector<T>& other);
explicit CuVector(const std::vector<T>& data);
CuVector& operator=(const CuVector<T>& other);
CuVector& operator=(T scalar);
explicit CuVector(const size_t numberOfElements);
CuVector(const T* dataOnHost, const size_t numberOfElements);
virtual ~CuVector();
T* data();
const T* data() const;
template <int BlockDimension>
void copyFromHost(const Dune::BlockVector<Dune::FieldVector<T, BlockDimension>>& bvector)
{
// TODO: [perf] vector.dim() can be replaced by bvector.N() * BlockDimension
if (detail::to_size_t(m_numberOfElements) != bvector.dim()) {
OPM_THROW(std::runtime_error,
fmt::format("Given incompatible vector size. CuVector has size {}, \n"
"however, BlockVector has N() = {}, and dim = {}.",
m_numberOfElements,
bvector.N(),
bvector.dim()));
}
const auto dataPointer = static_cast<const T*>(&(bvector[0][0]));
copyFromHost(dataPointer, m_numberOfElements);
}
template <int BlockDimension>
void copyToHost(Dune::BlockVector<Dune::FieldVector<T, BlockDimension>>& bvector) const
{
// TODO: [perf] vector.dim() can be replaced by bvector.N() * BlockDimension
if (detail::to_size_t(m_numberOfElements) != bvector.dim()) {
OPM_THROW(std::runtime_error,
fmt::format("Given incompatible vector size. CuVector has size {},\n however, the BlockVector "
"has has N() = {}, and dim() = {}.",
m_numberOfElements,
bvector.N(),
bvector.dim()));
}
const auto dataPointer = static_cast<T*>(&(bvector[0][0]));
copyToHost(dataPointer, m_numberOfElements);
}
void copyFromHost(const T* dataPointer, size_t numberOfElements);
void copyToHost(T* dataPointer, size_t numberOfElements) const;
void copyFromHost(const std::vector<T>& data);
void copyToHost(std::vector<T>& data) const;
CuVector<T>& operator*=(const T& scalar);
CuVector<T>& axpy(T alpha, const CuVector<T>& y);
CuVector<T>& operator+=(const CuVector<T>& other);
CuVector<T>& operator-=(const CuVector<T>& other);
T dot(const CuVector<T>& other) const;
T two_norm() const;
T dot(const CuVector<T>& other, const CuVector<int>& indexSet, CuVector<T>& buffer) const;
T two_norm(const CuVector<int>& indexSet, CuVector<T>& buffer) const;
T dot(const CuVector<T>& other, const CuVector<int>& indexSet) const;
T two_norm(const CuVector<int>& indexSet) const;
size_type dim() const;
std::vector<T> asStdVector() const;
template <int blockSize>
Dune::BlockVector<Dune::FieldVector<T, blockSize>> asDuneBlockVector() const
{
OPM_ERROR_IF(dim() % blockSize != 0,
fmt::format("blockSize is not a multiple of dim(). Given blockSize = {}, and dim() = {}",
blockSize,
dim()));
Dune::BlockVector<Dune::FieldVector<T, blockSize>> returnValue(dim() / blockSize);
copyToHost(returnValue);
return returnValue;
}
std::size_t to_size_t(int i) to_size_t converts a (on most relevant platforms) a 32 bit signed int to a 64 bits unsigned int Definition: safe_conversion.hpp:85 ◆ toDebugString()
References Opm::to_string(). |
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||